#include <G4QNeutronNuclearCrossSection.hh>
Inheritance diagram for G4QNeutronNuclearCrossSection:
Public Member Functions | |
~G4QNeutronNuclearCrossSection () | |
virtual G4double | GetCrossSection (G4bool fCS, G4double pMom, G4int tgZ, G4int tgN, G4int pPDG=2112) |
G4double | CalculateCrossSection (G4bool CS, G4int F, G4int I, G4int PDG, G4int Z, G4int N, G4double Momentum) |
Static Public Member Functions | |
static G4VQCrossSection * | GetPointer () |
Protected Member Functions | |
G4QNeutronNuclearCrossSection () |
Definition at line 48 of file G4QNeutronNuclearCrossSection.hh.
G4QNeutronNuclearCrossSection::G4QNeutronNuclearCrossSection | ( | ) | [inline, protected] |
G4QNeutronNuclearCrossSection::~G4QNeutronNuclearCrossSection | ( | ) |
G4double G4QNeutronNuclearCrossSection::CalculateCrossSection | ( | G4bool | CS, | |
G4int | F, | |||
G4int | I, | |||
G4int | PDG, | |||
G4int | Z, | |||
G4int | N, | |||
G4double | Momentum | |||
) | [virtual] |
Implements G4VQCrossSection.
Definition at line 232 of file G4QNeutronNuclearCrossSection.cc.
References G4VQCrossSection::EquLinearFit(), G4cerr, G4cout, G4endl, and CLHEP::detail::n.
Referenced by GetCrossSection().
00234 { 00235 static const G4double THmin=1.; // default minimum Momentum (MeV/c) Threshold 00236 static const G4double THmiG=THmin*.001; // minimum Momentum (GeV/c) Threshold 00237 static const G4double dP=10.; // step for the LEN (Low ENergy) table MeV/c 00238 static const G4double dPG=dP*.001; // step for the LEN (Low ENergy) table GeV/c 00239 static const G4int nL=105; // A#of LEN points in E (step 10 MeV/c) 00240 static const G4double Pmin=THmin+(nL-1)*dP; // minP for the HighE part with safety 00241 static const G4double Pmax=227000.; // maxP for the HEN (High ENergy) part 227 GeV 00242 static const G4int nH=224; // A#of HEN points in lnE 00243 static const G4double milP=std::log(Pmin);// Low logarithm energy for the HEN part 00244 static const G4double malP=std::log(Pmax);// High logarithm energy (each 2.75 percent) 00245 static const G4double dlP=(malP-milP)/(nH-1); // Step in log energy in the HEN part 00246 static const G4double milPG=std::log(.001*Pmin);// Low logarithmEnergy for HEN part GeV/c 00247 // 00248 // Associative memory for acceleration 00249 //static std::vector <G4double> spA; // shadowing coefficients (A-dependent) 00250 #ifdef debug 00251 G4cout<<"G4QNeutNCS::CalCS:N="<<targN<<",Z="<<targZ<<",P="<<Momentum<<">"<<THmin<<G4endl; 00252 #endif 00253 G4double sigma=0.; 00254 if(F&&I) sigma=0.; // @@ *!* Fake line *!* to use F & I !!!Temporary!!! 00255 //G4double A=targN+targZ; // A of the target 00256 #ifdef debug 00257 G4cout<<"G4QNeutNucCS::CalCS: A="<<A<<",F="<<F<<",I="<<I<<",nL="<<nL<<",nH="<<nH<<G4endl; 00258 #endif 00259 if(F<=0) // This isotope was not the last used isotop 00260 { 00261 if(F<0) // This isotope was found in DAMDB =-----=> RETRIEVE 00262 { 00263 G4int sync=LEN->size(); 00264 if(sync<=I) G4cerr<<"*!*G4QNetronNuclCS::CalcCrossSect:Sync="<<sync<<"<="<<I<<G4endl; 00265 lastLEN=(*LEN)[I]; // Pointer to prepared LowEnergy cross sections 00266 lastHEN=(*HEN)[I]; // Pointer to prepared High Energy cross sections 00267 } 00268 else // This isotope wasn't calculated before => CREATE 00269 { 00270 lastLEN = new G4double[nL]; // Allocate memory for the new LEN cross sections 00271 lastHEN = new G4double[nH]; // Allocate memory for the new HEN cross sections 00272 // --- Instead of making a separate function --- 00273 G4double P=THmiG; // Table threshold in GeV/c 00274 for(G4int n=0; n<nL; n++) 00275 { 00276 lastLEN[n] = CrossSectionLin(targZ, targN, P); 00277 #ifdef debug 00278 G4cout<<"-Li->G4QNeutNucCS::CalcCS: P="<<P<<", S="<<lastLEN[n]<<G4endl; 00279 #endif 00280 P+=dPG; 00281 } 00282 G4double lP=milPG; 00283 for(G4int n=0; n<nH; n++) 00284 { 00285 lastHEN[n] = CrossSectionLog(targZ, targN, lP); 00286 #ifdef debug 00287 G4cout<<"-Li->G4QNeutNucCS::CalcCS: lP="<<lP<<", S="<<lastHEN[n]<<G4endl; 00288 #endif 00289 lP+=dlP; 00290 } 00291 // --- End of possible separate function 00292 // *** The synchronization check *** 00293 G4int sync=LEN->size(); 00294 #ifdef debug 00295 G4cout<<"-*->G4QNeutNCS::CCS:Z="<<targZ<<",N="<<targN<<",I="<<I<<",S="<<sync<<G4endl; 00296 #endif 00297 if(sync!=I) 00298 { 00299 G4cerr<<"***G4QNetronNuclearCS::CalcCrossSect: Sync="<<sync<<"#"<<I<<", Z=" <<targZ 00300 <<", N="<<targN<<", F="<<F<<G4endl; 00301 //G4Exception("G4ProtonNuclearCS::CalculateCS:","39",FatalException,"overflow DB"); 00302 } 00303 LEN->push_back(lastLEN); // remember the Low Energy Table 00304 HEN->push_back(lastHEN); // remember the High Energy Table 00305 } // End of creation of the new set of parameters 00306 } // End of parameters udate 00307 // =------------------= NOW the Magic Formula =---------------------------= 00308 #ifdef debug 00309 G4cout<<"G4QNeutCS::CalcCS:lTH="<<lastTH<<",Pi="<<Pmin<<",dP="<<dP<<",dlP="<<dlP<<G4endl; 00310 #endif 00311 if (Momentum<lastTH) return 0.; // It must be already checked in the interface class 00312 else if (Momentum<Pmin) // High Energy region 00313 { 00314 #ifdef debug 00315 G4cout<<"G4QNeutCS::CalcCS:bLEN nL="<<nL<<",TH="<<THmin<<",dP="<<dP<<G4endl; 00316 #endif 00317 sigma=EquLinearFit(Momentum,nL,THmin,dP,lastLEN); 00318 #ifdef debugn 00319 if(sigma<0.) 00320 G4cout<<"G4QNeutCS::CalcCS: E="<<Momentum<<",T="<<THmin<<",dP="<<dP<<G4endl; 00321 #endif 00322 } 00323 else if (Momentum<Pmax) // High Energy region 00324 { 00325 G4double lP=std::log(Momentum); 00326 #ifdef debug 00327 G4cout<<"G4QProtNucCS::CalcCS: before HEN nH="<<nH<<",iE="<<milP<<",dlP="<<dlP<<G4endl; 00328 #endif 00329 sigma=EquLinearFit(lP,nH,milP,dlP,lastHEN); 00330 } 00331 else // UHE region (calculation, not frequent) 00332 { 00333 G4double P=0.001*Momentum; // Approximation formula is for P in GeV/c 00334 sigma=CrossSectionFormula(targZ, targN, P, std::log(P)); 00335 } 00336 #ifdef debug 00337 G4cout<<"G4QNeutronNuclearCrossSection::CalcCS: CS="<<sigma<<G4endl; 00338 #endif 00339 if(sigma<0.) return 0.; 00340 return sigma; 00341 }
G4double G4QNeutronNuclearCrossSection::GetCrossSection | ( | G4bool | fCS, | |
G4double | pMom, | |||
G4int | tgZ, | |||
G4int | tgN, | |||
G4int | pPDG = 2112 | |||
) | [virtual] |
Reimplemented from G4VQCrossSection.
Definition at line 82 of file G4QNeutronNuclearCrossSection.cc.
References CalculateCrossSection(), G4cout, G4endl, G4VQCrossSection::ThresholdEnergy(), and G4VQCrossSection::tolerance.
00084 { 00085 //A.R.23-Oct-2012 Shadowed variable static G4double tolerance=0.001; // Tolerance (0.1%) to consider as "the same mom" 00086 static G4int j; // A#0f Z/N-records already tested in AMDB 00087 static std::vector <G4int> colN; // Vector of N for calculated nuclei (isotops) 00088 static std::vector <G4int> colZ; // Vector of Z for calculated nuclei (isotops) 00089 static std::vector <G4double> colP; // Vector of last momenta for the reaction 00090 static std::vector <G4double> colTH; // Vector of energy thresholds for the reaction 00091 static std::vector <G4double> colCS; // Vector of last cross sections for the reaction 00092 // ***---*** End of the mandatory Static Definitions of the Associative Memory ***---*** 00093 #ifdef debug 00094 G4cout<<"G4QNeutCS::GetCS:>>> f="<<fCS<<", p="<<pMom<<", Z="<<tgZ<<"("<<lastZ<<") ,N=" 00095 <<tgN<<"("<<lastN<<"),PDG=2112, thresh="<<lastTH<<",Sz="<<colN.size()<<G4endl; 00096 #endif 00097 if(PDG!=2112) G4cout<<"-Warning-G4QNeutronCS::GetCS:**Not a neutron**,PDG="<<PDG<<G4endl; 00098 G4bool in=false; // By default the isotope must be found in the AMDB 00099 if(tgN!=lastN || tgZ!=lastZ) // The nucleus was not the last used isotope 00100 { 00101 in = false; // By default the isotope haven't be found in AMDB 00102 lastP = 0.; // New momentum history (nothing to compare with) 00103 lastN = tgN; // The last N of the calculated nucleus 00104 lastZ = tgZ; // The last Z of the calculated nucleus 00105 lastI = colN.size(); // Size of the Associative Memory DB in the heap 00106 j = 0; // A#0f records found in DB for this projectile 00107 #ifdef debug 00108 G4cout<<"G4QNeutCS::GetCS: the amount of records in the AMDB lastI="<<lastI<<G4endl; 00109 #endif 00110 if(lastI) for(G4int i=0; i<lastI; i++) // AMDB exists, try to find the (Z,N) isotope 00111 { 00112 if(colN[i]==tgN && colZ[i]==tgZ) // Try the record "i" in the AMDB 00113 { 00114 lastI=i; // Remember the index for future fast/last use 00115 lastTH =colTH[i]; // The last THreshold (A-dependent) 00116 #ifdef debug 00117 G4cout<<"G4QNCS::GetCS: *Found* P="<<pMom<<",Threshold="<<lastTH<<",j="<<j<<G4endl; 00118 #endif 00119 if(pMom<=lastTH) 00120 { 00121 #ifdef debug 00122 G4cout<<"G4QNCS::GetCS:Found,P="<<pMom<<" < Threshold="<<lastTH<<",CS=0"<<G4endl; 00123 #endif 00124 return 0.; // Energy is below the Threshold value 00125 } 00126 lastP =colP [i]; // Last Momentum (A-dependent) 00127 lastCS =colCS[i]; // Last CrossSect (A-dependent) 00128 if(std::fabs(lastP-pMom)<tolerance*pMom) 00129 //if(lastP==pMom) // VI do not use tolerance 00130 { 00131 #ifdef debug 00132 G4cout<<"...G4QNCS::GetCS:.DoNothing.P="<<pMom<<",CS="<<lastCS*millibarn<<G4endl; 00133 #endif 00134 //CalculateCrossSection(fCS,-1,j,2112,lastZ,lastN,pMom); // Update param's only 00135 return lastCS*millibarn; // Use theLastCS 00136 } 00137 in = true; // This is the case when the isotop is found in DB 00138 // Momentum pMom is in IU ! @@ Units 00139 #ifdef debug 00140 G4cout<<"G4QNCS::GX:UpdatDB P="<<pMom<<",f="<<fCS<<",lI="<<lastI<<",j="<<j<<G4endl; 00141 #endif 00142 lastCS=CalculateCrossSection(fCS,-1,j,2112,lastZ,lastN,pMom); // read & update 00143 #ifdef debug 00144 G4cout<<"G4QNeutCS::GetCrosSec: *****> New (inDB) Calculated CS="<<lastCS<<G4endl; 00145 #endif 00146 if(lastCS<=0. && pMom>lastTH) // Correct the threshold (@@ No intermediate Zeros) 00147 { 00148 #ifdef debug 00149 G4cout<<"G4QNeutCS::GetCS: New P="<<pMom<<"(CS=0) > Threshold="<<lastTH<<G4endl; 00150 #endif 00151 lastCS=0.; 00152 lastTH=pMom; 00153 } 00154 break; // Go out of the LOOP 00155 } 00156 #ifdef debug 00157 G4cout<<"-->G4QNeutronNuclearCrossSec::GetCrosSec: pPDG=2112, j="<<j<<", N="<<colN[i] 00158 <<",Z["<<i<<"]="<<colZ[i]<<G4endl; 00159 #endif 00160 j++; // Increment a#0f records found in DB 00161 } 00162 #ifdef debug 00163 G4cout<<"-?-G4QNeuCS::GetCS:RC Z="<<tgZ<<",N="<<tgN<<",in="<<in<<",j="<<j<<"?"<<G4endl; 00164 #endif 00165 if(!in) // This isotope has not been calculated previously 00166 { 00167 #ifdef debug 00168 G4cout<<"^^^G4QNeutCS::GetCS:CalcNew P="<<pMom<<",f="<<fCS<<",lastI="<<lastI<<G4endl; 00169 #endif 00171 lastCS=CalculateCrossSection(fCS,0,j,2112,lastZ,lastN,pMom); //calculate & create 00172 //if(lastCS>0.) // It means that the AMBD was initialized 00173 //{ 00174 00175 lastTH = ThresholdEnergy(tgZ, tgN); // The Threshold Energy which is now the last 00176 #ifdef debug 00177 G4cout<<"G4QPrCrossSection::GetCrossSect: NewThresh="<<lastTH<<",P="<<pMom<<G4endl; 00178 #endif 00179 colN.push_back(tgN); 00180 colZ.push_back(tgZ); 00181 colP.push_back(pMom); 00182 colTH.push_back(lastTH); 00183 colCS.push_back(lastCS); 00184 #ifdef debug 00185 G4cout<<"G4QNCS::GetCrosSec: recCS="<<lastCS<<",lZ="<<lastN<<",lN="<<lastZ<<G4endl; 00186 #endif 00187 //} // M.K. Presence of H1 with high threshold breaks the syncronization 00188 #ifdef pdebug 00189 G4cout<<"G4QNCS::GetCS:1st, P="<<pMom<<"(MeV),CS="<<lastCS*millibarn<<"(mb)"<<G4endl; 00190 #endif 00191 return lastCS*millibarn; 00192 } // End of creation of the new set of parameters 00193 else 00194 { 00195 #ifdef debug 00196 G4cout<<"G4QNeutronNuclCrossSections::GetCS: Update lastI="<<lastI<<",j="<<j<<G4endl; 00197 #endif 00198 colP[lastI]=pMom; 00199 colCS[lastI]=lastCS; 00200 } 00201 } // End of parameters udate 00202 else if(pMom<=lastTH) 00203 { 00204 #ifdef debug 00205 G4cout<<"G4QNeutCS::GetCS:Current P="<<pMom<<" < Threshold="<<lastTH<<", CS=0"<<G4endl; 00206 #endif 00207 return 0.; // Momentum is below the Threshold Value -> CS=0 00208 } 00209 else if(std::fabs(lastP-pMom)<tolerance*pMom) 00210 //else if(lastP==pMom) // VI do not use tolerance 00211 { 00212 #ifdef debug 00213 G4cout<<"..G4QNCS::GetCS:OldNZ&P="<<lastP<<"="<<pMom<<",CS="<<lastCS*millibarn<<G4endl; 00214 #endif 00215 return lastCS*millibarn; // Use theLastCS 00216 } 00217 else // It is the last used -> use the current tables 00218 { 00219 #ifdef debug 00220 G4cout<<"-!-G4QNCS::GetCS:UseCur P="<<pMom<<",f="<<fCS<<",I="<<lastI<<",j="<<j<<G4endl; 00221 #endif 00222 lastCS=CalculateCrossSection(fCS,1,j,2112,lastZ,lastN,pMom); // Only read and UpdateDB 00223 lastP=pMom; 00224 } 00225 #ifdef debug 00226 G4cout<<"==>G4QNCS::GetCroSec: P="<<pMom<<"(MeV), CS="<<lastCS*millibarn<<"(mb)"<<G4endl; 00227 #endif 00228 return lastCS*millibarn; 00229 }
G4VQCrossSection * G4QNeutronNuclearCrossSection::GetPointer | ( | ) | [static] |
Definition at line 64 of file G4QNeutronNuclearCrossSection.cc.
Referenced by G4QHadronInelasticDataSet::GetIsoCrossSection(), G4QNGamma::GetMeanFreePath(), and G4QInelastic::GetMeanFreePath().
00065 { 00066 static G4QNeutronNuclearCrossSection theCrossSection; //**Static body of Cross Section** 00067 return &theCrossSection; 00068 }