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

#include <G4ChipsNeutronElasticXS.hh>

Inheritance diagram for G4ChipsNeutronElasticXS:
Collaboration diagram for G4ChipsNeutronElasticXS:

Public Member Functions

 G4ChipsNeutronElasticXS ()
 
 ~G4ChipsNeutronElasticXS ()
 
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)
 
G4double GetHMaxT ()
 
- 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 G4ChipsNeutronElasticXS.hh.

Constructor & Destructor Documentation

G4ChipsNeutronElasticXS::G4ChipsNeutronElasticXS ( )

Definition at line 62 of file G4ChipsNeutronElasticXS.cc.

62  :G4VCrossSectionDataSet(Default_Name()), nPoints(128), nLast(nPoints-1)
63 {
64  lPMin=-8.; // Min tabulated log Momentum (D)
65  lPMax= 8.; // Max tabulated log Momentum (D)
66  dlnP=(lPMax-lPMin)/nLast;// Log step in table (D)
67  onlyCS=true;// Flag to calc only CS (not Si/Bi)(L)
68  lastSIG=0.; // Last calculated cross section (L)
69  lastLP=-10.;// Last log(momOfIncidentHadron) (L)
70  lastTM=0.; // Last t_maximum (L)
71  theSS=0.; // The Last sq.slope of 1st difMax (L)
72  theS1=0.; // The Last mantissa of 1st difMax (L)
73  theB1=0.; // The Last slope of 1st difr. Max (L)
74  theS2=0.; // The Last mantissa of 2nd difMax (L)
75  theB2=0.; // The Last slope of 2nd difr. Max (L)
76  theS3=0.; // The Last mantissa of 3d difrMax (L)
77  theB3=0.; // The Last slope of 3d difructMax (L)
78  theS4=0.; // The Last mantissa of 4th difMax (L)
79  theB4=0.; // The Last slope of 4th difr. Max (L)
80  lastTZ=0; // Last atomic number of the target
81  lastTN=0; // Last # of neutrons in the target
82  lastPIN=0.; // Last initialized max momentum
83  lastCST=0; // Elastic cross-section table
84  lastPAR=0; // Parameters of FunctionalCalculation
85  lastSST=0; // E-dep of sq.slope of the 1st difMax
86  lastS1T=0; // E-dep of mantissa of the 1st difMax
87  lastB1T=0; // E-dep of theSlope of the 1st difMax
88  lastS2T=0; // E-dep of mantissa of the 2nd difMax
89  lastB2T=0; // E-dep of theSlope of the 2nd difMax
90  lastS3T=0; // E-dep of mantissa of the 3d difrMax
91  lastB3T=0; // E-dep of the slope of the 3d difMax
92  lastS4T=0; // E-dep of mantissa of the 4th difMax
93  lastB4T=0; // E-dep of theSlope of the 4th difMax
94  lastN=0; // The last N of calculated nucleus
95  lastZ=0; // The last Z of calculated nucleus
96  lastP=0.; // Last used in cross section Momentum
97  lastTH=0.; // Last threshold momentum
98  lastCS=0.; // Last value of the Cross Section
99  lastI=0; // The last position in the DAMDB
100 
101  mNeut= G4Neutron::Neutron()->GetPDGMass()*.001;// MeV to GeV
102  mProt= G4Proton::Proton()->GetPDGMass()*.001;// MeV to GeV
103  mNeut2= mNeut*mNeut;
104 }
G4VCrossSectionDataSet(const G4String &nam="")
static const char * Default_Name()
static G4Proton * Proton()
Definition: G4Proton.cc:93
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
G4double GetPDGMass() const

Here is the call graph for this function:

G4ChipsNeutronElasticXS::~G4ChipsNeutronElasticXS ( )

Definition at line 106 of file G4ChipsNeutronElasticXS.cc.

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

Member Function Documentation

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

Reimplemented from G4VCrossSectionDataSet.

Definition at line 145 of file G4ChipsNeutronElasticXS.cc.

146 {
147  outFile << "G4ChipsNeutronElasticXS provides the elastic cross\n"
148  << "section for neutron nucleus scattering as a function of incident\n"
149  << "momentum. The cross section is calculated using M. Kossov's\n"
150  << "CHIPS parameterization of cross section data.\n";
151 }
static const char* G4ChipsNeutronElasticXS::Default_Name ( )
inlinestatic

Definition at line 55 of file G4ChipsNeutronElasticXS.hh.

55 {return "ChipsNeutronElasticXS";}

Here is the caller graph for this function:

G4double G4ChipsNeutronElasticXS::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 173 of file G4ChipsNeutronElasticXS.cc.

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

Definition at line 1864 of file G4ChipsNeutronElasticXS.cc.

1865 {
1866  static const G4double GeVSQ=gigaelectronvolt*gigaelectronvolt;
1867  static const G4double third=1./3.;
1868  static const G4double fifth=1./5.;
1869  static const G4double sevth=1./7.;
1870  if(PDG!=2112) G4cout<<"*Warning*G4ChipsNeutronElasticXS::GetExT:PDG="<<PDG<<G4endl;
1871  if(onlyCS) G4cout<<"*Warning*G4ChipsNeutronElasticXS::GetExchangeT:onCS=1"<<G4endl;
1872  if(lastLP<-4.3) return lastTM*GeVSQ*G4UniformRand();// S-wave for p<14 MeV/c (kinE<.1MeV)
1873  G4double q2=0.;
1874  if(tgZ==1 && tgN==0) // ===> n+p=n+p
1875  {
1876  G4double E1=lastTM*theB1;
1877  G4double R1=(1.-std::exp(-E1));
1878  G4double E2=lastTM*theB2;
1879  G4double R2=(1.-std::exp(-E2));
1880  G4double I1=R1*theS1;
1881  G4double I2=R2*theS2/theB2;
1882  //G4double I3=R3*theS3/theB3;
1883  G4double I12=I1+I2;
1884  //G4double rand=(I12+I3)*G4UniformRand();
1885  G4double rand=I12*G4UniformRand();
1886  if (rand<I1 )
1887  {
1888  G4double ran=R1*G4UniformRand();
1889  if(ran>1.) ran=1.;
1890  q2=-std::log(1.-ran)/theB1; // t-chan
1891  }
1892  else
1893  {
1894  G4double ran=R2*G4UniformRand();
1895  if(ran>1.) ran=1.;
1896  q2=lastTM+std::log(1.-ran)/theB2; // u-chan (ChEx)
1897  }
1898  }
1899  else
1900  {
1901  G4double a=tgZ+tgN;
1902  G4double E1=lastTM*(theB1+lastTM*theSS);
1903  G4double R1=(1.-std::exp(-E1));
1904  G4double tss=theSS+theSS; // for future solution of quadratic equation (imediate check)
1905  G4double tm2=lastTM*lastTM;
1906  G4double E2=lastTM*tm2*theB2; // power 3 for lowA, 5 for HighA (1st)
1907  if(a>6.5)E2*=tm2; // for heavy nuclei
1908  G4double R2=(1.-std::exp(-E2));
1909  G4double E3=lastTM*theB3;
1910  if(a>6.5)E3*=tm2*tm2*tm2; // power 1 for lowA, 7 (2nd) for HighA
1911  G4double R3=(1.-std::exp(-E3));
1912  G4double E4=lastTM*theB4;
1913  G4double R4=(1.-std::exp(-E4));
1914  G4double I1=R1*theS1;
1915  G4double I2=R2*theS2;
1916  G4double I3=R3*theS3;
1917  G4double I4=R4*theS4;
1918  G4double I12=I1+I2;
1919  G4double I13=I12+I3;
1920  G4double rand=(I13+I4)*G4UniformRand();
1921  if(rand<I1)
1922  {
1923  G4double ran=R1*G4UniformRand();
1924  if(ran>1.) ran=1.;
1925  q2=-std::log(1.-ran)/theB1;
1926  if(std::fabs(tss)>1.e-7) q2=(std::sqrt(theB1*(theB1+(tss+tss)*q2))-theB1)/tss;
1927  }
1928  else if(rand<I12)
1929  {
1930  G4double ran=R2*G4UniformRand();
1931  if(ran>1.) ran=1.;
1932  q2=-std::log(1.-ran)/theB2;
1933  if(q2<0.) q2=0.;
1934  if(a<6.5) q2=std::pow(q2,third);
1935  else q2=std::pow(q2,fifth);
1936  }
1937  else if(rand<I13)
1938  {
1939  G4double ran=R3*G4UniformRand();
1940  if(ran>1.) ran=1.;
1941  q2=-std::log(1.-ran)/theB3;
1942  if(q2<0.) q2=0.;
1943  if(a>6.5) q2=std::pow(q2,sevth);
1944  }
1945  else
1946  {
1947  G4double ran=R4*G4UniformRand();
1948  if(ran>1.) ran=1.;
1949  q2=-std::log(1.-ran)/theB4;
1950  if(a<6.5) q2=lastTM-q2; // u reduced for lightA (starts from 0)
1951  }
1952  }
1953  if(q2<0.) q2=0.;
1954  if(!(q2>=-1.||q2<=1.)) G4cout<<"*NAN*G4QNeutronElCroSect::GetExchangeT: -t="<<q2<<G4endl;
1955  if(q2>lastTM)
1956  {
1957  q2=lastTM;
1958  }
1959  return q2*GeVSQ;
1960 }
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
static constexpr double gigaelectronvolt
Definition: G4SIunits.hh:209
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4double G4ChipsNeutronElasticXS::GetHMaxT ( )

Definition at line 1986 of file G4ChipsNeutronElasticXS.cc.

1987 {
1988  static const G4double HGeVSQ=gigaelectronvolt*gigaelectronvolt/2.;
1989  return lastTM*HGeVSQ;
1990 }
static constexpr double gigaelectronvolt
Definition: G4SIunits.hh:209
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4double G4ChipsNeutronElasticXS::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 160 of file G4ChipsNeutronElasticXS.cc.

164 {
165  G4double pMom=Pt->GetTotalMomentum();
166  G4int tgN = A - tgZ;
167 
168  return GetChipsCrossSection(pMom, tgZ, tgN, 2112);
169 }
virtual G4double GetChipsCrossSection(G4double momentum, G4int Z, G4int N, G4int pdg)
int G4int
Definition: G4Types.hh:78
G4double GetTotalMomentum() const
double A(double temperature)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4bool G4ChipsNeutronElasticXS::IsIsoApplicable ( const G4DynamicParticle Pt,
G4int  Z,
G4int  A,
const G4Element elm,
const G4Material mat 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 153 of file G4ChipsNeutronElasticXS.cc.

156 {
157  return true;
158 }

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