G4QHyperonPlusElasticCrossSection Class Reference

#include <G4QHyperonPlusElasticCrossSection.hh>

Inheritance diagram for G4QHyperonPlusElasticCrossSection:

G4VQCrossSection

Public Member Functions

 ~G4QHyperonPlusElasticCrossSection ()
virtual G4double GetCrossSection (G4bool fCS, G4double pMom, G4int tgZ, G4int tgN, G4int pPDG=3222)
G4double CalculateCrossSection (G4bool CS, G4int F, G4int I, G4int pPDG, G4int Z, G4int N, G4double pP)
G4double GetSlope (G4int tZ, G4int tN, G4int pPDG)
G4double GetExchangeT (G4int tZ, G4int tN, G4int pPDG)
G4double GetHMaxT ()

Static Public Member Functions

static G4VQCrossSectionGetPointer ()

Protected Member Functions

 G4QHyperonPlusElasticCrossSection ()

Detailed Description

Definition at line 47 of file G4QHyperonPlusElasticCrossSection.hh.


Constructor & Destructor Documentation

G4QHyperonPlusElasticCrossSection::G4QHyperonPlusElasticCrossSection (  )  [protected]

Definition at line 100 of file G4QHyperonPlusElasticCrossSection.cc.

00101 {
00102 }

G4QHyperonPlusElasticCrossSection::~G4QHyperonPlusElasticCrossSection (  ) 

Definition at line 104 of file G4QHyperonPlusElasticCrossSection.cc.

00105 {
00106   std::vector<G4double*>::iterator pos;
00107   for (pos=CST.begin(); pos<CST.end(); pos++)
00108   { delete [] *pos; }
00109   CST.clear();
00110   for (pos=PAR.begin(); pos<PAR.end(); pos++)
00111   { delete [] *pos; }
00112   PAR.clear();
00113   for (pos=SST.begin(); pos<SST.end(); pos++)
00114   { delete [] *pos; }
00115   SST.clear();
00116   for (pos=S1T.begin(); pos<S1T.end(); pos++)
00117   { delete [] *pos; }
00118   S1T.clear();
00119   for (pos=B1T.begin(); pos<B1T.end(); pos++)
00120   { delete [] *pos; }
00121   B1T.clear();
00122   for (pos=S2T.begin(); pos<S2T.end(); pos++)
00123   { delete [] *pos; }
00124   S2T.clear();
00125   for (pos=B2T.begin(); pos<B2T.end(); pos++)
00126   { delete [] *pos; }
00127   B2T.clear();
00128   for (pos=S3T.begin(); pos<S3T.end(); pos++)
00129   { delete [] *pos; }
00130   S3T.clear();
00131   for (pos=B3T.begin(); pos<B3T.end(); pos++)
00132   { delete [] *pos; }
00133   B3T.clear();
00134   for (pos=S4T.begin(); pos<S4T.end(); pos++)
00135   { delete [] *pos; }
00136   S4T.clear();
00137   for (pos=B4T.begin(); pos<B4T.end(); pos++)
00138   { delete [] *pos; }
00139   B4T.clear();
00140 }


Member Function Documentation

G4double G4QHyperonPlusElasticCrossSection::CalculateCrossSection ( G4bool  CS,
G4int  F,
G4int  I,
G4int  pPDG,
G4int  Z,
G4int  N,
G4double  pP 
) [virtual]

Implements G4VQCrossSection.

Definition at line 286 of file G4QHyperonPlusElasticCrossSection.cc.

References G4cout, and G4endl.

Referenced by GetCrossSection().

00288 {
00289   // *** Begin of Associative Memory DB for acceleration of the cross section calculations
00290   static std::vector <G4double>  PIN;   // Vector of max initialized log(P) in the table
00291   // *** End of Static Definitions (Associative Memory Data Base) ***
00292   G4double pMom=pIU/GeV;                // All calculations are in GeV
00293   onlyCS=CS;                            // Flag to calculate only CS (not Si/Bi)
00294 #ifdef pdebug
00295   G4cout<<"G4QHyperonPlusElasticCS::CalcCS: onlyCS="<<onlyCS<<",F="<<F<<",p="<<pIU<<G4endl;
00296 #endif
00297   lastLP=std::log(pMom);                // Make a logarithm of the momentum for calculation
00298   if(F)                                 // This isotope was found in AMDB =>RETRIEVE/UPDATE
00299   {
00300     if(F<0)                             // the AMDB must be loded
00301     {
00302       lastPIN = PIN[I];                 // Max log(P) initialised for this table set
00303       lastPAR = PAR[I];                 // Pointer to the parameter set
00304       lastCST = CST[I];                 // Pointer to the total sross-section table
00305       lastSST = SST[I];                 // Pointer to the first squared slope
00306       lastS1T = S1T[I];                 // Pointer to the first mantissa
00307       lastB1T = B1T[I];                 // Pointer to the first slope
00308       lastS2T = S2T[I];                 // Pointer to the second mantissa
00309       lastB2T = B2T[I];                 // Pointer to the second slope
00310       lastS3T = S3T[I];                 // Pointer to the third mantissa
00311       lastB3T = B3T[I];                 // Pointer to the rhird slope
00312       lastS4T = S4T[I];                 // Pointer to the 4-th mantissa
00313       lastB4T = B4T[I];                 // Pointer to the 4-th slope
00314 #ifdef pdebug
00315       G4cout<<"G4QHyperonPlusElCS::CalCS:DB updated for I="<<I<<",*,PIN4="<<PIN[4]<<G4endl;
00316 #endif
00317     }
00318 #ifdef pdebug
00319     G4cout<<"G4QHyperonPlusElasticCS::CalcCS:*read*,LP="<<lastLP<<",PIN="<<lastPIN<<G4endl;
00320 #endif
00321     if(lastLP>lastPIN && lastLP<lPMax)
00322     {
00323       lastPIN=GetPTables(lastLP,lastPIN,PDG,tgZ,tgN);// Can update upper logP-Limit in tabs
00324 #ifdef pdebug
00325       G4cout<<"G4QHPElCS::CalcCS:updated(I),LP="<<lastLP<<"<IN["<<I<<"]="<<lastPIN<<G4endl;
00326 #endif
00327       PIN[I]=lastPIN;                   // Remember the new P-Limit of the tables
00328     }
00329   }
00330   else                                  // This isotope wasn't initialized => CREATE
00331   {
00332     lastPAR = new G4double[nPoints];    // Allocate memory for parameters of CS function
00333     lastPAR[nLast]=0;                   // Initialization for VALGRIND
00334     lastCST = new G4double[nPoints];    // Allocate memory for Tabulated CS function    
00335     lastSST = new G4double[nPoints];    // Allocate memory for Tabulated first sqaredSlope 
00336     lastS1T = new G4double[nPoints];    // Allocate memory for Tabulated first mantissa 
00337     lastB1T = new G4double[nPoints];    // Allocate memory for Tabulated first slope    
00338     lastS2T = new G4double[nPoints];    // Allocate memory for Tabulated second mantissa
00339     lastB2T = new G4double[nPoints];    // Allocate memory for Tabulated second slope   
00340     lastS3T = new G4double[nPoints];    // Allocate memory for Tabulated third mantissa 
00341     lastB3T = new G4double[nPoints];    // Allocate memory for Tabulated third slope    
00342     lastS4T = new G4double[nPoints];    // Allocate memory for Tabulated 4-th mantissa 
00343     lastB4T = new G4double[nPoints];    // Allocate memory for Tabulated 4-th slope    
00344 #ifdef pdebug
00345     G4cout<<"G4QHyperonPlusElasticCS::CalcCS:*ini*,lstLP="<<lastLP<<",min="<<lPMin<<G4endl;
00346 #endif
00347     lastPIN = GetPTables(lastLP,lPMin,PDG,tgZ,tgN); // Returns the new P-limit for tables
00348 #ifdef pdebug
00349     G4cout<<"G4QHyPlElCS::CCS:i,Z="<<tgZ<<",N="<<tgN<<",PDG="<<PDG<<",LP"<<lastPIN<<G4endl;
00350 #endif
00351     PIN.push_back(lastPIN);             // Fill parameters of CS function to AMDB
00352     PAR.push_back(lastPAR);             // Fill parameters of CS function to AMDB
00353     CST.push_back(lastCST);             // Fill Tabulated CS function to AMDB    
00354     SST.push_back(lastSST);             // Fill Tabulated first sq.slope to AMDB 
00355     S1T.push_back(lastS1T);             // Fill Tabulated first mantissa to AMDB 
00356     B1T.push_back(lastB1T);             // Fill Tabulated first slope to AMDB    
00357     S2T.push_back(lastS2T);             // Fill Tabulated second mantissa to AMDB 
00358     B2T.push_back(lastB2T);             // Fill Tabulated second slope to AMDB    
00359     S3T.push_back(lastS3T);             // Fill Tabulated third mantissa to AMDB 
00360     B3T.push_back(lastB3T);             // Fill Tabulated third slope to AMDB    
00361     S4T.push_back(lastS4T);             // Fill Tabulated 4-th mantissa to AMDB 
00362     B4T.push_back(lastB4T);             // Fill Tabulated 4-th slope to AMDB    
00363   } // End of creation/update of the new set of parameters and tables
00364   // =-------------= NOW Update (if necessary) and Calculate the Cross Section =----------=
00365 #ifdef pdebug
00366   G4cout<<"G4QHPElCS::CalcCS:?update?,LP="<<lastLP<<",IN="<<lastPIN<<",ML="<<lPMax<<G4endl;
00367 #endif
00368   if(lastLP>lastPIN && lastLP<lPMax)
00369   {
00370     lastPIN = GetPTables(lastLP,lastPIN,PDG,tgZ,tgN);
00371 #ifdef pdebug
00372     G4cout<<"G4QHyperonPlusElCS::CalcCS:*updated(O)*,LP="<<lastLP<<"<IN="<<lastPIN<<G4endl;
00373 #endif
00374   }
00375 #ifdef pdebug
00376   G4cout<<"G4QHyPlElCS::CalcCS:lastLP="<<lastLP<<",lPM="<<lPMin<<",lPIN="<<lastPIN<<G4endl;
00377 #endif
00378   if(!onlyCS) lastTM=GetQ2max(PDG, tgZ, tgN, pMom); // Calculate (-t)_max=Q2_max (GeV2)
00379 #ifdef pdebug
00380   G4cout<<"G4QHyPlElCrosSec::CalcCS: oCS="<<onlyCS<<",-t="<<lastTM<<", p="<<lastLP<<G4endl;
00381 #endif
00382   if(lastLP>lPMin && lastLP<=lastPIN)   // Linear fit is made using precalculated tables
00383   {
00384     if(lastLP==lastPIN)
00385     {
00386       G4double shift=(lastLP-lPMin)/dlnP+.000001; // Log distance from lPMin
00387       G4int    blast=static_cast<int>(shift); // this is a bin number of the lower edge (0)
00388       if(blast<0||blast>=nLast) G4cout<<"G4QHyPlElCS::CCS:b="<<blast<<",n="<<nLast<<G4endl;
00389       lastSIG = lastCST[blast];
00390       if(!onlyCS)                       // Skip the differential cross-section parameters
00391       {
00392         theSS  = lastSST[blast];
00393         theS1  = lastS1T[blast];
00394         theB1  = lastB1T[blast];
00395         theS2  = lastS2T[blast];
00396         theB2  = lastB2T[blast];
00397         theS3  = lastS3T[blast];
00398         theB3  = lastB3T[blast];
00399         theS4  = lastS4T[blast];
00400         theB4  = lastB4T[blast];
00401       }
00402 #ifdef pdebug
00403       G4cout<<"G4QHyperonPlusElasticCS::CalculateCS:(E) S1="<<theS1<<",B1="<<theB1<<G4endl;
00404 #endif
00405     }
00406     else
00407     {
00408       G4double shift=(lastLP-lPMin)/dlnP;        // a shift from the beginning of the table
00409       G4int    blast=static_cast<int>(shift);    // the lower bin number
00410       if(blast<0)   blast=0;
00411       if(blast>=nLast) blast=nLast-1;            // low edge of the last bin
00412       shift-=blast;                              // step inside the unit bin
00413       G4int lastL=blast+1;                       // the upper bin number
00414       G4double SIGL=lastCST[blast];              // the basic value of the cross-section
00415       lastSIG= SIGL+shift*(lastCST[lastL]-SIGL); // calculated total elastic cross-section
00416 #ifdef pdebug
00417       G4cout<<"G4QHyperonPlusElasticCrossSection::CalcCrossSection:Sig="<<lastSIG<<",P="
00418             <<pMom<<",Z="<<tgZ<<",N="<<tgN<<",PDG="<<PDG<<",onlyCS="<<onlyCS<<G4endl;
00419 #endif
00420       if(!onlyCS)                       // Skip the differential cross-section parameters
00421       {
00422         G4double SSTL=lastSST[blast];           // the low bin of the first squared slope
00423         theSS=SSTL+shift*(lastSST[lastL]-SSTL); // the basic value of the first sq.slope
00424         G4double S1TL=lastS1T[blast];           // the low bin of the first mantissa
00425         theS1=S1TL+shift*(lastS1T[lastL]-S1TL); // the basic value of the first mantissa
00426         G4double B1TL=lastB1T[blast];           // the low bin of the first slope
00427 #ifdef pdebug
00428         G4cout<<"G4QHyperonPlusElasticCS::CalcCrossSect: b="<<blast<<", ls="<<lastL<<",SL="
00429               <<S1TL<<",SU="<<lastS1T[lastL]<<",BL="<<B1TL<<",BU="<<lastB1T[lastL]<<G4endl;
00430 #endif
00431         theB1=B1TL+shift*(lastB1T[lastL]-B1TL); // the basic value of the first slope
00432         G4double S2TL=lastS2T[blast];           // the low bin of the second mantissa
00433         theS2=S2TL+shift*(lastS2T[lastL]-S2TL); // the basic value of the second mantissa
00434         G4double B2TL=lastB2T[blast];           // the low bin of the second slope
00435         theB2=B2TL+shift*(lastB2T[lastL]-B2TL); // the basic value of the second slope
00436         G4double S3TL=lastS3T[blast];           // the low bin of the third mantissa
00437         theS3=S3TL+shift*(lastS3T[lastL]-S3TL); // the basic value of the third mantissa
00438 #ifdef pdebug
00439         G4cout<<"G4QHyperonPlusElasticCrossSection::CalcCS: s3l="<<S3TL<<", sh3="<<shift
00440               <<", s3h="<<lastS3T[lastL]<<", b="<<blast<<", l="<<lastL<<G4endl;
00441 #endif
00442         G4double B3TL=lastB3T[blast];           // the low bin of the third slope
00443         theB3=B3TL+shift*(lastB3T[lastL]-B3TL); // the basic value of the third slope
00444         G4double S4TL=lastS4T[blast];           // the low bin of the 4-th mantissa
00445         theS4=S4TL+shift*(lastS4T[lastL]-S4TL); // the basic value of the 4-th mantissa
00446 #ifdef pdebug
00447         G4cout<<"G4QHyperonPlusElasticCrossSection::CalcCS: s4l="<<S4TL<<", sh4="<<shift
00448               <<",s4h="<<lastS4T[lastL]<<",b="<<blast<<",l="<<lastL<<G4endl;
00449 #endif
00450         G4double B4TL=lastB4T[blast];           // the low bin of the 4-th slope
00451         theB4=B4TL+shift*(lastB4T[lastL]-B4TL); // the basic value of the 4-th slope
00452       }
00453 #ifdef pdebug
00454       G4cout<<"G4QHyperonPlusElasticCS::CalculateCS:(I) S1="<<theS1<<",B1="<<theB1<<G4endl;
00455 #endif
00456     }
00457   }
00458   else lastSIG=GetTabValues(lastLP, PDG, tgZ, tgN); // Direct calculation beyond the table
00459   if(lastSIG<0.) lastSIG = 0.;                   // @@ a Warning print can be added
00460 #ifdef pdebug
00461   G4cout<<"G4QHyperonPlusElasticCrossSection::CalculateCS: END, onlyCS="<<onlyCS<<G4endl;
00462 #endif
00463   return lastSIG;
00464 }

G4double G4QHyperonPlusElasticCrossSection::GetCrossSection ( G4bool  fCS,
G4double  pMom,
G4int  tgZ,
G4int  tgN,
G4int  pPDG = 3222 
) [virtual]

Reimplemented from G4VQCrossSection.

Definition at line 151 of file G4QHyperonPlusElasticCrossSection.cc.

References CalculateCrossSection(), G4cout, G4endl, and G4VQCrossSection::ThresholdEnergy().

00153 {
00154   static std::vector <G4int>    colN;  // Vector of N for calculated nuclei (isotops)
00155   static std::vector <G4int>    colZ;  // Vector of Z for calculated nuclei (isotops)
00156   static std::vector <G4double> colP;  // Vector of last momenta for the reaction
00157   static std::vector <G4double> colTH; // Vector of energy thresholds for the reaction
00158   static std::vector <G4double> colCS; // Vector of last cross sections for the reaction
00159   // ***---*** End of the mandatory Static Definitions of the Associative Memory ***---***
00160   G4double pEn=pMom;
00161   onlyCS=fCS;
00162 #ifdef pdebug
00163   G4cout<<"G4QHypPlElCS::GetCS:>>> f="<<fCS<<", p="<<pMom<<", Z="<<tgZ<<"("<<lastZ<<") ,N="
00164         <<tgN<<"("<<lastN<<"), T="<<pEn<<"("<<lastTH<<")"<<",Sz="<<colN.size()<<G4endl;
00165   //CalculateCrossSection(fCS,-27,j,pPDG,lastZ,lastN,pMom); // DUMMY TEST
00166 #endif
00167   if(pPDG!=3222)
00168   {
00169     G4cout<<"*Warning*G4QHyperonPlusElaCS::GetCS:*>Found pPDG="<<pPDG<<" ==> CS=0"<<G4endl;
00170     //CalculateCrossSection(fCS,-27,j,pPDG,lastZ,lastN,pMom); // DUMMY TEST
00171     return 0.;                         // projectile PDG=0 is a mistake (?!) @@
00172   }
00173   G4bool in=false;                   // By default the isotope must be found in the AMDB
00174   lastP   = 0.;                      // New momentum history (nothing to compare with)
00175   lastN   = tgN;                     // The last N of the calculated nucleus
00176   lastZ   = tgZ;                     // The last Z of the calculated nucleus
00177   lastI   = colN.size();             // Size of the Associative Memory DB in the heap
00178   if(lastI) for(G4int i=0; i<lastI; i++) // Loop over proj/tgZ/tgN lines of DB
00179   {                                  // The nucleus with projPDG is found in AMDB
00180     if(colN[i]==tgN && colZ[i]==tgZ) // Isotope is foind in AMDB
00181     {
00182       lastI=i;
00183       lastTH =colTH[i];              // Last THreshold (A-dependent)
00184 #ifdef pdebug
00185       G4cout<<"G4QHPElCS::GetCS:*Found* P="<<pMom<<",Threshold="<<lastTH<<",i="<<i<<G4endl;
00186       //CalculateCrossSection(fCS,-27,i,pPDG,lastZ,lastN,pMom); // DUMMY TEST
00187 #endif
00188       if(pEn<=lastTH)
00189       {
00190 #ifdef pdebug
00191         G4cout<<"G4QHPElCS::GetCS:Found T="<<pEn<<" < Threshold="<<lastTH<<",CS=0"<<G4endl;
00192         //CalculateCrossSection(fCS,-27,i,pPDG,lastZ,lastN,pMom); // DUMMY TEST
00193 #endif
00194         return 0.;                   // Energy is below the Threshold value
00195       }
00196       lastP  =colP [i];              // Last Momentum  (A-dependent)
00197       lastCS =colCS[i];              // Last CrossSect (A-dependent)
00198       //  if(std::fabs(lastP/pMom-1.)<tolerance) //VI (do not use tolerance)
00199       if(lastP == pMom)              // Do not recalculate
00200       {
00201 #ifdef pdebug
00202         G4cout<<"G4QHyperonPlusElastCS::GetCS: P="<<pMom<<",CS="<<lastCS*millibarn<<G4endl;
00203 #endif
00204         CalculateCrossSection(fCS,-1,i,pPDG,lastZ,lastN,pMom); // Update param's only
00205         return lastCS*millibarn;     // Use theLastCS
00206       }
00207       in = true;                       // This is the case when the isotop is found in DB
00208       // Momentum pMom is in IU ! @@ Units
00209 #ifdef pdebug
00210       G4cout<<"G4QHPElCS::G:UpdateDB P="<<pMom<<",f="<<fCS<<",I="<<lastI<<",i="<<i<<G4endl;
00211 #endif
00212       lastCS=CalculateCrossSection(fCS,-1,i,pPDG,lastZ,lastN,pMom); // read & update
00213 #ifdef pdebug
00214       G4cout<<"G4QHyperonPlusElCS::GetCrosSec:**>New(inDB) Calculated CS="<<lastCS<<G4endl;
00215       //CalculateCrossSection(fCS,-27,i,pPDG,lastZ,lastN,pMom); // DUMMY TEST
00216 #endif
00217       if(lastCS<=0. && pEn>lastTH)    // Correct the threshold
00218       {
00219 #ifdef pdebug
00220         G4cout<<"G4QHyperonPlusElCS::GetCS:T="<<pEn<<"(CS=0) > Threshold="<<lastTH<<G4endl;
00221 #endif
00222         lastTH=pEn;
00223       }
00224       break;                           // Go out of the LOOP with found lastI
00225     }
00226 #ifdef pdebug
00227     G4cout<<"---G4QHyperonPlusElasticCrossSection::GetCrosSec:pPDG="<<pPDG<<",i="<<i<<",N="
00228           <<colN[i]<<",Z["<<i<<"]="<<colZ[i]<<G4endl;
00229     //CalculateCrossSection(fCS,-27,i,pPDG,lastZ,lastN,pMom); // DUMMY TEST
00230 #endif
00231   } // End of attampt to find the nucleus in DB
00232   if(!in)                            // This nucleus has not been calculated previously
00233   {
00234 #ifdef pdebug
00235     G4cout<<"G4QHyperonPlusElCS::GetCS:CNewP="<<pMom<<",f="<<fCS<<",lastI="<<lastI<<G4endl;
00236 #endif
00238     lastCS=CalculateCrossSection(fCS,0,lastI,pPDG,lastZ,lastN,pMom);//calculate&create
00239     if(lastCS<=0.)
00240     {
00241       lastTH = ThresholdEnergy(tgZ, tgN); // The Threshold Energy which is now the last
00242 #ifdef pdebug
00243       G4cout<<"G4QHyperonPlusElasticCrossSection::GCS:NewThr="<<lastTH<<",T="<<pEn<<G4endl;
00244 #endif
00245       if(pEn>lastTH)
00246       {
00247 #ifdef pdebug
00248         G4cout<<"G4QHyperonPlusElCS::GCS:1st T="<<pEn<<"(CS=0)>Threshold="<<lastTH<<G4endl;
00249 #endif
00250         lastTH=pEn;
00251       }
00252     }
00253 #ifdef pdebug
00254     G4cout<<"G4QHyperonPlusElCS::GetCS:NCS="<<lastCS<<",lZ="<<lastN<<",lN="<<lastZ<<G4endl;
00255     //CalculateCrossSection(fCS,-27,lastI,pPDG,lastZ,lastN,pMom); // DUMMY TEST
00256 #endif
00257     colN.push_back(tgN);
00258     colZ.push_back(tgZ);
00259     colP.push_back(pMom);
00260     colTH.push_back(lastTH);
00261     colCS.push_back(lastCS);
00262 #ifdef pdebug
00263     G4cout<<"G4QHPElCS::GetCS:1st,P="<<pMom<<"(MeV),CS="<<lastCS*millibarn<<"(mb)"<<G4endl;
00264     //CalculateCrossSection(fCS,-27,lastI,pPDG,lastZ,lastN,pMom); // DUMMY TEST
00265 #endif
00266     return lastCS*millibarn;
00267   } // End of creation of the new set of parameters
00268   else
00269   {
00270 #ifdef pdebug
00271     G4cout<<"G4QHyperonPlusElasticCrossSection::GetCS: Update lastI="<<lastI<<G4endl;
00272 #endif
00273     colP[lastI]=pMom;
00274     colCS[lastI]=lastCS;
00275   }
00276 #ifdef pdebug
00277   G4cout<<"G4QHPElCS::GetCSec:End,P="<<pMom<<"(MeV),CS="<<lastCS*millibarn<<"(mb)"<<G4endl;
00278   //CalculateCrossSection(fCS,-27,lastI,pPDG,lastZ,lastN,pMom); // DUMMY TEST
00279   G4cout<<"G4QHyperonPlusElasticCrossSection::GetCrosSec:**End**, onlyCS="<<onlyCS<<G4endl;
00280 #endif
00281   return lastCS*millibarn;
00282 }

G4double G4QHyperonPlusElasticCrossSection::GetExchangeT ( G4int  tZ,
G4int  tN,
G4int  pPDG 
) [virtual]

Reimplemented from G4VQCrossSection.

Definition at line 734 of file G4QHyperonPlusElasticCrossSection.cc.

References G4cout, G4endl, and G4UniformRand.

00735 {
00736   static const G4double GeVSQ=gigaelectronvolt*gigaelectronvolt;
00737   static const G4double third=1./3.;
00738   static const G4double fifth=1./5.;
00739   static const G4double sevth=1./7.;
00740 #ifdef tdebug
00741   G4cout<<"G4QHyPlElCS::GetExchT:F="<<onlyCS<<",Z="<<tgZ<<",N="<<tgN<<",PDG="<<PDG<<G4endl;
00742 #endif
00743   if(PDG!=3222)G4cout<<"*Warning*G4QHyperonPlusElasticCrossSection::GET:PDG="<<PDG<<G4endl;
00744   if(onlyCS)G4cout<<"*Warning*G4QHyperonPlusElasticCrossSection::GetExT: onlyCS=1"<<G4endl;
00745   if(lastLP<-4.3) return lastTM*GeVSQ*G4UniformRand();// S-wave for p<14 MeV/c (kinE<.1MeV)
00746   G4double q2=0.;
00747   if(tgZ==1 && tgN==0)                // ===> p+p=p+p
00748   {
00749 #ifdef tdebug
00750     G4cout<<"G4QHyperonPlusElCS::GetExT: TM="<<lastTM<<",S1="<<theS1<<",B1="<<theB1<<",S2="
00751           <<theS2<<",B2="<<theB2<<",S3="<<theS3<<",B3="<<theB3<<",GeV2="<<GeVSQ<<G4endl;
00752 #endif
00753     G4double E1=lastTM*theB1;
00754     G4double R1=(1.-std::exp(-E1));
00755 #ifdef tdebug
00756     G4double ts1=-std::log(1.-R1)/theB1;
00757     G4double ds1=std::fabs(ts1-lastTM)/lastTM;
00758     if(ds1>.0001)
00759       G4cout<<"*Warning*G4QHyperonPlusElasticCrossSection::GetExT:1p "<<ts1<<"#"<<lastTM
00760             <<",d="<<ds1<<",R1="<<R1<<",E1="<<E1<<G4endl;
00761 #endif
00762     G4double E2=lastTM*theB2;
00763     G4double R2=(1.-std::exp(-E2*E2*E2));
00764 #ifdef tdebug
00765     G4double ts2=std::pow(-std::log(1.-R2),.333333333)/theB2;
00766     G4double ds2=std::fabs(ts2-lastTM)/lastTM;
00767     if(ds2>.0001)
00768       G4cout<<"*Warning*G4QHyperonPlusElasticCrossSection::GetExT:2p "<<ts2<<"#"<<lastTM
00769             <<",d="<<ds2<<",R2="<<R2<<",E2="<<E2<<G4endl;
00770 #endif
00771     G4double E3=lastTM*theB3;
00772     G4double R3=(1.-std::exp(-E3));
00773 #ifdef tdebug
00774     G4double ts3=-std::log(1.-R3)/theB3;
00775     G4double ds3=std::fabs(ts3-lastTM)/lastTM;
00776     if(ds3>.0001)
00777       G4cout<<"*Warning*G4QHyperonPlusElasticCrossSection::GetExT:3p "<<ts3<<"#"<<lastTM
00778             <<",d="<<ds3<<",R3="<<R1<<",E3="<<E3<<G4endl;
00779 #endif
00780     G4double I1=R1*theS1/theB1;
00781     G4double I2=R2*theS2;
00782     G4double I3=R3*theS3;
00783     G4double I12=I1+I2;
00784     G4double rand=(I12+I3)*G4UniformRand();
00785     if     (rand<I1 )
00786     {
00787       G4double ran=R1*G4UniformRand();
00788       if(ran>1.) ran=1.;
00789       q2=-std::log(1.-ran)/theB1;
00790     }
00791     else if(rand<I12)
00792     {
00793       G4double ran=R2*G4UniformRand();
00794       if(ran>1.) ran=1.;
00795       q2=-std::log(1.-ran);
00796       if(q2<0.) q2=0.;
00797       q2=std::pow(q2,third)/theB2;
00798     }
00799     else
00800     {
00801       G4double ran=R3*G4UniformRand();
00802       if(ran>1.) ran=1.;
00803       q2=-std::log(1.-ran)/theB3;
00804     }
00805   }
00806   else
00807   {
00808     G4double a=tgZ+tgN;
00809 #ifdef tdebug
00810     G4cout<<"G4QHyperonPlusElasticCrossSection::GetExT:a="<<a<<",t="<<lastTM<<",S1="<<theS1
00811           <<",B1="<<theB1<<",SS="<<theSS<<",S2="<<theS2<<",B2="<<theB2<<",S3="<<theS3
00812           <<",B3="<<theB3<<",S4="<<theS4<<",B4="<<theB4<<G4endl;
00813 #endif
00814     G4double E1=lastTM*(theB1+lastTM*theSS);
00815     G4double R1=(1.-std::exp(-E1));
00816     G4double tss=theSS+theSS; // for future solution of quadratic equation (imediate check)
00817 #ifdef tdebug
00818     G4double ts1=-std::log(1.-R1)/theB1;
00819     if(std::fabs(tss)>1.e-7) ts1=(std::sqrt(theB1*(theB1+(tss+tss)*ts1))-theB1)/tss;
00820     G4double ds1=(ts1-lastTM)/lastTM;
00821     if(ds1>.0001)
00822       G4cout<<"*Warning*G4QHyperonPlusElasticCrossSection::GetExT:1a "<<ts1<<"#"<<lastTM
00823             <<",d="<<ds1<<",R1="<<R1<<",E1="<<E1<<G4endl;
00824 #endif
00825     G4double tm2=lastTM*lastTM;
00826     G4double E2=lastTM*tm2*theB2;                   // power 3 for lowA, 5 for HighA (1st)
00827     if(a>6.5)E2*=tm2;                               // for heavy nuclei
00828     G4double R2=(1.-std::exp(-E2));
00829 #ifdef tdebug
00830     G4double ts2=-std::log(1.-R2)/theB2;
00831     if(a<6.5)ts2=std::pow(ts2,third);
00832     else     ts2=std::pow(ts2,fifth);
00833     G4double ds2=std::fabs(ts2-lastTM)/lastTM;
00834     if(ds2>.0001)
00835       G4cout<<"*Warning*G4QHyperonPlusElasticCrossSection::GetExT:2a "<<ts2<<"#"<<lastTM
00836             <<",d="<<ds2<<",R2="<<R2<<",E2="<<E2<<G4endl;
00837 #endif
00838     G4double E3=lastTM*theB3;
00839     if(a>6.5)E3*=tm2*tm2*tm2;                       // power 1 for lowA, 7 (2nd) for HighA
00840     G4double R3=(1.-std::exp(-E3));
00841 #ifdef tdebug
00842     G4double ts3=-std::log(1.-R3)/theB3;
00843     if(a>6.5)ts3=std::pow(ts3,sevth);
00844     G4double ds3=std::fabs(ts3-lastTM)/lastTM;
00845     if(ds3>.0001)
00846       G4cout<<"*Warning*G4QHyperonPlusElasticCrossSection::GetExT:3a "<<ts3<<"#"<<lastTM
00847             <<",d="<<ds3<<",R3="<<R3<<",E3="<<E3<<G4endl;
00848 #endif
00849     G4double E4=lastTM*theB4;
00850     G4double R4=(1.-std::exp(-E4));
00851 #ifdef tdebug
00852     G4double ts4=-std::log(1.-R4)/theB4;
00853     G4double ds4=std::fabs(ts4-lastTM)/lastTM;
00854     if(ds4>.0001)
00855       G4cout<<"*Warning*G4QHyperonPlusElasticCrossSection::GetExT:4a "<<ts4<<"#"<<lastTM
00856             <<",d="<<ds4<<",R4="<<R4<<",E4="<<E4<<G4endl;
00857 #endif
00858     G4double I1=R1*theS1;
00859     G4double I2=R2*theS2;
00860     G4double I3=R3*theS3;
00861     G4double I4=R4*theS4;
00862     G4double I12=I1+I2;
00863     G4double I13=I12+I3;
00864     G4double rand=(I13+I4)*G4UniformRand();
00865 #ifdef tdebug
00866     G4cout<<"G4QHPElCS::GExT:1="<<I1<<",2="<<I2<<",3="<<I3<<",4="<<I4<<",r="<<rand<<G4endl;
00867 #endif
00868     if(rand<I1)
00869     {
00870       G4double ran=R1*G4UniformRand();
00871       if(ran>1.) ran=1.;
00872       q2=-std::log(1.-ran)/theB1;
00873       if(std::fabs(tss)>1.e-7) q2=(std::sqrt(theB1*(theB1+(tss+tss)*q2))-theB1)/tss;
00874 #ifdef tdebug
00875       G4cout<<"G4QHPElCS::GExT:Q2="<<q2<<",ss="<<tss/2<<",b1="<<theB1<<",t1="<<ts1<<G4endl;
00876 #endif
00877     }
00878     else if(rand<I12)
00879     {
00880       G4double ran=R2*G4UniformRand();
00881       if(ran>1.) ran=1.;
00882       q2=-std::log(1.-ran)/theB2;
00883       if(q2<0.) q2=0.;
00884       if(a<6.5) q2=std::pow(q2,third);
00885       else      q2=std::pow(q2,fifth);
00886 #ifdef tdebug
00887       G4cout<<"G4QHPElCS::GetExT: Q2="<<q2<<",r2="<<R2<<",b2="<<theB2<<",t2="<<ts2<<G4endl;
00888 #endif
00889     }
00890     else if(rand<I13)
00891     {
00892       G4double ran=R3*G4UniformRand();
00893       if(ran>1.) ran=1.;
00894       q2=-std::log(1.-ran)/theB3;
00895       if(q2<0.) q2=0.;
00896       if(a>6.5) q2=std::pow(q2,sevth);
00897 #ifdef tdebug
00898       G4cout<<"G4QHPElCS::GetExT: Q2="<<q2<<",r3="<<R2<<",b3="<<theB2<<",t3="<<ts2<<G4endl;
00899 #endif
00900     }
00901     else
00902     {
00903       G4double ran=R4*G4UniformRand();
00904       if(ran>1.) ran=1.;
00905       q2=-std::log(1.-ran)/theB4;
00906       if(a<6.5) q2=lastTM-q2;                    // u reduced for lightA (starts from 0)
00907 #ifdef tdebug
00908       G4cout<<"G4QHPElCS::GExT:Q2="<<q2<<",m="<<lastTM<<",b4="<<theB3<<",t4="<<ts3<<G4endl;
00909 #endif
00910     }
00911   }
00912   if(q2<0.) q2=0.;
00913   if(!(q2>=-1.||q2<=1.)) G4cout<<"*NAN*G4QHyperonPlusElasticCS::GetExchT: -t="<<q2<<G4endl;
00914   if(q2>lastTM)
00915   {
00916 #ifdef tdebug
00917     G4cout<<"*Warning*G4QHyperonPlusElasticCrossSection::GET:-t="<<q2<<">"<<lastTM<<G4endl;
00918 #endif
00919     q2=lastTM;
00920   }
00921   return q2*GeVSQ;
00922 }

G4double G4QHyperonPlusElasticCrossSection::GetHMaxT (  )  [virtual]

Reimplemented from G4VQCrossSection.

Definition at line 950 of file G4QHyperonPlusElasticCrossSection.cc.

00951 {
00952   static const G4double HGeVSQ=gigaelectronvolt*gigaelectronvolt/2.;
00953   return lastTM*HGeVSQ;
00954 }

G4VQCrossSection * G4QHyperonPlusElasticCrossSection::GetPointer (  )  [static]

Definition at line 143 of file G4QHyperonPlusElasticCrossSection.cc.

Referenced by G4QElastic::GetMeanFreePath(), and G4QElastic::PostStepDoIt().

00144 {
00145   static G4QHyperonPlusElasticCrossSection theCrossSection;//StaticBody of the QElCrossSect
00146   return &theCrossSection;
00147 }

G4double G4QHyperonPlusElasticCrossSection::GetSlope ( G4int  tZ,
G4int  tN,
G4int  pPDG 
) [virtual]

Reimplemented from G4VQCrossSection.

Definition at line 925 of file G4QHyperonPlusElasticCrossSection.cc.

References FatalException, G4cout, G4endl, and G4Exception().

00926 {
00927   static const G4double GeVSQ=gigaelectronvolt*gigaelectronvolt;
00928 #ifdef tdebug
00929   G4cout<<"G4QHyperonPlusECS::GetS:"<<onlyCS<<",Z="<<tgZ<<",N="<<tgN<<",PDG="<<PDG<<G4endl;
00930 #endif
00931   if(onlyCS)G4cout<<"*Warning*G4QHyperonPlusElasticCrossSection::GetSl:onlCS=true"<<G4endl;
00932   if(lastLP<-4.3) return 0.;          // S-wave for p<14 MeV/c (kinE<.1MeV)
00933   if(PDG != 3222)
00934   {
00935     // G4cout<<"*Error*G4QHyperonPlusElasticCrossSection::GetSlope: PDG="<<PDG<<", Z="<<tgZ
00936     //       <<", N="<<tgN<<", while it is defined only for PDG=3222"<<G4endl;
00937     // throw G4QException("G4QHyperonPlusElasticCrossSection::GetSlope:SigP are implemented");
00938     G4ExceptionDescription ed;
00939     ed << "PDG = " << PDG << ", Z = " << tgZ << ", N = " << tgN
00940        << ", while it is defined only for PDG=3222" << G4endl;
00941     G4Exception("G4QHyperonPlusElasticCrossSection::GetSlope()", "HAD_CHPS_0000",
00942                 FatalException, ed);
00943   }
00944   if(theB1<0.) theB1=0.;
00945   if(!(theB1>=-1.||theB1<=1.)) G4cout<<"*NAN*G4QHyperonPlusElCS::GetSlope:"<<theB1<<G4endl;
00946   return theB1/GeVSQ;
00947 }


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:53:08 2013 for Geant4 by  doxygen 1.4.7