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

#include <G4ChipsPionMinusElasticXS.hh>

Inheritance diagram for G4ChipsPionMinusElasticXS:
Collaboration diagram for G4ChipsPionMinusElasticXS:

Public Member Functions

 G4ChipsPionMinusElasticXS ()
 
 ~G4ChipsPionMinusElasticXS ()
 
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 47 of file G4ChipsPionMinusElasticXS.hh.

Constructor & Destructor Documentation

G4ChipsPionMinusElasticXS::G4ChipsPionMinusElasticXS ( )

Definition at line 59 of file G4ChipsPionMinusElasticXS.cc.

59  :G4VCrossSectionDataSet(Default_Name()), nPoints(128), nLast(nPoints-1)
60 {
61  lPMin=-8.; //Min tabulated logarithmMomentum(D)
62  lPMax= 8.; //Max tabulated logarithmMomentum(D)
63  dlnP=(lPMax-lPMin)/nLast;//LogStep inTheTable(D)
64  onlyCS=true;//Flag toCalcul OnlyCS(not Si/Bi)(L)
65  lastSIG=0.; //Last calculated cross section (L)
66  lastLP=-10.;//Last log(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 difr.Max(L)
70  theB1=0.; //TheLastSlope of 1st difruct.Max(L)
71  theS2=0.; //TheLastMantissa of 2nd difr.Max(L)
72  theB2=0.; //TheLastSlope of 2nd difruct.Max(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 difr.Max(L)
76  theB4=0.; //TheLastSlope of 4th difruct.Max(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; // Parameters ForFunctionCalculation
82  lastSST=0; // E-dep of SqardSlope 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 in CrossSection 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="")
G4ChipsPionMinusElasticXS::~G4ChipsPionMinusElasticXS ( )

Definition at line 99 of file G4ChipsPionMinusElasticXS.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 G4ChipsPionMinusElasticXS::CrossSectionDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 138 of file G4ChipsPionMinusElasticXS.cc.

139 {
140  outFile << "G4ChipsPionMinusElasticXS provides the elastic cross\n"
141  << "section for pion- 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* G4ChipsPionMinusElasticXS::Default_Name ( )
inlinestatic

Definition at line 55 of file G4ChipsPionMinusElasticXS.hh.

55 {return "ChipsPionMinusElasticXS";}

Here is the caller graph for this function:

G4double G4ChipsPionMinusElasticXS::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 167 of file G4ChipsPionMinusElasticXS.cc.

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

Definition at line 604 of file G4ChipsPionMinusElasticXS.cc.

605 {
606  static const G4double GeVSQ=gigaelectronvolt*gigaelectronvolt;
607  static const G4double third=1./3.;
608  static const G4double fifth=1./5.;
609  static const G4double sevth=1./7.;
610  if(PDG!=-211)G4cout<<"Warning*G4ChipsPionMinusElasticXS::GetExT:PDG="<<PDG<<G4endl;
611  if(onlyCS)G4cout<<"Warning*G4ChipsPionMinusElasticXS::GetExchanT:onlyCS=1"<<G4endl;
612  if(lastLP<-4.3) return lastTM*GeVSQ*G4UniformRand();// S-wave for p<14 MeV/c (kinE<.1MeV)
613  G4double q2=0.;
614  if(tgZ==1 && tgN==0) // ===> p+p=p+p
615  {
616  G4double E1=lastTM*theB1;
617  G4double R1=(1.-G4Exp(-E1));
618  G4double E2=lastTM*theB2;
619  G4double R2=(1.-G4Exp(-E2*E2*E2));
620  G4double E3=lastTM*theB3;
621  G4double R3=(1.-G4Exp(-E3));
622  G4double I1=R1*theS1/theB1;
623  G4double I2=R2*theS2;
624  G4double I3=R3*theS3;
625  G4double I12=I1+I2;
626  G4double rand=(I12+I3)*G4UniformRand();
627  if (rand<I1 )
628  {
629  G4double ran=R1*G4UniformRand();
630  if(ran>1.) ran=1.;
631  q2=-G4Log(1.-ran)/theB1;
632  }
633  else if(rand<I12)
634  {
635  G4double ran=R2*G4UniformRand();
636  if(ran>1.) ran=1.;
637  q2=-G4Log(1.-ran);
638  if(q2<0.) q2=0.;
639  q2=G4Pow::GetInstance()->powA(q2,third)/theB2;
640  }
641  else
642  {
643  G4double ran=R3*G4UniformRand();
644  if(ran>1.) ran=1.;
645  q2=-G4Log(1.-ran)/theB3;
646  }
647  }
648  else
649  {
650  G4double a=tgZ+tgN;
651  G4double E1=lastTM*(theB1+lastTM*theSS);
652  G4double R1=(1.-G4Exp(-E1));
653  G4double tss=theSS+theSS; // for future solution of quadratic equation (imediate check)
654  G4double tm2=lastTM*lastTM;
655  G4double E2=lastTM*tm2*theB2; // power 3 for lowA, 5 for HighA (1st)
656  if(a>6.5)E2*=tm2; // for heavy nuclei
657  G4double R2=(1.-G4Exp(-E2));
658  G4double E3=lastTM*theB3;
659  if(a>6.5)E3*=tm2*tm2*tm2; // power 1 for lowA, 7 (2nd) for HighA
660  G4double R3=(1.-G4Exp(-E3));
661  G4double E4=lastTM*theB4;
662  G4double R4=(1.-G4Exp(-E4));
663  G4double I1=R1*theS1;
664  G4double I2=R2*theS2;
665  G4double I3=R3*theS3;
666  G4double I4=R4*theS4;
667  G4double I12=I1+I2;
668  G4double I13=I12+I3;
669  G4double rand=(I13+I4)*G4UniformRand();
670  if(rand<I1)
671  {
672  G4double ran=R1*G4UniformRand();
673  if(ran>1.) ran=1.;
674  q2=-G4Log(1.-ran)/theB1;
675  if(std::fabs(tss)>1.e-7) q2=(std::sqrt(theB1*(theB1+(tss+tss)*q2))-theB1)/tss;
676  }
677  else if(rand<I12)
678  {
679  G4double ran=R2*G4UniformRand();
680  if(ran>1.) ran=1.;
681  q2=-G4Log(1.-ran)/theB2;
682  if(q2<0.) q2=0.;
683  if(a<6.5) q2=G4Pow::GetInstance()->powA(q2,third);
684  else q2=G4Pow::GetInstance()->powA(q2,fifth);
685  }
686  else if(rand<I13)
687  {
688  G4double ran=R3*G4UniformRand();
689  if(ran>1.) ran=1.;
690  q2=-G4Log(1.-ran)/theB3;
691  if(q2<0.) q2=0.;
692  if(a>6.5) q2=G4Pow::GetInstance()->powA(q2,sevth);
693  }
694  else
695  {
696  G4double ran=R4*G4UniformRand();
697  if(ran>1.) ran=1.;
698  q2=-G4Log(1.-ran)/theB4;
699  if(a<6.5) q2=lastTM-q2; // u reduced for lightA (starts from 0)
700  }
701  }
702  if(q2<0.) q2=0.;
703  if(!(q2>=-1.||q2<=1.)) G4cout<<"*NAN*G4QElasticCrossSect::GetExchangeT: -t="<<q2<<G4endl;
704  if(q2>lastTM)
705  {
706  q2=lastTM;
707  }
708  return q2*GeVSQ;
709 }
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 G4ChipsPionMinusElasticXS::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 156 of file G4ChipsPionMinusElasticXS.cc.

160 {
161  G4double pMom=Pt->GetTotalMomentum();
162  G4int tgN = A - tgZ;
163 
164  return GetChipsCrossSection(pMom, tgZ, tgN, -212);
165 }
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 G4ChipsPionMinusElasticXS::IsIsoApplicable ( const G4DynamicParticle Pt,
G4int  Z,
G4int  A,
const G4Element elm,
const G4Material mat 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 147 of file G4ChipsPionMinusElasticXS.cc.

150 {
151  return true;
152 }

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