#include <G4ChipsPionPlusElasticXS.hh>
Inheritance diagram for G4ChipsPionPlusElasticXS:
Public Member Functions | |
G4ChipsPionPlusElasticXS () | |
~G4ChipsPionPlusElasticXS () | |
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) |
Static Public Member Functions | |
static const char * | Default_Name () |
Definition at line 47 of file G4ChipsPionPlusElasticXS.hh.
G4ChipsPionPlusElasticXS::G4ChipsPionPlusElasticXS | ( | ) |
Definition at line 54 of file G4ChipsPionPlusElasticXS.cc.
00054 :G4VCrossSectionDataSet(Default_Name()), nPoints(128), nLast(nPoints-1) 00055 { 00056 lPMin=-8.; // Min tabulated logarithmMomentum(D) 00057 lPMax= 8.; // Max tabulated logarithmMomentum(D) 00058 dlnP=(lPMax-lPMin)/nLast;// LogStep inTheTable(D) 00059 onlyCS=true;// Flag toCalcul OnlyCS(not Si/Bi)(L) 00060 lastSIG=0.; // Last calculated cross section (L) 00061 lastLP=-10.;// Last log(mom_of IncidentHadron)(L) 00062 lastTM=0.; // Last t_maximum (L) 00063 theSS=0.; // TheLastSqSlope of 1st difr.Max(L) 00064 theS1=0.; // TheLastMantissa of 1st difr.Max(L) 00065 theB1=0.; // TheLastSlope of 1st difruct.Max(L) 00066 theS2=0.; // TheLastMantissa of 2nd difr.Max(L) 00067 theB2=0.; // TheLastSlope of 2nd difruct.Max(L) 00068 theS3=0.; // TheLastMantissa of 3d difr. Max(L) 00069 theB3=0.; // TheLastSlope of 3d difruct. Max(L) 00070 theS4=0.; // TheLastMantissa of 4th difr.Max(L) 00071 theB4=0.; // TheLastSlope of 4th difruct.Max(L) 00072 lastTZ=0; // Last atomic number of the target 00073 lastTN=0; // Last # of neutrons in the target 00074 lastPIN=0.; // Last initialized max momentum 00075 lastCST=0; // Elastic cross-section table 00076 lastPAR=0; // ParametersForFunctionalCalculation 00077 lastSST=0; // E-dep of SqaredSlope of 1st difMax 00078 lastS1T=0; // E-dep of mantissa of 1st dif.Max 00079 lastB1T=0; // E-dep of the slope of 1st difMax 00080 lastS2T=0; // E-dep of mantissa of 2nd difrMax 00081 lastB2T=0; // E-dep of the slope of 2nd difMax 00082 lastS3T=0; // E-dep of mantissa of 3d difr.Max 00083 lastB3T=0; // E-dep of the slope of 3d difrMax 00084 lastS4T=0; // E-dep of mantissa of 4th difrMax 00085 lastB4T=0; // E-dep of the slope of 4th difMax 00086 lastN=0; // The last N of calculated nucleus 00087 lastZ=0; // The last Z of calculated nucleus 00088 lastP=0.; // LastUsed in cross section Momentum 00089 lastTH=0.; // Last threshold momentum 00090 lastCS=0.; // Last value of the Cross Section 00091 lastI=0; // The last position in the DAMDB 00092 }
G4ChipsPionPlusElasticXS::~G4ChipsPionPlusElasticXS | ( | ) |
Definition at line 94 of file G4ChipsPionPlusElasticXS.cc.
00095 { 00096 std::vector<G4double*>::iterator pos; 00097 for (pos=CST.begin(); pos<CST.end(); pos++) 00098 { delete [] *pos; } 00099 CST.clear(); 00100 for (pos=PAR.begin(); pos<PAR.end(); pos++) 00101 { delete [] *pos; } 00102 PAR.clear(); 00103 for (pos=SST.begin(); pos<SST.end(); pos++) 00104 { delete [] *pos; } 00105 SST.clear(); 00106 for (pos=S1T.begin(); pos<S1T.end(); pos++) 00107 { delete [] *pos; } 00108 S1T.clear(); 00109 for (pos=B1T.begin(); pos<B1T.end(); pos++) 00110 { delete [] *pos; } 00111 B1T.clear(); 00112 for (pos=S2T.begin(); pos<S2T.end(); pos++) 00113 { delete [] *pos; } 00114 S2T.clear(); 00115 for (pos=B2T.begin(); pos<B2T.end(); pos++) 00116 { delete [] *pos; } 00117 B2T.clear(); 00118 for (pos=S3T.begin(); pos<S3T.end(); pos++) 00119 { delete [] *pos; } 00120 S3T.clear(); 00121 for (pos=B3T.begin(); pos<B3T.end(); pos++) 00122 { delete [] *pos; } 00123 B3T.clear(); 00124 for (pos=S4T.begin(); pos<S4T.end(); pos++) 00125 { delete [] *pos; } 00126 S4T.clear(); 00127 for (pos=B4T.begin(); pos<B4T.end(); pos++) 00128 { delete [] *pos; } 00129 B4T.clear(); 00130 }
static const char* G4ChipsPionPlusElasticXS::Default_Name | ( | ) | [inline, static] |
Definition at line 55 of file G4ChipsPionPlusElasticXS.hh.
Referenced by G4ChipsComponentXS::G4ChipsComponentXS(), and G4ChipsElasticModel::G4ChipsElasticModel().
G4double G4ChipsPionPlusElasticXS::GetChipsCrossSection | ( | G4double | momentum, | |
G4int | Z, | |||
G4int | N, | |||
G4int | pdg | |||
) | [virtual] |
Definition at line 155 of file G4ChipsPionPlusElasticXS.cc.
Referenced by GetIsoCrossSection(), and G4ChipsElasticModel::SampleInvariantT().
00156 { 00157 static std::vector <G4int> colN; // Vector of N for calculated nuclei (isotops) 00158 static std::vector <G4int> colZ; // Vector of Z for calculated nuclei (isotops) 00159 static std::vector <G4double> colP; // Vector of last momenta for the reaction 00160 static std::vector <G4double> colTH; // Vector of energy thresholds for the reaction 00161 static std::vector <G4double> colCS; // Vector of last cross sections for the reaction 00162 // ***---*** End of the mandatory Static Definitions of the Associative Memory ***---*** 00163 00164 G4double pEn=pMom; 00165 G4bool fCS = false; 00166 onlyCS=fCS; 00167 00168 G4bool in=false; // By default the isotope must be found in the AMDB 00169 lastP = 0.; // New momentum history (nothing to compare with) 00170 lastN = tgN; // The last N of the calculated nucleus 00171 lastZ = tgZ; // The last Z of the calculated nucleus 00172 lastI = colN.size(); // Size of the Associative Memory DB in the heap 00173 if(lastI) for(G4int i=0; i<lastI; i++) // Loop over proj/tgZ/tgN lines of DB 00174 { // The nucleus with projPDG is found in AMDB 00175 if(colN[i]==tgN && colZ[i]==tgZ) // Isotope is foind in AMDB 00176 { 00177 lastI=i; 00178 lastTH =colTH[i]; // Last THreshold (A-dependent) 00179 if(pEn<=lastTH) 00180 { 00181 return 0.; // Energy is below the Threshold value 00182 } 00183 lastP =colP [i]; // Last Momentum (A-dependent) 00184 lastCS =colCS[i]; // Last CrossSect (A-dependent) 00185 // if(std::fabs(lastP/pMom-1.)<tolerance) //VI (do not use tolerance) 00186 if(lastP == pMom) // Do not recalculate 00187 { 00188 CalculateCrossSection(fCS,-1,i,211,lastZ,lastN,pMom); // Update param's only 00189 return lastCS*millibarn; // Use theLastCS 00190 } 00191 in = true; // This is the case when the isotop is found in DB 00192 // Momentum pMom is in IU ! @@ Units 00193 lastCS=CalculateCrossSection(fCS,-1,i,211,lastZ,lastN,pMom); // read & update 00194 if(lastCS<=0. && pEn>lastTH) // Correct the threshold 00195 { 00196 lastTH=pEn; 00197 } 00198 break; // Go out of the LOOP with found lastI 00199 } 00200 } // End of attampt to find the nucleus in DB 00201 if(!in) // This nucleus has not been calculated previously 00202 { 00204 lastCS=CalculateCrossSection(fCS,0,lastI,211,lastZ,lastN,pMom);//calculate&create 00205 if(lastCS<=0.) 00206 { 00207 lastTH = 0; //ThresholdEnergy(tgZ, tgN); // The Threshold Energy which is now the last 00208 if(pEn>lastTH) 00209 { 00210 lastTH=pEn; 00211 } 00212 } 00213 colN.push_back(tgN); 00214 colZ.push_back(tgZ); 00215 colP.push_back(pMom); 00216 colTH.push_back(lastTH); 00217 colCS.push_back(lastCS); 00218 return lastCS*millibarn; 00219 } // End of creation of the new set of parameters 00220 else 00221 { 00222 colP[lastI]=pMom; 00223 colCS[lastI]=lastCS; 00224 } 00225 return lastCS*millibarn; 00226 }
Definition at line 601 of file G4ChipsPionPlusElasticXS.cc.
References G4cout, and G4UniformRand.
Referenced by G4ChipsElasticModel::SampleInvariantT().
00602 { 00603 static const G4double GeVSQ=gigaelectronvolt*gigaelectronvolt; 00604 static const G4double third=1./3.; 00605 static const G4double fifth=1./5.; 00606 static const G4double sevth=1./7.; 00607 if(PDG!= 211)G4cout<<"*Warning*G4ChipsPionPlusElasticXS::GetExT:PDG="<<PDG<<G4endl; 00608 if(onlyCS)G4cout<<"*Warning*G4ChipsPionPlusElasticXS::GetExchanT:onlyCS=1"<<G4endl; 00609 if(lastLP<-4.3) return lastTM*GeVSQ*G4UniformRand();// S-wave for p<14 MeV/c (kinE<.1MeV) 00610 G4double q2=0.; 00611 if(tgZ==1 && tgN==0) // ===> p+p=p+p 00612 { 00613 G4double E1=lastTM*theB1; 00614 G4double R1=(1.-std::exp(-E1)); 00615 G4double E2=lastTM*theB2; 00616 G4double R2=(1.-std::exp(-E2*E2*E2)); 00617 G4double E3=lastTM*theB3; 00618 G4double R3=(1.-std::exp(-E3)); 00619 G4double I1=R1*theS1/theB1; 00620 G4double I2=R2*theS2; 00621 G4double I3=R3*theS3; 00622 G4double I12=I1+I2; 00623 G4double rand=(I12+I3)*G4UniformRand(); 00624 if (rand<I1 ) 00625 { 00626 G4double ran=R1*G4UniformRand(); 00627 if(ran>1.) ran=1.; 00628 q2=-std::log(1.-ran)/theB1; 00629 } 00630 else if(rand<I12) 00631 { 00632 G4double ran=R2*G4UniformRand(); 00633 if(ran>1.) ran=1.; 00634 q2=-std::log(1.-ran); 00635 if(q2<0.) q2=0.; 00636 q2=std::pow(q2,third)/theB2; 00637 } 00638 else 00639 { 00640 G4double ran=R3*G4UniformRand(); 00641 if(ran>1.) ran=1.; 00642 q2=-std::log(1.-ran)/theB3; 00643 } 00644 } 00645 else 00646 { 00647 G4double a=tgZ+tgN; 00648 G4double E1=lastTM*(theB1+lastTM*theSS); 00649 G4double R1=(1.-std::exp(-E1)); 00650 G4double tss=theSS+theSS; // for future solution of quadratic equation (imediate check) 00651 G4double tm2=lastTM*lastTM; 00652 G4double E2=lastTM*tm2*theB2; // power 3 for lowA, 5 for HighA (1st) 00653 if(a>6.5)E2*=tm2; // for heavy nuclei 00654 G4double R2=(1.-std::exp(-E2)); 00655 G4double E3=lastTM*theB3; 00656 if(a>6.5)E3*=tm2*tm2*tm2; // power 1 for lowA, 7 (2nd) for HighA 00657 G4double R3=(1.-std::exp(-E3)); 00658 G4double E4=lastTM*theB4; 00659 G4double R4=(1.-std::exp(-E4)); 00660 G4double I1=R1*theS1; 00661 G4double I2=R2*theS2; 00662 G4double I3=R3*theS3; 00663 G4double I4=R4*theS4; 00664 G4double I12=I1+I2; 00665 G4double I13=I12+I3; 00666 G4double rand=(I13+I4)*G4UniformRand(); 00667 if(rand<I1) 00668 { 00669 G4double ran=R1*G4UniformRand(); 00670 if(ran>1.) ran=1.; 00671 q2=-std::log(1.-ran)/theB1; 00672 if(std::fabs(tss)>1.e-7) q2=(std::sqrt(theB1*(theB1+(tss+tss)*q2))-theB1)/tss; 00673 } 00674 else if(rand<I12) 00675 { 00676 G4double ran=R2*G4UniformRand(); 00677 if(ran>1.) ran=1.; 00678 q2=-std::log(1.-ran)/theB2; 00679 if(q2<0.) q2=0.; 00680 if(a<6.5) q2=std::pow(q2,third); 00681 else q2=std::pow(q2,fifth); 00682 } 00683 else if(rand<I13) 00684 { 00685 G4double ran=R3*G4UniformRand(); 00686 if(ran>1.) ran=1.; 00687 q2=-std::log(1.-ran)/theB3; 00688 if(q2<0.) q2=0.; 00689 if(a>6.5) q2=std::pow(q2,sevth); 00690 } 00691 else 00692 { 00693 G4double ran=R4*G4UniformRand(); 00694 if(ran>1.) ran=1.; 00695 q2=-std::log(1.-ran)/theB4; 00696 if(a<6.5) q2=lastTM-q2; // u reduced for lightA (starts from 0) 00697 } 00698 } 00699 if(q2<0.) q2=0.; 00700 if(!(q2>=-1.||q2<=1.)) G4cout<<"*NAN*G4QElasticCrossSect::GetExchangeT: -t="<<q2<<G4endl; 00701 if(q2>lastTM) 00702 { 00703 q2=lastTM; 00704 } 00705 return q2*GeVSQ; 00706 }
G4double G4ChipsPionPlusElasticXS::GetIsoCrossSection | ( | const G4DynamicParticle * | , | |
G4int | tgZ, | |||
G4int | A, | |||
const G4Isotope * | iso = 0 , |
|||
const G4Element * | elm = 0 , |
|||
const G4Material * | mat = 0 | |||
) | [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 143 of file G4ChipsPionPlusElasticXS.cc.
References GetChipsCrossSection(), and G4DynamicParticle::GetTotalMomentum().
00148 { 00149 G4double pMom=Pt->GetTotalMomentum(); 00150 G4int tgN = A - tgZ; 00151 00152 return GetChipsCrossSection(pMom, tgZ, tgN, 211); 00153 }
G4bool G4ChipsPionPlusElasticXS::IsIsoApplicable | ( | const G4DynamicParticle * | Pt, | |
G4int | Z, | |||
G4int | A, | |||
const G4Element * | elm, | |||
const G4Material * | mat | |||
) | [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 132 of file G4ChipsPionPlusElasticXS.cc.
References G4DynamicParticle::GetDefinition(), and G4PionPlus::PionPlus().
00135 { 00136 G4ParticleDefinition* particle = Pt->GetDefinition(); 00137 if (particle == G4PionPlus::PionPlus() ) return true; 00138 return false; 00139 }