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

#include <G4ChipsAntiBaryonElasticXS.hh>

Inheritance diagram for G4ChipsAntiBaryonElasticXS:
Collaboration diagram for G4ChipsAntiBaryonElasticXS:

Public Member Functions

 G4ChipsAntiBaryonElasticXS ()
 
 ~G4ChipsAntiBaryonElasticXS ()
 
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 G4ChipsAntiBaryonElasticXS.hh.

Constructor & Destructor Documentation

G4ChipsAntiBaryonElasticXS::G4ChipsAntiBaryonElasticXS ( )

Definition at line 59 of file G4ChipsAntiBaryonElasticXS.cc.

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

Definition at line 99 of file G4ChipsAntiBaryonElasticXS.cc.

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

Member Function Documentation

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

Reimplemented from G4VCrossSectionDataSet.

Definition at line 138 of file G4ChipsAntiBaryonElasticXS.cc.

139 {
140  outFile << "G4ChipsAntiBaryonElasticXS provides the elastic cross\n"
141  << "section for anti-baryon nucleus scattering as a function of incident\n"
142  << "momentum. The cross section is calculated using M. Kossov's\n"
143  << "CHIPS parameterization of cross section data.\n";
144 }
static const char* G4ChipsAntiBaryonElasticXS::Default_Name ( )
inlinestatic

Definition at line 55 of file G4ChipsAntiBaryonElasticXS.hh.

55 {return "ChipsAntiBaryonElasticXS";}

Here is the caller graph for this function:

G4double G4ChipsAntiBaryonElasticXS::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 206 of file G4ChipsAntiBaryonElasticXS.cc.

207 {
208  G4bool fCS = false;
209 
210  G4double pEn=pMom;
211  onlyCS=fCS;
212 
213  G4bool in=false; // By default the isotope must be found in the AMDB
214  lastP = 0.; // New momentum history (nothing to compare with)
215  lastN = tgN; // The last N of the calculated nucleus
216  lastZ = tgZ; // The last Z of the calculated nucleus
217  lastI = colN.size(); // Size of the Associative Memory DB in the heap
218  if(lastI) for(G4int i=0; i<lastI; i++) // Loop over proj/tgZ/tgN lines of DB
219  { // The nucleus with projPDG is found in AMDB
220  if(colN[i]==tgN && colZ[i]==tgZ) // Isotope is foind in AMDB
221  {
222  lastI=i;
223  lastTH =colTH[i]; // Last THreshold (A-dependent)
224  if(pEn<=lastTH)
225  {
226  return 0.; // Energy is below the Threshold value
227  }
228  lastP =colP [i]; // Last Momentum (A-dependent)
229  lastCS =colCS[i]; // Last CrossSect (A-dependent)
230  // if(std::fabs(lastP/pMom-1.)<tolerance) //VI (do not use tolerance)
231  if(lastP == pMom) // Do not recalculate
232  {
233  CalculateCrossSection(fCS,-1,i,pPDG,lastZ,lastN,pMom); // Update param's only
234  return lastCS*millibarn; // Use theLastCS
235  }
236  in = true; // This is the case when the isotop is found in DB
237  // Momentum pMom is in IU ! @@ Units
238  lastCS=CalculateCrossSection(fCS,-1,i,pPDG,lastZ,lastN,pMom); // read & update
239  if(lastCS<=0. && pEn>lastTH) // Correct the threshold
240  {
241  lastTH=pEn;
242  }
243  break; // Go out of the LOOP with found lastI
244  }
245  } // End of attampt to find the nucleus in DB
246  if(!in) // This nucleus has not been calculated previously
247  {
249  lastCS=CalculateCrossSection(fCS,0,lastI,pPDG,lastZ,lastN,pMom);//calculate&create
250  if(lastCS<=0.)
251  {
252  lastTH = 0; // ThresholdEnergy(tgZ, tgN); // The Threshold Energy which is now the last
253  if(pEn>lastTH)
254  {
255  lastTH=pEn;
256  }
257  }
258  colN.push_back(tgN);
259  colZ.push_back(tgZ);
260  colP.push_back(pMom);
261  colTH.push_back(lastTH);
262  colCS.push_back(lastCS);
263  return lastCS*millibarn;
264  } // End of creation of the new set of parameters
265  else
266  {
267  colP[lastI]=pMom;
268  colCS[lastI]=lastCS;
269  }
270  return lastCS*millibarn;
271 }
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 G4ChipsAntiBaryonElasticXS::GetExchangeT ( G4int  tZ,
G4int  tN,
G4int  pPDG 
)

Definition at line 644 of file G4ChipsAntiBaryonElasticXS.cc.

645 {
646  static const G4double GeVSQ=gigaelectronvolt*gigaelectronvolt;
647  static const G4double third=1./3.;
648  static const G4double fifth=1./5.;
649  static const G4double sevth=1./7.;
650 
651  if(PDG<-3334 || PDG>-1111)G4cout<<"*Warning*G4QAntiBaryonElCS::GetExT:PDG="<<PDG<<G4endl;
652  if(onlyCS)G4cout<<"WarningG4ChipsAntiBaryonElasticXS::GetExchanT:onlyCS=1"<<G4endl;
653  if(lastLP<-4.3) return lastTM*GeVSQ*G4UniformRand();// S-wave for p<14 MeV/c (kinE<.1MeV)
654  G4double q2=0.;
655  if(tgZ==1 && tgN==0) // ===> p+p=p+p
656  {
657  G4double E1=lastTM*theB1;
658  G4double R1=(1.-G4Exp(-E1));
659  G4double E2=lastTM*theB2;
660  G4double R2=(1.-G4Exp(-E2*E2*E2));
661  G4double E3=lastTM*theB3;
662  G4double R3=(1.-G4Exp(-E3));
663  G4double I1=R1*theS1/theB1;
664  G4double I2=R2*theS2;
665  G4double I3=R3*theS3;
666  G4double I12=I1+I2;
667  G4double rand=(I12+I3)*G4UniformRand();
668  if (rand<I1 )
669  {
670  G4double ran=R1*G4UniformRand();
671  if(ran>1.) ran=1.;
672  q2=-G4Log(1.-ran)/theB1;
673  }
674  else if(rand<I12)
675  {
676  G4double ran=R2*G4UniformRand();
677  if(ran>1.) ran=1.;
678  q2=-G4Log(1.-ran);
679  if(q2<0.) q2=0.;
680  q2=G4Pow::GetInstance()->powA(q2,third)/theB2;
681  }
682  else
683  {
684  G4double ran=R3*G4UniformRand();
685  if(ran>1.) ran=1.;
686  q2=-G4Log(1.-ran)/theB3;
687  }
688  }
689  else
690  {
691  G4double a=tgZ+tgN;
692  G4double E1=lastTM*(theB1+lastTM*theSS);
693  G4double R1=(1.-G4Exp(-E1));
694  G4double tss=theSS+theSS; // for future solution of quadratic equation (imediate check)
695  G4double tm2=lastTM*lastTM;
696  G4double E2=lastTM*tm2*theB2; // power 3 for lowA, 5 for HighA (1st)
697  if(a>6.5)E2*=tm2; // for heavy nuclei
698  G4double R2=(1.-G4Exp(-E2));
699  G4double E3=lastTM*theB3;
700  if(a>6.5)E3*=tm2*tm2*tm2; // power 1 for lowA, 7 (2nd) for HighA
701  G4double R3=(1.-G4Exp(-E3));
702  G4double E4=lastTM*theB4;
703  G4double R4=(1.-G4Exp(-E4));
704  G4double I1=R1*theS1;
705  G4double I2=R2*theS2;
706  G4double I3=R3*theS3;
707  G4double I4=R4*theS4;
708  G4double I12=I1+I2;
709  G4double I13=I12+I3;
710  G4double rand=(I13+I4)*G4UniformRand();
711  if(rand<I1)
712  {
713  G4double ran=R1*G4UniformRand();
714  if(ran>1.) ran=1.;
715  q2=-G4Log(1.-ran)/theB1;
716  if(std::fabs(tss)>1.e-7) q2=(std::sqrt(theB1*(theB1+(tss+tss)*q2))-theB1)/tss;
717  }
718  else if(rand<I12)
719  {
720  G4double ran=R2*G4UniformRand();
721  if(ran>1.) ran=1.;
722  q2=-G4Log(1.-ran)/theB2;
723  if(q2<0.) q2=0.;
724  if(a<6.5) q2=G4Pow::GetInstance()->powA(q2,third);
725  else q2=G4Pow::GetInstance()->powA(q2,fifth);
726  }
727  else if(rand<I13)
728  {
729  G4double ran=R3*G4UniformRand();
730  if(ran>1.) ran=1.;
731  q2=-G4Log(1.-ran)/theB3;
732  if(q2<0.) q2=0.;
733  if(a>6.5) q2=G4Pow::GetInstance()->powA(q2,sevth);
734  }
735  else
736  {
737  G4double ran=R4*G4UniformRand();
738  if(ran>1.) ran=1.;
739  q2=-G4Log(1.-ran)/theB4;
740  if(a<6.5) q2=lastTM-q2; // u reduced for lightA (starts from 0)
741  }
742  }
743  if(q2<0.) q2=0.;
744  if(!(q2>=-1.||q2<=1.))G4cout<<"*NAN*G4QaBElasticCrossSect::GetExchangeT:-t="<<q2<<G4endl;
745  if(q2>lastTM)
746  {
747  q2=lastTM;
748  }
749  return q2*GeVSQ;
750 }
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:

Here is the caller graph for this function:

G4double G4ChipsAntiBaryonElasticXS::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 194 of file G4ChipsAntiBaryonElasticXS.cc.

198 {
199  G4double pMom=Pt->GetTotalMomentum();
200  G4int tgN = A - tgZ;
201  G4int pdg = Pt->GetDefinition()->GetPDGEncoding();
202 
203  return GetChipsCrossSection(pMom, tgZ, tgN, pdg);
204 }
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 G4ChipsAntiBaryonElasticXS::IsIsoApplicable ( const G4DynamicParticle Pt,
G4int  Z,
G4int  A,
const G4Element elm,
const G4Material mat 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 146 of file G4ChipsAntiBaryonElasticXS.cc.

149 {
150 
151  /*
152  if(particle == G4AntiNeutron::AntiNeutron())
153  {
154  return true;
155  }
156  else if(particle == G4AntiProton::AntiProton())
157  {
158  return true;
159  }
160  else if(particle == G4AntiLambda::AntiLambda())
161  {
162  return true;
163  }
164  else if(particle == G4AntiSigmaPlus::AntiSigmaPlus())
165  {
166  return true;
167  }
168  else if(particle == G4AntiSigmaMinus::AntiSigmaMinus())
169  {
170  return true;
171  }
172  else if(particle == G4AntiSigmaZero::AntiSigmaZero())
173  {
174  return true;
175  }
176  else if(particle == G4AntiXiMinus::AntiXiMinus())
177  {
178  return true;
179  }
180  else if(particle == G4AntiXiZero::AntiXiZero())
181  {
182  return true;
183  }
184  else if(particle == G4AntiOmegaMinus::AntiOmegaMinus())
185  {
186  return true;
187  }
188  */
189  return true;
190 }

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