Geant4  10.02.p02
G4ScreenedNuclearRecoil.hh
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 //
28 //
29 // $Id: G4ScreenedNuclearRecoil.hh 91266 2015-06-29 06:48:42Z gcosmo $
30 //
31 //
32 // G4ScreenedNuclearRecoil.hh,v 1.24 2008/05/01 19:58:59 marcus Exp
33 // GEANT4 tag
34 //
35 //
36 //
37 // Class Description
38 // Process for screened electromagnetic nuclear elastic scattering;
39 // Physics comes from:
40 // Marcus H. Mendenhall and Robert A. Weller,
41 // "Algorithms for the rapid computation of classical cross sections
42 // for screened Coulomb collisions "
43 // Nuclear Instruments and Methods in Physics Research B58 (1991) 11-17
44 // The only input required is a screening function phi(r/a) which is the ratio
45 // of the actual interatomic potential for two atoms with atomic numbers Z1 and
46 // Z2,
47 // to the unscreened potential Z1*Z2*e^2/r where e^2 is elm_coupling in Geant4
48 // units
49 // the actual screening tables are computed externally in a python module
50 // "screened_scattering.py"
51 // to allow very specific screening functions to be added if desired, without
52 // messing with the insides of this code.
53 //
54 // First version, April 2004, Marcus H. Mendenhall, Vanderbilt University
55 // May 1, 2008 -- Added code to allow process to have zero cross section above
56 // max energy, to coordinate with G4MSC. -- mhm
57 //
58 // Class Description - End
59 
60 
61 #ifndef G4ScreenedNuclearRecoil_h
62 #define G4ScreenedNuclearRecoil_h 1
63 
64 #include "globals.hh"
65 #include "G4VDiscreteProcess.hh"
66 #include "G4ParticleChange.hh"
67 #include "c2_function.hh"
68 
69 #include "G4PhysicalConstants.hh"
70 #include "G4SystemOfUnits.hh"
71 
72 #include <map>
73 #include <vector>
74 
76 
80 
81 typedef struct G4ScreeningTables {
85 
86 // A class for loading ScreenedCoulombCrossSections
88 {
89 public:
92 
93  static const char* CVSHeaderVers() { return
94  "G4ScreenedNuclearRecoil.hh,v 1.24 2008/05/01 19:58:59 marcus Exp GEANT4 tag ";
95  }
96  static const char* CVSFileVers();
97 };
98 
99 // A class for loading ScreenedCoulombCrossSections
101 {
102 public:
103 
108 
109  typedef std::map<G4int, G4ScreeningTables> ScreeningMap;
110 
111  // a local, fast-access mapping of a particle's Z to its full definition
112  typedef std::map<G4int, class G4ParticleDefinition *> ParticleCache;
113 
114  // LoadData is called by G4ScreenedNuclearRecoil::GetMeanFreePath
115  // It loads the data tables, builds the elemental cross-section tables.
116  virtual void LoadData(G4String screeningKey, G4int z1, G4double m1,
117  G4double recoilCutoff) = 0;
118 
119  // BuildMFPTables is called by G4ScreenedNuclearRecoil::GetMeanFreePath
120  //to build the MFP tables for each material
121  void BuildMFPTables(void); // scan the MaterialsTable and construct MFP
122  //tables
123 
124  virtual G4ScreenedCoulombCrossSection *create() = 0;
125  // a 'virtual constructor' which clones the class
127  { return &(screeningData[Z]); }
128  void SetVerbosity(G4int v) { verbosity=v; }
129 
130  // this process needs element selection weighted only by number density
132  (const G4MaterialCutsCouple* couple);
133 
134  enum { nMassMapElements=116 };
135 
137  { return z1 <= nMassMapElements ? massmap[z1] : 2.5*z1; }
138 
139  // get the mean-free-path table for the indexed material
140  const G4_c2_function * operator [] (G4int materialIndex) {
141  return MFPTables.find(materialIndex)!=MFPTables.end() ?
142  &(MFPTables[materialIndex].get()) : (G4_c2_function *)0;
143  }
144 
145 protected:
146  ScreeningMap screeningData; // screening tables for each element
147  ParticleCache targetMap;
149  std::map<G4int, G4_c2_const_ptr > sigmaMap;
150  // total cross section for each element
151  std::map<G4int, G4_c2_const_ptr > MFPTables; // MFP for each material
152 
153 private:
155 
156 };
157 
158 typedef struct G4CoulombKinematicsInfo {
165 
167 public:
168  virtual void DoCollisionStep(class G4ScreenedNuclearRecoil *master,
169  const class G4Track& aTrack, const class G4Step& aStep)=0;
171 };
172 
175 
176 public:
178  virtual void DoCollisionStep(class G4ScreenedNuclearRecoil *master,
179  const class G4Track& aTrack, const class G4Step& aStep);
180 
182  const G4ScreeningTables *screen, G4double eps, G4double beta);
183 
185 
186 protected:
187  // the c2_functions we need to do the work.
191 };
192 
194  public G4ScreenedCollisionStage {
195 
196 public:
198  virtual void DoCollisionStep(class G4ScreenedNuclearRecoil *master,
199  const class G4Track& aTrack, const class G4Step& aStep);
200  virtual ~G4SingleScatter() {}
201 };
202 
209  public G4VDiscreteProcess
210 {
211 public:
212 
214 
234 
235  G4ScreenedNuclearRecoil(const G4String& processName = "ScreenedElastic",
236  const G4String &ScreeningKey="zbl", G4bool GenerateRecoils=1,
237  G4double RecoilCutoff=100.0*eV, G4double PhysicsCutoff=10.0*eV);
238 
240  virtual ~G4ScreenedNuclearRecoil();
241 
243  virtual G4double GetMeanFreePath(const G4Track&, G4double,
244  G4ForceCondition* );
245 
247  virtual G4VParticleChange* PostStepDoIt(const G4Track& aTrack,
248  const G4Step& aStep);
251  virtual G4bool IsApplicable(const G4ParticleDefinition& aParticleType);
254  virtual void BuildPhysicsTable(const G4ParticleDefinition& aParticleType);
257  virtual void DumpPhysicsTable(const G4ParticleDefinition& aParticleType);
263 
265  G4double apsis); // return true if hard collision
266 
268 
270  G4double GetNIEL() const { return NIEL; }
271 
273  void ResetTables();
274  // clear all data tables to allow changing energy cutoff, materials, etc.
275 
286 
288 
290  std::string GetScreeningKey() const { return screeningKey; }
291 
295 
297 
300 
307 
308  void EnableRecoils(G4bool flag) { generateRecoils=flag; }
309 
312 
318 
319  void SetMFPScaling(G4double scale) { MFPScale=scale; }
320 
322  G4double GetMFPScaling() const { return MFPScale; }
323 
330 
332 
335 
340 
342 
345 
349 
351  ResetTables(); }
352 
355 
358 
359  void SetNIELPartitionFunction(const G4VNIELPartition *part);
360 
367 
368  void SetCrossSectionHardening(G4double fraction, G4double HardeningFactor){
369  hardeningFraction=fraction;
370  hardeningFactor=HardeningFactor;
371  }
372 
375 
378 
381  return currentInteractionLength; }
382 
386 
389  }
390 
392  G4int GetVerboseLevel() const { return verboseLevel; }
393 
394  std::map<G4int, G4ScreenedCoulombCrossSection*> &GetCrossSectionHandlers()
395  { return crossSectionHandlers; }
396 
397  void ClearStages(void);
398 
400  collisionStages.push_back(stage); }
401 
403 
405 
407 
411  { return static_cast<G4ParticleChange &>(*pParticleChange); }
412 
415 
416  void DepositEnergy(G4int z1, G4double a1, const G4Material *material,
417  G4double energy);
418 
419 protected:
422 
425 
436 
438  std::vector<G4ScreenedCollisionStage *> collisionStages;
439 
440  std::map<G4int, G4ScreenedCoulombCrossSection*> crossSectionHandlers;
441 
445 };
446 
447 // A customized G4CrossSectionHandler which gets its data from
448 // an external program
449 
451 {
452 public:
454 
458 
460  const G4ScreenedCoulombCrossSection &src)
462 
464 
465  virtual void LoadData(G4String screeningKey, G4int z1, G4double m1,
466  G4double recoilCutoff);
467 
469  { return new G4NativeScreenedCoulombCrossSection(*this); }
470 
471  // get a list of available keys
472  std::vector<G4String> GetScreeningKeys() const;
473 
474  typedef G4_c2_function &(*ScreeningFunc)(G4int z1, G4int z2,
475  size_t nPoints, G4double rMax, G4double *au);
476 
478  phiMap[name]=fn;
479  }
480 
481 private:
482  // this is a map used to look up screening function generators
483  std::map<std::string, ScreeningFunc> phiMap;
484 };
485 
486 G4_c2_function &ZBLScreening(G4int z1, G4int z2, size_t npoints,
487  G4double rMax, G4double *auval);
488 
489 G4_c2_function &MoliereScreening(G4int z1, G4int z2, size_t npoints,
490  G4double rMax, G4double *auval);
491 
492 G4_c2_function &LJScreening(G4int z1, G4int z2, size_t npoints,
493  G4double rMax, G4double *auval);
494 
495 G4_c2_function &LJZBLScreening(G4int z1, G4int z2, size_t npoints,
496  G4double rMax, G4double *auval);
497 
498 #endif
G4NativeScreenedCoulombCrossSection(const G4ScreenedCoulombCrossSection &src)
std::vector< G4String > GetScreeningKeys() const
const G4_c2_function * operator[](G4int materialIndex)
G4int verboseLevel
Definition: G4VProcess.hh:368
G4CoulombKinematicsInfo kinematics
G4double GetCurrentInteractionLength() const
the the interaciton length used in the last scattering.
void AddStage(G4ScreenedCollisionStage *stage)
virtual void LoadData(G4String screeningKey, G4int z1, G4double m1, G4double recoilCutoff)=0
struct G4ScreeningTables G4ScreeningTables
G4double highEnergyLimit
the energy per nucleon above which the MFP is constant
virtual void DoCollisionStep(class G4ScreenedNuclearRecoil *master, const class G4Track &aTrack, const class G4Step &aStep)
G4ParticleDefinition * recoilIon
G4bool GetAvoidNuclearReactions() const
get the flag indicating whether hadronic collisions are ignored.
G4CoulombKinematicsInfo & GetKinematics()
G4ScreenedCoulombCrossSection(const G4ScreenedCoulombCrossSection &src)
G4String name
Definition: TRTMaterials.hh:40
static const G4double a1
virtual void DumpPhysicsTable(const G4ParticleDefinition &aParticleType)
Export physics tables for persistency. Not Implemented.
G4double GetRecoilCutoff() const
get the recoil cutoff
const G4VNIELPartition * NIELPartitionFunction
static const G4double massmap[nMassMapElements+1]
static const G4double eps
virtual G4double GetMeanFreePath(const G4Track &, G4double, G4ForceCondition *)
used internally by Geant4 machinery
void DepositEnergy(G4int z1, G4double a1, const G4Material *material, G4double energy)
take the given energy, and use the material information to partition it into NIEL and ionizing energy...
struct G4CoulombKinematicsInfo G4CoulombKinematicsInfo
int G4int
Definition: G4Types.hh:78
virtual void BuildPhysicsTable(const G4ParticleDefinition &aParticleType)
Build physics tables in advance. Not Implemented.
std::map< G4int, G4ScreeningTables > ScreeningMap
Provides the headers for the general c2_function algebra which supports fast, flexible operations on ...
std::vector< G4ScreenedCollisionStage * > collisionStages
virtual void DoCollisionStep(class G4ScreenedNuclearRecoil *master, const class G4Track &aTrack, const class G4Step &aStep)=0
G4_c2_function & LJScreening(G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
std::map< std::string, ScreeningFunc > phiMap
void SetPhysicsCutoff(G4double energy)
set the energy to which screening tables are computed. Typically, this is 10 eV or so...
void AllowEnergyDeposition(G4bool flag)
enable or disable all energy deposition by this process
c2_const_plugin_function_p< G4double > & phifunc
G4ScreenedCoulombCrossSection * externalCrossSectionConstructor
virtual G4ScreenedCoulombCrossSection * create()
double A(double temperature)
G4ScreenedCoulombCrossSection * crossSection
G4double GetHardeningFactor() const
get the boost factor in use.
void SetRecoilCutoff(G4double energy)
set the minimum energy (per nucleon) at which recoils can be generated, and the energy (per nucleon) ...
std::map< G4int, G4_c2_const_ptr > MFPTables
std::map< G4int, class G4ParticleDefinition * > ParticleCache
void SetMFPScaling(G4double scale)
set the mean free path scaling as specified
bool G4bool
Definition: G4Types.hh:79
std::map< G4int, G4ScreenedCoulombCrossSection * > crossSectionHandlers
G4double GetMFPScaling() const
get the MFPScaling parameter
G4double currentInteractionLength
Definition: G4VProcess.hh:297
const G4ScreeningTables * GetScreening(G4int Z)
G4ParticleDefinition * SelectRandomUnweightedTarget(const G4MaterialCutsCouple *couple)
G4bool GetAllowEnergyDeposition() const
get flag indicating whether deposition is enabled
G4double lowEnergyLimit
the energy per nucleon below which the MFP is zero
std::map< G4int, G4_c2_const_ptr > sigmaMap
G4_c2_function & ZBLScreening(G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
void SetCrossSectionHardening(G4double fraction, G4double HardeningFactor)
set the cross section boost to provide faster computation of backscattering
Definition: G4Step.hh:76
G4bool GetEnableRecoils() const
find out if generation of recoils is enabled.
void SetMaxEnergyForScattering(G4double energy)
set the upper energy beyond which this process has no section
void AvoidNuclearReactions(G4bool flag)
enable or disable whether this process will skip collisions which are close enough they need hadronic...
G4ScreenedNuclearRecoil(const G4String &processName="ScreenedElastic", const G4String &ScreeningKey="zbl", G4bool GenerateRecoils=1, G4double RecoilCutoff=100.0 *eV, G4double PhysicsCutoff=10.0 *eV)
Construct the process and set some physics parameters for it.
virtual G4bool IsApplicable(const G4ParticleDefinition &aParticleType)
test if a prticle of type aParticleType can use this process
void SetExternalCrossSectionHandler(G4ScreenedCoulombCrossSection *cs)
set a function to compute screening tables, if the user needs non-standard behavior.
void SetNIELPartitionFunction(const G4VNIELPartition *part)
set the pointer to a class for paritioning energy into NIEL
virtual G4ScreenedCoulombCrossSection * create()=0
static const double eV
Definition: G4SIunits.hh:212
c2_function< G4double > G4_c2_function
G4double energy(const ThreeVector &p, const G4double m)
G4_c2_function & LJZBLScreening(G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
std::string GetScreeningKey() const
find out what screening function we are using
c2_const_ptr< G4double > G4_c2_const_ptr
virtual ~G4ScreenedNuclearRecoil()
destructor
virtual G4VParticleChange * PostStepDoIt(const G4Track &aTrack, const G4Step &aStep)
used internally by Geant4 machinery
virtual G4bool CheckNuclearCollision(G4double A, G4double A1, G4double apsis)
deterine if the moving particle is within the strong force of the selected nucleus ...
void AddScreeningFunction(G4String name, ScreeningFunc fn)
c2_ptr< G4double > G4_c2_ptr
virtual G4ScreenedCoulombCrossSection * GetNewCrossSectionHandler(void)
G4double GetHardeningFraction() const
get the fraction of particles which will have boosted scattering
G4double GetPhysicsCutoff() const
get the physics cutoff energy.
virtual void DoCollisionStep(class G4ScreenedNuclearRecoil *master, const class G4Track &aTrack, const class G4Step &aStep)
double G4double
Definition: G4Types.hh:76
G4int GetVerboseLevel() const
get the verbosity.
G4_c2_function &(* ScreeningFunc)(G4int z1, G4int z2, size_t nPoints, G4double rMax, G4double *au)
G4double GetNIEL() const
Get non-ionizing energy loss for last step.
G4ForceCondition
class G4ParticleChange & GetParticleChange()
get the pointer to our ParticleChange object. for internal use, primarily.
void EnableRecoils(G4bool flag)
enable or disable the generation of recoils. recoils are disabled, the energy they would have receiv...
G4bool DoScreeningComputation(class G4ScreenedNuclearRecoil *master, const G4ScreeningTables *screen, G4double eps, G4double beta)
G4_c2_function & MoliereScreening(G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
void ResetTables()
clear precomputed screening tables
A process which handles screened Coulomb collisions between nuclei.
virtual void LoadData(G4String screeningKey, G4int z1, G4double m1, G4double recoilCutoff)
G4NativeScreenedCoulombCrossSection(const G4NativeScreenedCoulombCrossSection &src)
G4double processMaxEnergy
the energy per nucleon beyond which the cross section is zero, to cross over to G4MSC ...
std::map< G4int, G4ScreenedCoulombCrossSection * > & GetCrossSectionHandlers()