Geant4_10
G4IonisParamMat.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4IonisParamMat.cc 76883 2013-11-18 12:50:08Z gcosmo $
27 //
28 //
29 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
30 
31 // 09-07-98, data moved from G4Material, M.Maire
32 // 18-07-98, bug corrected in ComputeDensityEffect() for gas
33 // 16-01-01, bug corrected in ComputeDensityEffect() E100eV (L.Urban)
34 // 08-02-01, fShellCorrectionVector correctly handled (mma)
35 // 28-10-02, add setMeanExcitationEnergy (V.Ivanchenko)
36 // 06-09-04, factor 2 to shell correction term (V.Ivanchenko)
37 // 10-05-05, add a missing coma in FindMeanExcitationEnergy() - Bug#746 (mma)
38 // 27-09-07, add computation of parameters for ions (V.Ivanchenko)
39 // 04-03-08, remove reference to G4NistManager. Add fBirks constant (mma)
40 // 30-10-09, add G4DensityEffectData class and density effect computation (VI)
41 
42 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
43 
44 #include "G4IonisParamMat.hh"
45 #include "G4Material.hh"
46 #include "G4DensityEffectData.hh"
47 #include "G4NistManager.hh"
48 #include "G4Pow.hh"
49 #include "G4PhysicalConstants.hh"
50 #include "G4SystemOfUnits.hh"
51 #include "G4Log.hh"
52 #include "G4Exp.hh"
53 
54 G4DensityEffectData* G4IonisParamMat::fDensityData = 0;
55 
56 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
57 
59  : fMaterial(material)
60 {
61  fBirks = 0.;
62  fMeanEnergyPerIon = 0.0;
63  twoln10 = 2.*G4Pow::GetInstance()->logZ(10);
64 
65  // minimal set of default parameters for density effect
66  fCdensity = 0.0;
67  fD0density = 0.0;
68  fAdjustmentFactor = 1.0;
69  if(!fDensityData) { fDensityData = new G4DensityEffectData(); }
70 
71  // compute parameters
72  ComputeMeanParameters();
73  ComputeDensityEffect();
74  ComputeFluctModel();
75  ComputeIonParameters();
76 }
77 
78 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
79 
80 // Fake default constructor - sets only member data and allocates memory
81 // for usage restricted to object persistency
82 
84  : fMaterial(0), fShellCorrectionVector(0)
85 {
86  fMeanExcitationEnergy = 0.0;
87  fLogMeanExcEnergy = 0.0;
88  fTaul = 0.0;
89  fCdensity = 0.0;
90  fMdensity = 0.0;
91  fAdensity = 0.0;
92  fX0density = 0.0;
93  fX1density = 0.0;
94  fD0density = 0.0;
95  fPlasmaEnergy = 0.0;
96  fAdjustmentFactor = 0.0;
97  fF1fluct = 0.0;
98  fF2fluct = 0.0;
99  fEnergy1fluct = 0.0;
100  fLogEnergy1fluct = 0.0;
101  fEnergy2fluct = 0.0;
102  fLogEnergy2fluct = 0.0;
103  fEnergy0fluct = 0.0;
104  fRateionexcfluct = 0.0;
105  fZeff = 0.0;
106  fFermiEnergy = 0.0;
107  fLfactor = 0.0;
108  fInvA23 = 0.0;
109  fBirks = 0.0;
110  fMeanEnergyPerIon = 0.0;
111  twoln10 = 2.*G4Pow::GetInstance()->logZ(10);
112 }
113 
114 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
115 
117 {
118  if (fShellCorrectionVector) { delete [] fShellCorrectionVector; }
119  if (fDensityData) { delete fDensityData; }
120  fDensityData = 0;
121  fShellCorrectionVector = 0;
122 }
123 
124 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
125 
126 void G4IonisParamMat::ComputeMeanParameters()
127 {
128  // compute mean excitation energy and shell correction vector
129  fTaul = (*(fMaterial->GetElementVector()))[0]->GetIonisation()->GetTaul();
130 
131  fMeanExcitationEnergy = 0.;
132  fLogMeanExcEnergy = 0.;
133 
134  size_t nElements = fMaterial->GetNumberOfElements();
135  const G4ElementVector* elmVector = fMaterial->GetElementVector();
136  const G4double* nAtomsPerVolume = fMaterial->GetVecNbOfAtomsPerVolume();
137 
138  G4String ch = fMaterial->GetChemicalFormula();
139 
140  if(ch != "") { fMeanExcitationEnergy = FindMeanExcitationEnergy(ch); }
141 
142  // Chemical formula defines mean excitation energy
143  if(fMeanExcitationEnergy > 0.0) {
144  fLogMeanExcEnergy = G4Log(fMeanExcitationEnergy);
145 
146  // Compute average
147  } else {
148  for (size_t i=0; i < nElements; i++) {
149  const G4Element* elm = (*elmVector)[i];
150  fLogMeanExcEnergy += nAtomsPerVolume[i]*elm->GetZ()
152  }
153  fLogMeanExcEnergy /= fMaterial->GetTotNbOfElectPerVolume();
154  fMeanExcitationEnergy = G4Exp(fLogMeanExcEnergy);
155  }
156 
157  fShellCorrectionVector = new G4double[3];
158 
159  for (G4int j=0; j<=2; j++)
160  {
161  fShellCorrectionVector[j] = 0.;
162 
163  for (size_t k=0; k<nElements; k++) {
164  fShellCorrectionVector[j] += nAtomsPerVolume[k]
165  *(((*elmVector)[k])->GetIonisation()->GetShellCorrectionVector())[j];
166  }
167  fShellCorrectionVector[j] *= 2.0/fMaterial->GetTotNbOfElectPerVolume();
168  }
169 }
170 
171 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
172 
174 {
175  return fDensityData;
176 }
177 
178 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
179 
180 void G4IonisParamMat::ComputeDensityEffect()
181 {
182  G4State State = fMaterial->GetState();
183 
184  // Check if density effect data exist in the table
185  // R.M. Sternheimer, Atomic Data and Nuclear Data Tables, 30: 261 (1984)
186  G4int idx = fDensityData->GetIndex(fMaterial->GetName());
187  G4int nelm= fMaterial->GetNumberOfElements();
188  G4int Z0 = G4lrint((*(fMaterial->GetElementVector()))[0]->GetZ());
189  if(idx < 0 && 1 == nelm) {
190  idx = fDensityData->GetElementIndex(Z0, fMaterial->GetState());
191  }
192 
193  //G4cout<<"DensityEffect for "<<fMaterial->GetName()<<" "<< idx << G4endl;
194 
195  if(idx >= 0) {
196 
197  // Take parameters for the density effect correction from
198  // R.M. Sternheimer et al. Density Effect For The Ionization Loss
199  // of Charged Particles in Various Substances.
200  // Atom. Data Nucl. Data Tabl. 30 (1984) 261-271.
201 
202  fCdensity = fDensityData->GetCdensity(idx);
203  fMdensity = fDensityData->GetMdensity(idx);
204  fAdensity = fDensityData->GetAdensity(idx);
205  fX0density = fDensityData->GetX0density(idx);
206  fX1density = fDensityData->GetX1density(idx);
207  fD0density = fDensityData->GetDelta0density(idx);
208  fPlasmaEnergy = fDensityData->GetPlasmaEnergy(idx);
209  fAdjustmentFactor = fDensityData->GetAdjustmentFactor(idx);
210 
211  // Correction for base material
212  const G4Material* bmat = fMaterial->GetBaseMaterial();
213  if(bmat) {
214  G4double corr = G4Log(bmat->GetDensity()/fMaterial->GetDensity());
215  fCdensity += corr;
216  fX0density += corr/twoln10;
217  fX1density += corr/twoln10;
218  }
219 
220  } else {
221 
222  static const G4double Cd2 = 4*pi*hbarc_squared*classic_electr_radius;
223  fPlasmaEnergy = std::sqrt(Cd2*fMaterial->GetTotNbOfElectPerVolume());
224 
225  // Compute parameters for the density effect correction in DE/Dx formula.
226  // The parametrization is from R.M. Sternheimer, Phys. Rev.B,3:3681 (1971)
227  G4int icase;
228 
229  fCdensity = 1. + 2*G4Log(fMeanExcitationEnergy/fPlasmaEnergy);
230  //
231  // condensed materials
232  //
233  if ((State == kStateSolid)||(State == kStateLiquid)) {
234 
235  static const G4double E100eV = 100.*eV;
236  static const G4double ClimiS[] = {3.681 , 5.215 };
237  static const G4double X0valS[] = {1.0 , 1.5 };
238  static const G4double X1valS[] = {2.0 , 3.0 };
239 
240  if(fMeanExcitationEnergy < E100eV) { icase = 0; }
241  else { icase = 1; }
242 
243  if(fCdensity < ClimiS[icase]) { fX0density = 0.2; }
244  else { fX0density = 0.326*fCdensity - X0valS[icase]; }
245 
246  fX1density = X1valS[icase]; fMdensity = 3.0;
247 
248  //special: Hydrogen
249  if (1 == nelm && 1 == Z0) {
250  fX0density = 0.425; fX1density = 2.0; fMdensity = 5.949;
251  }
252  }
253  //
254  // gases
255  //
256  if (State == kStateGas) {
257 
258  fMdensity = 3.;
259  fX1density = 4.0;
260  //static const G4double ClimiG[] = { 10. , 10.5 , 11. , 11.5 , 12.25 , 13.804};
261  //static const G4double X0valG[] = { 1.6 , 1.7 , 1.8 , 1.9 , 2.0 , 2.0 };
262  //static const G4double X1valG[] = { 4.0 , 4.0 , 4.0 , 4.0 , 4.0 , 5.0 };
263 
264  if(fCdensity < 10.) {
265  fX0density = 1.6;
266  } else if(fCdensity < 11.5) {
267  fX0density = 1.6 + 0.2*(fCdensity - 10.);
268  } else if(fCdensity < 12.25) {
269  fX0density = 1.9 + (fCdensity - 11.5)/7.5;
270  } else if(fCdensity < 13.804) {
271  fX0density = 2.0;
272  fX1density = 4.0 + (fCdensity - 12.25)/1.554;
273  } else {
274  fX0density = 0.326*fCdensity-2.5; fX1density = 5.0;
275  }
276 
277  //special: Hydrogen
278  if (1 == nelm && 1 == Z0) {
279  fX0density = 1.837; fX1density = 3.0; fMdensity = 4.754;
280  }
281 
282  //special: Helium
283  if (1 == nelm && 2 == Z0) {
284  fX0density = 2.191; fX1density = 3.0; fMdensity = 3.297;
285  }
286  }
287  }
288 
289  // change parameters if the gas is not in STP.
290  // For the correction the density(STP) is needed.
291  // Density(STP) is calculated here :
292 
293 
294  if (State == kStateGas) {
295  G4double Density = fMaterial->GetDensity();
296  G4double Pressure = fMaterial->GetPressure();
297  G4double Temp = fMaterial->GetTemperature();
298 
299  G4double DensitySTP = Density*STP_Pressure*Temp/(Pressure*STP_Temperature);
300 
301  G4double ParCorr = G4Log(Density/DensitySTP);
302 
303  fCdensity -= ParCorr;
304  fX0density -= ParCorr/twoln10;
305  fX1density -= ParCorr/twoln10;
306  }
307 
308  // fAdensity parameter can be fixed for not conductive materials
309  if(0.0 == fD0density) {
310  G4double Xa = fCdensity/twoln10;
311  fAdensity = twoln10*(Xa-fX0density)
312  /std::pow((fX1density-fX0density),fMdensity);
313  }
314  /*
315  G4cout << "G4IonisParamMat: density effect data for <" << fMaterial->GetName()
316  << "> " << G4endl;
317  G4cout << "Eplasma(eV)= " << fPlasmaEnergy/eV
318  << " rho= " << fAdjustmentFactor
319  << " -C= " << fCdensity
320  << " x0= " << fX0density
321  << " x1= " << fX1density
322  << " a= " << fAdensity
323  << " m= " << fMdensity
324  << G4endl;
325  */
326 }
327 
328 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
329 
330 void G4IonisParamMat::ComputeFluctModel()
331 {
332  // compute parameters for the energy loss fluctuation model
333  // needs an 'effective Z'
334  G4double Zeff = 0.;
335  for (size_t i=0;i<fMaterial->GetNumberOfElements();i++) {
336  Zeff += (fMaterial->GetFractionVector())[i]
337  *((*(fMaterial->GetElementVector()))[i]->GetZ());
338  }
339  if (Zeff > 2.) { fF2fluct = 2./Zeff; }
340  else { fF2fluct = 0.; }
341 
342  fF1fluct = 1. - fF2fluct;
343  fEnergy2fluct = 10.*Zeff*Zeff*eV;
344  fLogEnergy2fluct = G4Log(fEnergy2fluct);
345  fLogEnergy1fluct = (fLogMeanExcEnergy - fF2fluct*fLogEnergy2fluct)
346  /fF1fluct;
347  fEnergy1fluct = G4Exp(fLogEnergy1fluct);
348  fEnergy0fluct = 10.*eV;
349  fRateionexcfluct = 0.4;
350 }
351 
352 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
353 
354 void G4IonisParamMat::ComputeIonParameters()
355 {
356  // get elements in the actual material,
357  const G4ElementVector* theElementVector = fMaterial->GetElementVector() ;
358  const G4double* theAtomicNumDensityVector =
359  fMaterial->GetAtomicNumDensityVector() ;
360  const G4int NumberOfElements = fMaterial->GetNumberOfElements() ;
361 
362  // loop for the elements in the material
363  // to find out average values Z, vF, lF
364  G4double z(0.0), vF(0.0), lF(0.0), norm(0.0), a23(0.0);
365 
366  if( 1 == NumberOfElements ) {
367  const G4Element* element = (*theElementVector)[0];
368  z = element->GetZ();
369  vF= element->GetIonisation()->GetFermiVelocity();
370  lF= element->GetIonisation()->GetLFactor();
371  a23 = 1.0/G4Pow::GetInstance()->Z23(G4int(element->GetN()));
372 
373  } else {
374  for (G4int iel=0; iel<NumberOfElements; iel++)
375  {
376  const G4Element* element = (*theElementVector)[iel];
377  const G4double weight = theAtomicNumDensityVector[iel];
378  norm += weight ;
379  z += element->GetZ() * weight;
380  vF += element->GetIonisation()->GetFermiVelocity() * weight;
381  lF += element->GetIonisation()->GetLFactor() * weight;
382  a23 += weight/G4Pow::GetInstance()->Z23(G4int(element->GetN()));
383  }
384  z /= norm;
385  vF /= norm;
386  lF /= norm;
387  a23 /= norm;
388  }
389  fZeff = z;
390  fLfactor = lF;
391  fFermiEnergy = 25.*keV*vF*vF;
392  fInvA23 = a23;
393 }
394 
395 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
396 
398 {
399  if(value == fMeanExcitationEnergy || value <= 0.0) { return; }
400  if (G4NistManager::Instance()->GetVerbose() > 1) {
401  G4cout << "G4Material: Mean excitation energy is changed for "
402  << fMaterial->GetName()
403  << " Iold= " << fMeanExcitationEnergy/eV
404  << "eV; Inew= " << value/eV << " eV;"
405  << G4endl;
406  }
407 
408  fMeanExcitationEnergy = value;
409 
410  // add corrections to density effect
411  G4double newlog = G4Log(value);
412  G4double corr = 2*(newlog - fLogMeanExcEnergy);
413  fCdensity += corr;
414  fX0density += corr/twoln10;
415  fX1density += corr/twoln10;
416 
417  // recompute parameters of fluctuation model
418  fLogMeanExcEnergy = newlog;
419  ComputeFluctModel();
420 }
421 
422 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
423 
425 {
426  // The data on mean excitation energy for compaunds
427  // from "Stopping Powers for Electrons and Positrons"
428  // ICRU Report N#37, 1984 (energy in eV)
429 
430  size_t numberOfMolecula = 54;
431  static const G4String name[54] = {
432  // gas 0 - 12
433  "NH_3", "C_4H_10", "CO_2", "C_2H_6", "C_7H_16-Gas",
434  // "G4_AMMONIA", "G4_BUTANE","G4_CARBON_DIOXIDE","G4_ETHANE", "G4_N-HEPTANE"
435  "C_6H_14-Gas", "CH_4", "NO", "N_2O", "C_8H_18-Gas",
436  // "G4_N-HEXANE" , "G4_METHANE", "x", "G4_NITROUS_OXIDE", "G4_OCTANE"
437  "C_5H_12-Gas", "C_3H_8", "H_2O-Gas",
438  // "G4_N-PENTANE", "G4_PROPANE", "G4_WATER_VAPOR"
439 
440  // liquid 13 - 39
441  "C_3H_6O", "C_6H_5NH_2", "C_6H_6", "C_4H_9OH", "CCl_4",
442  //"G4_ACETONE","G4_ANILINE","G4_BENZENE","G4_N-BUTYL_ALCOHOL","G4_CARBON_TETRACHLORIDE"
443  "C_6H_5Cl", "CHCl_3", "C_6H_12", "C_6H_4Cl_2", "C_4Cl_2H_8O",
444  //"G4_CHLOROBENZENE","G4_CHLOROFORM","G4_CYCLOHEXANE","G4_1,2-DICHLOROBENZENE","G4_DICHLORODIETHYL_ETHER"
445  "C_2Cl_2H_4", "(C_2H_5)_2O", "C_2H_5OH", "C_3H_5(OH)_3","C_7H_16",
446  //"G4_1,2-DICHLOROETHANE","G4_DIETHYL_ETHER","G4_ETHYL_ALCOHOL","G4_GLYCEROL","G4_N-HEPTANE"
447  "C_6H_14", "CH_3OH", "C_6H_5NO_2","C_5H_12", "C_3H_7OH",
448  //"G4_N-HEXANE","G4_METHANOL","G4_NITROBENZENE","G4_N-PENTANE","G4_N-PROPYL_ALCOHOL",
449  "C_5H_5N", "C_8H_8", "C_2Cl_4", "C_7H_8", "C_2Cl_3H",
450  //"G4_PYRIDINE","G4_POLYSTYRENE","G4_TETRACHLOROETHYLENE","G4_TOLUENE","G4_TRICHLOROETHYLENE"
451  "H_2O", "C_8H_10",
452  // "G4_WATER", "G4_XYLENE"
453 
454  // solid 40 - 53
455  "C_5H_5N_5", "C_5H_5N_5O", "(C_6H_11NO)-nylon", "C_25H_52",
456  // "G4_ADENINE", "G4_GUANINE", "G4_NYLON-6-6", "G4_PARAFFIN"
457  "(C_2H_4)-Polyethylene", "(C_5H_8O_2)-Polymethil_Methacrylate",
458  // "G4_ETHYLENE", "G4_PLEXIGLASS"
459  "(C_8H_8)-Polystyrene", "A-150-tissue", "Al_2O_3", "CaF_2",
460  // "G4_POLYSTYRENE", "G4_A-150_TISSUE", "G4_ALUMINUM_OXIDE", "G4_CALCIUM_FLUORIDE"
461  "LiF", "Photo_Emulsion", "(C_2F_4)-Teflon", "SiO_2"
462  // "G4_LITHIUM_FLUORIDE", "G4_PHOTO_EMULSION", "G4_TEFLON", "G4_SILICON_DIOXIDE"
463  } ;
464 
465  static const G4double meanExcitation[54] = {
466 
467  53.7, 48.3, 85.0, 45.4, 49.2,
468  49.1, 41.7, 87.8, 84.9, 49.5,
469  48.2, 47.1, 71.6,
470 
471  64.2, 66.2, 63.4, 59.9, 166.3,
472  89.1, 156.0, 56.4, 106.5, 103.3,
473  111.9, 60.0, 62.9, 72.6, 54.4,
474  54.0, 67.6, 75.8, 53.6, 61.1,
475  66.2, 64.0, 159.2, 62.5, 148.1,
476  75.0, 61.8,
477 
478  71.4, 75.0, 63.9, 48.3, 57.4,
479  74.0, 68.7, 65.1, 145.2, 166.,
480  94.0, 331.0, 99.1, 139.2
481  };
482 
483  G4double x = fMeanExcitationEnergy;
484 
485  for(size_t i=0; i<numberOfMolecula; i++) {
486  if(chFormula == name[i]) {
487  x = meanExcitation[i]*eV;
488  break;
489  }
490  }
491  return x;
492 }
493 
494 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
495 
497 {
498  fShellCorrectionVector = 0;
499  fMaterial = 0;
500  *this = right;
501 }
502 
503 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
504 
506 {
507  if (this != &right)
508  {
509  fMaterial = right.fMaterial;
510  fMeanExcitationEnergy = right.fMeanExcitationEnergy;
511  fLogMeanExcEnergy = right.fLogMeanExcEnergy;
512  if(fShellCorrectionVector){ delete [] fShellCorrectionVector; }
513  fShellCorrectionVector = new G4double[3];
514  fShellCorrectionVector[0] = right.fShellCorrectionVector[0];
515  fShellCorrectionVector[1] = right.fShellCorrectionVector[1];
516  fShellCorrectionVector[2] = right.fShellCorrectionVector[2];
517  fTaul = right.fTaul;
518  fCdensity = right.fCdensity;
519  fMdensity = right.fMdensity;
520  fAdensity = right.fAdensity;
521  fX0density = right.fX0density;
522  fX1density = right.fX1density;
523  fD0density = right.fD0density;
524  fPlasmaEnergy = right.fPlasmaEnergy;
525  fAdjustmentFactor = right.fAdjustmentFactor;
526  fF1fluct = right.fF1fluct;
527  fF2fluct = right.fF2fluct;
528  fEnergy1fluct = right.fEnergy1fluct;
529  fLogEnergy1fluct = right.fLogEnergy1fluct;
530  fEnergy2fluct = right.fEnergy2fluct;
531  fLogEnergy2fluct = right.fLogEnergy2fluct;
532  fEnergy0fluct = right.fEnergy0fluct;
533  fRateionexcfluct = right.fRateionexcfluct;
534  fZeff = right.fZeff;
535  fFermiEnergy = right.fFermiEnergy;
536  fLfactor = right.fLfactor;
537  fInvA23 = right.fInvA23;
538  fBirks = right.fBirks;
539  fMeanEnergyPerIon = right.fMeanEnergyPerIon;
540  fDensityData = right.fDensityData;
541  twoln10 = right.twoln10;
542  }
543  return *this;
544 }
545 
546 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
547 
549 {
550  return (this == (G4IonisParamMat*) &right);
551 }
552 
553 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
554 
556 {
557  return (this != (G4IonisParamMat*) &right);
558 }
559 
560 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... ....oooOO0OOooo....
561 
G4double GetPressure() const
Definition: G4Material.hh:181
static G4Pow * GetInstance()
Definition: G4Pow.cc:53
G4IonisParamMat & operator=(const G4IonisParamMat &)
G4double GetTotNbOfElectPerVolume() const
Definition: G4Material.hh:210
std::vector< G4Element * > G4ElementVector
G4double GetN() const
Definition: G4Element.hh:134
G4int GetElementIndex(G4int Z, G4State mState)
G4State
Definition: G4Material.hh:114
G4double FindMeanExcitationEnergy(const G4String &chFormula)
const G4String & GetChemicalFormula() const
Definition: G4Material.hh:177
void SetMeanExcitationEnergy(G4double value)
const G4String & GetName() const
Definition: G4Material.hh:176
Float_t norm
G4double GetZ() const
Definition: G4Element.hh:131
G4double GetDensity() const
Definition: G4Material.hh:178
const XML_Char * name
Definition: expat.h:151
tuple x
Definition: test.py:50
float STP_Temperature
Definition: hepunit.py:302
G4double GetLFactor() const
double weight
Definition: plottest35.C:25
#define State(theXInfo)
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:188
int G4int
Definition: G4Types.hh:78
static G4NistManager * Instance()
string material
Definition: eplot.py:19
G4double logZ(G4int Z) const
Definition: G4Pow.hh:165
G4double GetFermiVelocity() const
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:204
G4double GetMeanExcitationEnergy() const
G4GLOB_DLL std::ostream G4cout
G4double GetX1density(G4int idx)
G4int GetIndex(const G4String &matName)
virtual ~G4IonisParamMat()
G4double GetPlasmaEnergy(G4int idx)
const G4double * GetAtomicNumDensityVector() const
Definition: G4Material.hh:214
static G4DensityEffectData * GetDensityEffectData()
G4double G4Log(G4double x)
Definition: G4Log.hh:227
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:180
G4double GetAdensity(G4int idx)
int G4lrint(double ad)
Definition: templates.hh:163
G4IonisParamMat(G4Material *)
G4IonisParamElm * GetIonisation() const
Definition: G4Element.hh:198
const G4Material * GetBaseMaterial() const
Definition: G4Material.hh:231
tuple z
Definition: test.py:28
G4double Z23(G4int Z) const
Definition: G4Pow.hh:153
G4int operator!=(const G4IonisParamMat &) const
G4double GetCdensity(G4int idx)
G4double GetTemperature() const
Definition: G4Material.hh:180
const XML_Char int const XML_Char * value
Definition: expat.h:331
#define G4endl
Definition: G4ios.hh:61
G4double GetMdensity(G4int idx)
size_t GetNumberOfElements() const
Definition: G4Material.hh:184
G4State GetState() const
Definition: G4Material.hh:179
double G4double
Definition: G4Types.hh:76
G4double GetDelta0density(G4int idx)
const G4double * GetFractionVector() const
Definition: G4Material.hh:192
G4int operator==(const G4IonisParamMat &) const
G4double GetAdjustmentFactor(G4int idx)
G4double GetX0density(G4int idx)
int STP_Pressure
Definition: hepunit.py:303