Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4ChipsHyperonElasticXS Class Reference

#include <G4ChipsHyperonElasticXS.hh>

Inheritance diagram for G4ChipsHyperonElasticXS:
Collaboration diagram for G4ChipsHyperonElasticXS:

Public Member Functions

 G4ChipsHyperonElasticXS ()
 
 ~G4ChipsHyperonElasticXS ()
 
virtual void CrossSectionDescription (std::ostream &) const
 
virtual G4bool IsIsoApplicable (const G4DynamicParticle *Pt, G4int Z, G4int A, const G4Element *elm, const G4Material *mat)
 
virtual G4double GetIsoCrossSection (const G4DynamicParticle *, G4int tgZ, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
 
virtual G4double GetChipsCrossSection (G4double momentum, G4int Z, G4int N, G4int pdg)
 
G4double GetExchangeT (G4int tZ, G4int tN, G4int pPDG)
 
- Public Member Functions inherited from G4VCrossSectionDataSet
 G4VCrossSectionDataSet (const G4String &nam="")
 
virtual ~G4VCrossSectionDataSet ()
 
virtual G4bool IsElementApplicable (const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
 
G4double GetCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=0)
 
G4double ComputeCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=0)
 
virtual G4double GetElementCrossSection (const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
 
virtual G4IsotopeSelectIsotope (const G4Element *, G4double kinEnergy)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void DumpPhysicsTable (const G4ParticleDefinition &)
 
virtual G4int GetVerboseLevel () const
 
virtual void SetVerboseLevel (G4int value)
 
G4double GetMinKinEnergy () const
 
void SetMinKinEnergy (G4double value)
 
G4double GetMaxKinEnergy () const
 
void SetMaxKinEnergy (G4double value)
 
const G4StringGetName () const
 

Static Public Member Functions

static const char * Default_Name ()
 

Additional Inherited Members

- Protected Member Functions inherited from G4VCrossSectionDataSet
void SetName (const G4String &)
 
- Protected Attributes inherited from G4VCrossSectionDataSet
G4int verboseLevel
 

Detailed Description

Definition at line 46 of file G4ChipsHyperonElasticXS.hh.

Constructor & Destructor Documentation

G4ChipsHyperonElasticXS::G4ChipsHyperonElasticXS ( )

Definition at line 64 of file G4ChipsHyperonElasticXS.cc.

64  :G4VCrossSectionDataSet(Default_Name()), nPoints(128), nLast(nPoints-1)
65 {
66  lPMin=-8.; //Min tabulatedLogarithmMomentum(D)
67  lPMax= 8.; //Max tabulatedLogarithmMomentum(D)
68  dlnP=(lPMax-lPMin)/nLast;// LogStep inTable (D)
69  onlyCS=true;//Flag toCalculOnlyCS(not Si/Bi)(L)
70  lastSIG=0.; //Last calculated cross section (L)
71  lastLP=-10.;//LastLog(mom_of IncidentHadron)(L)
72  lastTM=0.; //Last t_maximum (L)
73  theSS=0.; //TheLastSqSlope of 1st difr.Max(L)
74  theS1=0.; //TheLastMantissa of 1st difrMax(L)
75  theB1=0.; //TheLastSlope of 1st difructMax(L)
76  theS2=0.; //TheLastMantissa of 2nd difrMax(L)
77  theB2=0.; //TheLastSlope of 2nd difructMax(L)
78  theS3=0.; //TheLastMantissa of 3d difr.Max(L)
79  theB3=0.; //TheLastSlope of 3d difruct.Max(L)
80  theS4=0.; //TheLastMantissa of 4th difrMax(L)
81  theB4=0.; //TheLastSlope of 4th difructMax(L)
82  lastTZ=0; // Last atomic number of the target
83  lastTN=0; // Last # of neutrons in the target
84  lastPIN=0.; // Last initialized max momentum
85  lastCST=0; // Elastic cross-section table
86  lastPAR=0; // ParametersForFunctionCalculation
87  lastSST=0; // E-dep ofSqardSlope of 1st difMax
88  lastS1T=0; // E-dep of mantissa of 1st dif.Max
89  lastB1T=0; // E-dep of the slope of 1st difMax
90  lastS2T=0; // E-dep of mantissa of 2nd difrMax
91  lastB2T=0; // E-dep of the slope of 2nd difMax
92  lastS3T=0; // E-dep of mantissa of 3d difr.Max
93  lastB3T=0; // E-dep of the slope of 3d difrMax
94  lastS4T=0; // E-dep of mantissa of 4th difrMax
95  lastB4T=0; // E-dep of the slope of 4th difMax
96  lastN=0; // The last N of calculated nucleus
97  lastZ=0; // The last Z of calculated nucleus
98  lastP=0.; // LastUsed inCrossSection Momentum
99  lastTH=0.; // Last threshold momentum
100  lastCS=0.; // Last value of the Cross Section
101  lastI=0; // The last position in the DAMDB
102 }
G4VCrossSectionDataSet(const G4String &nam="")
static const char * Default_Name()
G4ChipsHyperonElasticXS::~G4ChipsHyperonElasticXS ( )

Definition at line 104 of file G4ChipsHyperonElasticXS.cc.

105 {
106  std::vector<G4double*>::iterator pos;
107  for (pos=CST.begin(); pos<CST.end(); pos++)
108  { delete [] *pos; }
109  CST.clear();
110  for (pos=PAR.begin(); pos<PAR.end(); pos++)
111  { delete [] *pos; }
112  PAR.clear();
113  for (pos=SST.begin(); pos<SST.end(); pos++)
114  { delete [] *pos; }
115  SST.clear();
116  for (pos=S1T.begin(); pos<S1T.end(); pos++)
117  { delete [] *pos; }
118  S1T.clear();
119  for (pos=B1T.begin(); pos<B1T.end(); pos++)
120  { delete [] *pos; }
121  B1T.clear();
122  for (pos=S2T.begin(); pos<S2T.end(); pos++)
123  { delete [] *pos; }
124  S2T.clear();
125  for (pos=B2T.begin(); pos<B2T.end(); pos++)
126  { delete [] *pos; }
127  B2T.clear();
128  for (pos=S3T.begin(); pos<S3T.end(); pos++)
129  { delete [] *pos; }
130  S3T.clear();
131  for (pos=B3T.begin(); pos<B3T.end(); pos++)
132  { delete [] *pos; }
133  B3T.clear();
134  for (pos=S4T.begin(); pos<S4T.end(); pos++)
135  { delete [] *pos; }
136  S4T.clear();
137  for (pos=B4T.begin(); pos<B4T.end(); pos++)
138  { delete [] *pos; }
139  B4T.clear();
140 }
static const G4double pos

Member Function Documentation

void G4ChipsHyperonElasticXS::CrossSectionDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 143 of file G4ChipsHyperonElasticXS.cc.

144 {
145  outFile << "G4ChipsHyperonElasticXS provides the elastic cross\n"
146  << "section for hyperon nucleus scattering as a function of incident\n"
147  << "momentum. The cross section is calculated using M. Kossov's\n"
148  << "CHIPS parameterization of cross section data.\n";
149 }
static const char* G4ChipsHyperonElasticXS::Default_Name ( )
inlinestatic

Definition at line 54 of file G4ChipsHyperonElasticXS.hh.

54 {return "ChipsHyperonElasticXS";}

Here is the caller graph for this function:

G4double G4ChipsHyperonElasticXS::GetChipsCrossSection ( G4double  momentum,
G4int  Z,
G4int  N,
G4int  pdg 
)
virtual

!The slave functions must provide cross-sections in millibarns (mb) !! (not in IU)

Definition at line 172 of file G4ChipsHyperonElasticXS.cc.

173 {
174 
175  G4bool fCS = false;
176  G4double pEn=pMom;
177 
178  onlyCS=fCS;
179 
180  G4bool in=false; // By default the isotope must be found in the AMDB
181  lastP = 0.; // New momentum history (nothing to compare with)
182  lastN = tgN; // The last N of the calculated nucleus
183  lastZ = tgZ; // The last Z of the calculated nucleus
184  lastI = colN.size(); // Size of the Associative Memory DB in the heap
185  if(lastI) for(G4int i=0; i<lastI; i++) // Loop over proj/tgZ/tgN lines of DB
186  { // The nucleus with projPDG is found in AMDB
187  if(colN[i]==tgN && colZ[i]==tgZ) // Isotope is foind in AMDB
188  {
189  lastI=i;
190  lastTH =colTH[i]; // Last THreshold (A-dependent)
191  if(pEn<=lastTH)
192  {
193  return 0.; // Energy is below the Threshold value
194  }
195  lastP =colP [i]; // Last Momentum (A-dependent)
196  lastCS =colCS[i]; // Last CrossSect (A-dependent)
197  // if(std::fabs(lastP/pMom-1.)<tolerance) //VI (do not use tolerance)
198  if(lastP == pMom) // Do not recalculate
199  {
200  CalculateCrossSection(fCS,-1,i,pPDG,lastZ,lastN,pMom); // Update param's only
201  return lastCS*millibarn; // Use theLastCS
202  }
203  in = true; // This is the case when the isotop is found in DB
204  // Momentum pMom is in IU ! @@ Units
205  lastCS=CalculateCrossSection(fCS,-1,i,pPDG,lastZ,lastN,pMom); // read & update
206  if(lastCS<=0. && pEn>lastTH) // Correct the threshold
207  {
208  lastTH=pEn;
209  }
210  break; // Go out of the LOOP with found lastI
211  }
212  } // End of attampt to find the nucleus in DB
213  if(!in) // This nucleus has not been calculated previously
214  {
216  lastCS=CalculateCrossSection(fCS,0,lastI,pPDG,lastZ,lastN,pMom);//calculate&create
217  if(lastCS<=0.)
218  {
219  lastTH = 0; //ThresholdEnergy(tgZ, tgN); // The Threshold Energy which is now the last
220  if(pEn>lastTH)
221  {
222  lastTH=pEn;
223  }
224  }
225  colN.push_back(tgN);
226  colZ.push_back(tgZ);
227  colP.push_back(pMom);
228  colTH.push_back(lastTH);
229  colCS.push_back(lastCS);
230  return lastCS*millibarn;
231  } // End of creation of the new set of parameters
232  else
233  {
234  colP[lastI]=pMom;
235  colCS[lastI]=lastCS;
236  }
237  return lastCS*millibarn;
238 }
int G4int
Definition: G4Types.hh:78
bool G4bool
Definition: G4Types.hh:79
double G4double
Definition: G4Types.hh:76
static constexpr double millibarn
Definition: G4SIunits.hh:106

Here is the caller graph for this function:

G4double G4ChipsHyperonElasticXS::GetExchangeT ( G4int  tZ,
G4int  tN,
G4int  pPDG 
)

Definition at line 610 of file G4ChipsHyperonElasticXS.cc.

611 {
612  static const G4double GeVSQ=gigaelectronvolt*gigaelectronvolt;
613  static const G4double third=1./3.;
614  static const G4double fifth=1./5.;
615  static const G4double sevth=1./7.;
616  //AR-04Jun2014 if(PDG==3222 || PDG<3000 || PDG>3334)G4cout<<"*Warning*G4QHyElCS::GET:PDG="<<PDG<<G4endl;
617  if(PDG<3000 || PDG>3334)G4cout<<"*Warning*G4QHyElCS::GET:PDG="<<PDG<<G4endl;
618  if(onlyCS)G4cout<<"*Warning*G4ChipsHyperonElasticXS::GetExchanT: onlyCS=1"<<G4endl;
619  if(lastLP<-4.3) return lastTM*GeVSQ*G4UniformRand();// S-wave for p<14 MeV/c (kinE<.1MeV)
620  G4double q2=0.;
621  if(tgZ==1 && tgN==0) // ===> p+p=p+p
622  {
623  G4double E1=lastTM*theB1;
624  G4double R1=(1.-G4Exp(-E1));
625  G4double E2=lastTM*theB2;
626  G4double R2=(1.-G4Exp(-E2*E2*E2));
627  G4double E3=lastTM*theB3;
628  G4double R3=(1.-G4Exp(-E3));
629  G4double I1=R1*theS1/theB1;
630  G4double I2=R2*theS2;
631  G4double I3=R3*theS3;
632  G4double I12=I1+I2;
633  G4double rand=(I12+I3)*G4UniformRand();
634  if (rand<I1 )
635  {
636  G4double ran=R1*G4UniformRand();
637  if(ran>1.) ran=1.;
638  q2=-G4Log(1.-ran)/theB1;
639  }
640  else if(rand<I12)
641  {
642  G4double ran=R2*G4UniformRand();
643  if(ran>1.) ran=1.;
644  q2=-G4Log(1.-ran);
645  if(q2<0.) q2=0.;
646  q2=G4Pow::GetInstance()->powA(q2,third)/theB2;
647  }
648  else
649  {
650  G4double ran=R3*G4UniformRand();
651  if(ran>1.) ran=1.;
652  q2=-G4Log(1.-ran)/theB3;
653  }
654  }
655  else
656  {
657  G4double a=tgZ+tgN;
658  G4double E1=lastTM*(theB1+lastTM*theSS);
659  G4double R1=(1.-G4Exp(-E1));
660  G4double tss=theSS+theSS; // for future solution of quadratic equation (imediate check)
661  G4double tm2=lastTM*lastTM;
662  G4double E2=lastTM*tm2*theB2; // power 3 for lowA, 5 for HighA (1st)
663  if(a>6.5)E2*=tm2; // for heavy nuclei
664  G4double R2=(1.-G4Exp(-E2));
665  G4double E3=lastTM*theB3;
666  if(a>6.5)E3*=tm2*tm2*tm2; // power 1 for lowA, 7 (2nd) for HighA
667  G4double R3=(1.-G4Exp(-E3));
668  G4double E4=lastTM*theB4;
669  G4double R4=(1.-G4Exp(-E4));
670  G4double I1=R1*theS1;
671  G4double I2=R2*theS2;
672  G4double I3=R3*theS3;
673  G4double I4=R4*theS4;
674  G4double I12=I1+I2;
675  G4double I13=I12+I3;
676  G4double rand=(I13+I4)*G4UniformRand();
677  if(rand<I1)
678  {
679  G4double ran=R1*G4UniformRand();
680  if(ran>1.) ran=1.;
681  q2=-G4Log(1.-ran)/theB1;
682  if(std::fabs(tss)>1.e-7) q2=(std::sqrt(theB1*(theB1+(tss+tss)*q2))-theB1)/tss;
683  }
684  else if(rand<I12)
685  {
686  G4double ran=R2*G4UniformRand();
687  if(ran>1.) ran=1.;
688  q2=-G4Log(1.-ran)/theB2;
689  if(q2<0.) q2=0.;
690  if(a<6.5) q2=G4Pow::GetInstance()->powA(q2,third);
691  else q2=G4Pow::GetInstance()->powA(q2,fifth);
692  }
693  else if(rand<I13)
694  {
695  G4double ran=R3*G4UniformRand();
696  if(ran>1.) ran=1.;
697  q2=-G4Log(1.-ran)/theB3;
698  if(q2<0.) q2=0.;
699  if(a>6.5) q2=G4Pow::GetInstance()->powA(q2,sevth);
700  }
701  else
702  {
703  G4double ran=R4*G4UniformRand();
704  if(ran>1.) ran=1.;
705  q2=-G4Log(1.-ran)/theB4;
706  if(a<6.5) q2=lastTM-q2; // u reduced for lightA (starts from 0)
707  }
708  }
709  if(q2<0.) q2=0.;
710  if(!(q2>=-1.||q2<=1.))G4cout<<"*NAN*G4QHyElasticCrossSect::GetExchangeT:-t="<<q2<<G4endl;
711  if(q2>lastTM)
712  {
713  q2=lastTM;
714  }
715  return q2*GeVSQ;
716 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
static constexpr double gigaelectronvolt
Definition: G4SIunits.hh:209
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4ChipsHyperonElasticXS::GetIsoCrossSection ( const G4DynamicParticle Pt,
G4int  tgZ,
G4int  A,
const G4Isotope iso = 0,
const G4Element elm = 0,
const G4Material mat = 0 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 160 of file G4ChipsHyperonElasticXS.cc.

164 {
165  G4double pMom=Pt->GetTotalMomentum();
166  G4int tgN = A - tgZ;
167  G4int pdg = Pt->GetDefinition()->GetPDGEncoding();
168 
169  return GetChipsCrossSection(pMom, tgZ, tgN, pdg);
170 }
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
G4double GetTotalMomentum() const
double A(double temperature)
virtual G4double GetChipsCrossSection(G4double momentum, G4int Z, G4int N, G4int pdg)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4bool G4ChipsHyperonElasticXS::IsIsoApplicable ( const G4DynamicParticle Pt,
G4int  Z,
G4int  A,
const G4Element elm,
const G4Material mat 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 151 of file G4ChipsHyperonElasticXS.cc.

154 {
155  return true;
156 }

The documentation for this class was generated from the following files: