Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4RadioactiveDecayBase.hh
이 파일의 문서화 페이지로 가기
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 //
27 // //
28 // File: G4RadioactiveDecayBase.hh //
29 // Author: D.H. Wright (SLAC) //
30 // Date: 9 August 2017 //
31 // Description: version the G4RadioactiveDecay process by F. Lei and //
32 // P.R. Truscott with biasing and activation calculations //
33 // removed to a derived class. It performs alpha, beta, //
34 // electron capture and isomeric transition decays of //
35 // radioactive nuclei. //
36 // //
38 
39 #ifndef G4RadioactiveDecayBase_h
40 #define G4RadioactiveDecayBase_h 1
41 
42 
43 #include <vector>
44 #include <map>
46 
47 #include "G4ios.hh"
48 #include "globals.hh"
51 
52 #include "G4NucleusLimits.hh"
53 #include "G4ThreeVector.hh"
54 #include "G4Threading.hh"
55 
56 class G4Fragment;
59 
60 typedef std::map<G4String, G4DecayTable*> DecayTableMap;
61 
62 
64 {
65  // class description
66 
67  // Implementation of the radioactive decay process which simulates the
68  // decays of radioactive nuclei. These nuclei are submitted to RDM as
69  // G4Ions. The required half-lives and decay schemes are retrieved from
70  // the Radioactivity database which was derived from ENSDF.
71  // All decay products are submitted back to the particle tracking process
72  // through the G4ParticleChangeForRadDecay object.
73  // class description - end
74 
75  public: // with description
76 
77  G4RadioactiveDecayBase(const G4String& processName="RadioactiveDecayBase");
79 
80  virtual void ProcessDescription(std::ostream& outFile) const;
81 
83  // Return true if the specified isotope is
84  // 1) defined as "nucleus" and
85  // 2) it is within theNucleusLimit
86 
87  // Return decay table if it exists, if not, load it from file
89 
90  // Select a logical volume in which RDM applies
91  void SelectAVolume(const G4String aVolume);
92 
93  // Remove a logical volume from the RDM applied list
94  void DeselectAVolume(const G4String aVolume);
95 
96  // Select all logical volumes for the application of RDM
97  void SelectAllVolumes();
98 
99  // Remove all logical volumes from RDM applications
100  void DeselectAllVolumes();
101 
102  // Enable/disable ICM
103  void SetICM(G4bool icm) {applyICM = icm;}
104 
105  // Enable/disable ARM
106  void SetARM(G4bool arm) {applyARM = arm;}
107 
108  G4DecayTable* LoadDecayTable(const G4ParticleDefinition& theParentNucleus);
109  // Load the decay data of isotope theParentNucleus
110 
111  void AddUserDecayDataFile(G4int Z, G4int A,G4String filename);
112  // Allow the user to replace the radio-active decay data provided in Geant4
113  // by its own data file for a given isotope
114 
116  // Sets the VerboseLevel which controls duggering display
117 
118  inline G4int GetVerboseLevel() const {return verboseLevel;}
119  // Returns the VerboseLevel which controls level of debugging output
120 
121  inline void SetNucleusLimits(G4NucleusLimits theNucleusLimits1)
122  {theNucleusLimits = theNucleusLimits1 ;}
123  // Sets theNucleusLimits which specifies the range of isotopes
124  // the G4RadioactiveDecay applies.
125 
126  // Returns theNucleusLimits which specifies the range of isotopes used
127  // by G4RadioactiveDecay
129 
130  // Controls whether G4RadioactiveDecay uses fast beta simulation mode
131  // Currently does nothing - kept for backward compatibility
132  inline void SetFBeta (G4bool r ) { FBeta = r; }
133 
134  inline void SetDecayDirection(const G4ThreeVector& theDir) {
135  forceDecayDirection = theDir.unit();
136  }
137 
138  inline const G4ThreeVector& GetDecayDirection() const {
139  return forceDecayDirection;
140  }
141 
142  inline void SetDecayHalfAngle(G4double halfAngle=0.*CLHEP::deg) {
144  }
145 
147 
148  // Force direction (random within half-angle) for "visible" daughters
149  // (applies to electrons, positrons, gammas, neutrons, protons or alphas)
150  inline void SetDecayCollimation(const G4ThreeVector& theDir,
151  G4double halfAngle = 0.*CLHEP::deg) {
152  SetDecayDirection(theDir);
153  SetDecayHalfAngle(halfAngle);
154  }
155 
157 
158  G4VParticleChange* DecayIt(const G4Track& theTrack,
159  const G4Step& theStep);
160 
161  protected:
162 
163  G4DecayProducts* DoDecay(const G4ParticleDefinition& theParticleDef);
164 
165  // Apply directional bias for "visible" daughters (e+-, gamma, n, p, alpha)
166  void CollimateDecay(G4DecayProducts* products);
169 
170  G4double GetMeanFreePath(const G4Track& theTrack, G4double previousStepSize,
172 
173  G4double GetMeanLifeTime(const G4Track& theTrack,
175 
176 // G4double GetDecayTime();
177 
178  // ParticleChange for decay process
180 
183 
184  std::vector<G4String> ValidVolumes;
186 
187  static const G4double levelTolerance;
188 
189  private:
190 
193 
195 
198 
201 
202  // Parameters for pre-collimated (biased) decay products
205  static const G4ThreeVector origin; // (0,0,0) for convenience
206 
207  // Radioactive decay database directory path
209 
210  //User define radioactive decay data files replacing some files in the G4RADECAY database
211  std::map<G4int, G4String> theUserRadioactiveDataFiles;
212 
213  // Library of decay tables
215 #ifdef G4MULTITHREADED
216  static DecayTableMap* master_dkmap;
217 #endif
218 
219  // Remainder of life time at rest
222 
223 
224  // inline implementations
225  inline
228  {
231  return fRemainderLifeTime;
232  }
233 
234  inline
236  const G4Step& theStep)
237  {return DecayIt(theTrack, theStep);}
238 
239  inline
241  const G4Step& theStep)
242  {return DecayIt(theTrack, theStep);}
243 
244 #ifdef G4MULTITHREADED
245  public:
246  static G4Mutex radioactiveDecayMutex;
247 #endif
248 };
249 
250 #endif
251 
virtual void ProcessDescription(std::ostream &outFile) const
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4DecayTable * GetDecayTable(const G4ParticleDefinition *)
void BuildPhysicsTable(const G4ParticleDefinition &)
G4ThreeVector ChooseCollimationDirection() const
void SetDecayCollimation(const G4ThreeVector &theDir, G4double halfAngle=0.*CLHEP::deg)
G4RadioactiveDecayBaseMessenger * theRadioactiveDecayBaseMessenger
G4PhotonEvaporation * photonEvaporation
void SetDecayHalfAngle(G4double halfAngle=0.*CLHEP::deg)
G4VParticleChange * DecayIt(const G4Track &theTrack, const G4Step &theStep)
G4double condition(const G4ErrorSymMatrix &m)
G4double GetMeanLifeTime(const G4Track &theTrack, G4ForceCondition *condition)
void SetDecayDirection(const G4ThreeVector &theDir)
G4double AtRestGetPhysicalInteractionLength(const G4Track &track, G4ForceCondition *condition)
G4RadioactiveDecayBase & operator=(const G4RadioactiveDecayBase &right)
std::map< G4int, G4String > theUserRadioactiveDataFiles
Float_t Z
void SelectAVolume(const G4String aVolume)
static const G4ThreeVector origin
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
const XML_Char int const XML_Char * value
Definition: expat.h:331
void CollimateDecay(G4DecayProducts *products)
G4NucleusLimits GetNucleusLimits() const
std::vector< G4String > ValidVolumes
void SetNucleusLimits(G4NucleusLimits theNucleusLimits1)
G4DecayProducts * DoDecay(const G4ParticleDefinition &theParticleDef)
double A(double temperature)
G4VParticleChange * AtRestDoIt(const G4Track &theTrack, const G4Step &theStep)
Definition: G4Step.hh:76
G4ParticleChangeForRadDecay fParticleChangeForRadDecay
G4bool IsApplicable(const G4ParticleDefinition &)
void DeselectAVolume(const G4String aVolume)
Hep3Vector unit() const
static const G4double levelTolerance
G4VParticleChange * PostStepDoIt(const G4Track &theTrack, const G4Step &theStep)
int G4int
Definition: G4Types.hh:78
G4double GetDecayHalfAngle() const
void AddUserDecayDataFile(G4int Z, G4int A, G4String filename)
G4ForceCondition
void CollimateDecayProduct(G4DynamicParticle *product)
static constexpr double deg
virtual G4double AtRestGetPhysicalInteractionLength(const G4Track &, G4ForceCondition *)
std::map< G4String, G4DecayTable * > DecayTableMap
G4DecayTable * LoadDecayTable(const G4ParticleDefinition &theParentNucleus)
G4RadioactiveDecayBase(const G4String &processName="RadioactiveDecayBase")
G4double GetMeanFreePath(const G4Track &theTrack, G4double previousStepSize, G4ForceCondition *condition)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
const G4ThreeVector & GetDecayDirection() const
std::mutex G4Mutex
Definition: G4Threading.hh:84