Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4RadioactiveDecay.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 //
26 #ifndef G4RadioactiveDecay_h
27 #define G4RadioactiveDecay_h 1
28 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
29 //
30 // MODULE: G4RadioactiveDecay.hh
31 //
32 // Version: 0.b.4
33 // Date: 14/04/00
34 // Author: F Lei & P R Truscott
35 // Organisation: DERA UK
36 // Customer: ESA/ESTEC, NOORDWIJK
37 // Contract: 12115/96/JG/NL Work Order No. 3
38 //
39 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
40 //
41 // CHANGE HISTORY
42 // --------------
43 // 17 October 2011, L Desorgher - Add the method AddUserDecayDataFile
44 //
45 // 01 June 2011, M. Kelsey -- Add directional biasing interface to allow for
46 // "collimation" of decay daughters.
47 //
48 // 29 February 2000, P R Truscott, DERA UK
49 // 0.b.3 release.
50 //
51 // 13 April 2000, F Lei, DERA UK
52 // 0.b.4 release. No change to this file
53 //
54 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
56 
57 #include <vector>
58 #include <map>
60 
61 #include "G4ios.hh"
62 #include "globals.hh"
65 
66 #include "G4NucleusLimits.hh"
69 #include "G4RadioactivityTable.hh"
70 #include "G4ThreeVector.hh"
71 #include "G4Threading.hh"
72 
73 class G4Fragment;
76 
77 typedef std::vector<G4RadioactiveDecayChainsFromParent> G4RadioactiveDecayParentChainTable;
78 typedef std::vector<G4RadioactiveDecayRatesToDaughter> G4RadioactiveDecayRates;
79 typedef std::map<G4String, G4DecayTable*> DecayTableMap;
80 
81 
83 {
84  // class description
85 
86  // Implementation of the radioactive decay process which simulates the
87  // decays of radioactive nuclei. These nuclei are submitted to RDM as
88  // G4Ions. The required half-lives and decay schemes are retrieved from
89  // the Radioactivity database which was derived from ENSDF.
90  // All decay products are submitted back to the particle tracking process
91  // through the G4ParticleChangeForRadDecay object.
92  // class description - end
93 
94  public: // with description
95 
96  G4RadioactiveDecay(const G4String& processName="RadioactiveDecay");
98 
99  virtual void ProcessDescription(std::ostream& outFile) const;
100 
101  // Return true if the specified isotope is
102  // 1) defined as "nucleus" and
103  // 2) it is within theNucleusLimit
105 
106  // Return decay table if it exists, if not, load it from file
108 
109  // Select a logical volume in which RDM applies
110  void SelectAVolume(const G4String aVolume);
111 
112  // Remove a logical volume from the RDM applied list
113  void DeselectAVolume(const G4String aVolume);
114 
115  // Select all logical volumes for the application of RDM
116  void SelectAllVolumes();
117 
118  // Remove all logical volumes from RDM applications
119  void DeselectAllVolumes();
120 
121  // Set the decay biasing scheme using the data in "filename"
122  void SetDecayBias(G4String filename);
123 
124  // Set the half-life threshold for isomer production
126 
127  // Enable/disable ICM
128  void SetICM(G4bool icm) {applyICM = icm;}
129 
130  // Enable/disable ARM
131  void SetARM(G4bool arm) {applyARM = arm;}
132 
133  // Set source exposure function using histograms in "filename"
134  void SetSourceTimeProfile(G4String filename);
135 
137  // Returns true if the coefficient and decay time table for all the
138  // descendants of the specified isotope are ready.
139  // used in VR decay mode only
140 
142  // Calculates the coefficient and decay time table for all the descendents
143  // of the specified isotope. Adds the calculated table to the private data
144  // member "theParentChainTable".
145  // used in VR decay mode only
146 
148  // Used to retrieve the coefficient and decay time table for all the
149  // descendants of the specified isotope from "theParentChainTable"
150  // and place it in "chainsFromParent".
151  // used in VR decay mode only
152 
153  void SetDecayRate(G4int,G4int,G4double, G4int, std::vector<G4double>,
154  std::vector<G4double>);
155  // Sets "theDecayRate" with data supplied in the arguements.
156  // used in VR decay mode only
157 
158  std::vector<G4RadioactivityTable*> GetTheRadioactivityTables()
159  {return theRadioactivityTables;}
160  // Return vector of G4Radioactivity map - should be used in VR mode only
161 
162  G4DecayTable* LoadDecayTable(const G4ParticleDefinition& theParentNucleus);
163  // Load the decay data of isotope theParentNucleus
164 
165  void AddUserDecayDataFile(G4int Z, G4int A,G4String filename);
166  // Allow the user to replace the radio-active decay data provided in Geant4
167  // by its own data file for a given isotope
168 
170  // Sets the VerboseLevel which controls duggering display
171 
172  inline G4int GetVerboseLevel() const {return verboseLevel;}
173  // Returns the VerboseLevel which controls level of debugging output
174 
175  inline void SetNucleusLimits(G4NucleusLimits theNucleusLimits1)
176  {theNucleusLimits = theNucleusLimits1 ;}
177  // Sets theNucleusLimits which specifies the range of isotopes
178  // the G4RadioactiveDecay applies.
179 
181  {return theNucleusLimits;}
182  // Returns theNucleusLimits which specifies the range of isotopes
183  // the G4RadioactiveDecay applies
184 
185  // Controls whether G4RadioactiveDecay runs in analogue mode or
186  // variance reduction mode. SetBRBias, SetSplitNuclei and
187  // SetSourceTimeProfile all turn off analogue mode and use VR mode
188  inline void SetAnalogueMonteCarlo (G4bool r ) {
189  AnalogueMC = r;
191  }
192 
193  // Controls whether G4RadioactiveDecay uses fast beta simulation mode
194  // Currently does nothing - kept for backward compatibility
195  inline void SetFBeta (G4bool r ) { FBeta = r; }
196 
197  // Returns true if the simulation is an analogue Monte Carlo, and false if
198  // any of the biassing schemes have been selected.
200 
201  // Sets whether branching ration bias scheme applies.
202  inline void SetBRBias(G4bool r) {
203  BRBias = r;
205  }
206 
207  // Sets the number of times a nucleus will decay when biased
208  inline void SetSplitNuclei(G4int r) {
209  NSplit = r;
211  }
212 
213  // Returns the nuclear splitting number
214  inline G4int GetSplitNuclei () {return NSplit;}
215 
216  inline void SetDecayDirection(const G4ThreeVector& theDir) {
217  forceDecayDirection = theDir.unit();
218  }
219 
220  inline const G4ThreeVector& GetDecayDirection() const {
221  return forceDecayDirection;
222  }
223 
224  inline void SetDecayHalfAngle(G4double halfAngle=0.*CLHEP::deg) {
226  }
227 
229 
230  // Force direction (random within half-angle) for "visible" daughters
231  // (applies to electrons, positrons, gammas, neutrons, protons or alphas)
232  inline void SetDecayCollimation(const G4ThreeVector& theDir,
233  G4double halfAngle = 0.*CLHEP::deg) {
234  SetDecayDirection(theDir);
235  SetDecayHalfAngle(halfAngle);
236  }
237 
239 
240  G4VParticleChange* DecayIt(const G4Track& theTrack,
241  const G4Step& theStep);
242 
243  protected:
244 
245  G4DecayProducts* DoDecay(const G4ParticleDefinition& theParticleDef);
246 
247  // Apply directional bias for "visible" daughters (e+-, gamma, n, p, alpha)
248  void CollimateDecay(G4DecayProducts* products);
251 
252  G4double GetMeanFreePath(const G4Track& theTrack, G4double previousStepSize,
254 
255  G4double GetMeanLifeTime(const G4Track& theTrack,
257 
259 
261 
262  G4int GetDecayTimeBin(const G4double aDecayTime);
263 
264  //Add gamma,Xray,conversion,and auger electrons for bias mode
267  G4double currenTime,
268  std::vector<double>& weights_v,
269  std::vector<double>& times_v,
270  std::vector<G4DynamicParticle*>& secondaries_v);
271 
272  private:
273 
274  void StreamInfo(std::ostream& os, const G4String& endline);
275 
278 
281 
283 
289 
293 
294  // Parameters for pre-collimated (biased) decay products
297  static const G4ThreeVector origin; // (0,0,0) for convenience
298 
305 
306  std::vector<G4String> ValidVolumes;
308 
313 
314  // for the radioactivity tables
315  std::vector<G4RadioactivityTable*> theRadioactivityTables;
317  static const G4double levelTolerance;
318 
319  // Radioactive decay database directory path
321 
322  // User-defined radioactive decay data files
323  std::map<G4int, G4String> theUserRadioactiveDataFiles;
324 
325  // Library of decay tables
327 #ifdef G4MULTITHREADED
328  static DecayTableMap* master_dkmap;
329 #endif
330 
331  // Remainder of life time at rest
334 
335 
336  // ParticleChange for decay process
338 
339  // inline implementations
340  inline
343  {
346  return fRemainderLifeTime;
347  }
348 
349  inline
351  const G4Step& theStep)
352  {return DecayIt(theTrack, theStep);}
353 
354  inline
356  const G4Step& theStep)
357  {return DecayIt(theTrack, theStep);}
358 
359 #ifdef G4MULTITHREADED
360  public:
361  static G4Mutex radioactiveDecayMutex;
362 #endif
363 };
364 
365 #endif
366 
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4ParticleChangeForRadDecay fParticleChangeForRadDecay
void SetDecayRate(G4int, G4int, G4double, G4int, std::vector< G4double >, std::vector< G4double >)
G4double GetMeanLifeTime(const G4Track &theTrack, G4ForceCondition *condition)
G4RadioactiveDecaymessenger * theRadioactiveDecaymessenger
static const G4ThreeVector origin
G4int GetVerboseLevel() const
G4bool IsRateTableReady(const G4ParticleDefinition &)
G4NucleusLimits theNucleusLimits
G4PhotonEvaporation * photonEvaporation
G4double ConvolveSourceTimeProfile(const G4double, const G4double)
void SetBRBias(G4bool r)
std::map< G4int, G4String > theUserRadioactiveDataFiles
G4double GetMeanFreePath(const G4Track &theTrack, G4double previousStepSize, G4ForceCondition *condition)
std::vector< G4RadioactivityTable * > theRadioactivityTables
std::vector< G4String > ValidVolumes
G4double condition(const G4ErrorSymMatrix &m)
void SelectAVolume(const G4String aVolume)
G4bool IsApplicable(const G4ParticleDefinition &)
void AddDeexcitationSpectrumForBiasMode(G4ParticleDefinition *apartDef, G4double weight, G4double currenTime, std::vector< double > &weights_v, std::vector< double > &times_v, std::vector< G4DynamicParticle * > &secondaries_v)
void StreamInfo(std::ostream &os, const G4String &endline)
G4int GetDecayTimeBin(const G4double aDecayTime)
void SetVerboseLevel(G4int value)
std::vector< G4RadioactiveDecayChainsFromParent > G4RadioactiveDecayParentChainTable
Float_t Z
G4double AtRestGetPhysicalInteractionLength(const G4Track &track, G4ForceCondition *condition)
void SetNucleusLimits(G4NucleusLimits theNucleusLimits1)
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
G4NucleusLimits GetNucleusLimits() const
std::vector< G4RadioactiveDecayRatesToDaughter > G4RadioactiveDecayRates
void GetChainsFromParent(const G4ParticleDefinition &)
void BuildPhysicsTable(const G4ParticleDefinition &)
G4RadioactiveDecayRates theDecayRateVector
const XML_Char int const XML_Char * value
Definition: expat.h:331
void SetHLThreshold(G4double hl)
void SetDecayBias(G4String filename)
G4RadioactiveDecay(const G4String &processName="RadioactiveDecay")
static constexpr double s
G4RadioactiveDecay & operator=(const G4RadioactiveDecay &right)
double A(double temperature)
const G4ThreeVector & GetDecayDirection() const
void SetDecayCollimation(const G4ThreeVector &theDir, G4double halfAngle=0.*CLHEP::deg)
G4VParticleChange * AtRestDoIt(const G4Track &theTrack, const G4Step &theStep)
Definition: G4Step.hh:76
G4VParticleChange * PostStepDoIt(const G4Track &theTrack, const G4Step &theStep)
G4RadioactiveDecayParentChainTable theParentChainTable
void DeselectAVolume(const G4String aVolume)
std::vector< G4RadioactivityTable * > GetTheRadioactivityTables()
void SetDecayDirection(const G4ThreeVector &theDir)
void SetARM(G4bool arm)
G4double GetDecayHalfAngle() const
Hep3Vector unit() const
void CalculateChainsFromParent(const G4ParticleDefinition &)
void SetICM(G4bool icm)
virtual void ProcessDescription(std::ostream &outFile) const
G4RadioactiveDecayChainsFromParent chainsFromParent
double weight
Definition: plottest35.C:25
G4RadioactiveDecayRatesToDaughter ratesToDaughter
int G4int
Definition: G4Types.hh:78
void CollimateDecayProduct(G4DynamicParticle *product)
G4DecayTable * GetDecayTable(const G4ParticleDefinition *)
G4ThreeVector ChooseCollimationDirection() const
void AddUserDecayDataFile(G4int Z, G4int A, G4String filename)
G4ForceCondition
static constexpr double deg
virtual G4double AtRestGetPhysicalInteractionLength(const G4Track &, G4ForceCondition *)
void SetSplitNuclei(G4int r)
std::map< G4String, G4DecayTable * > DecayTableMap
void SetAnalogueMonteCarlo(G4bool r)
void CollimateDecay(G4DecayProducts *products)
static const G4double levelTolerance
void SetSourceTimeProfile(G4String filename)
void SetDecayHalfAngle(G4double halfAngle=0.*CLHEP::deg)
G4DecayTable * LoadDecayTable(const G4ParticleDefinition &theParentNucleus)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4ThreeVector forceDecayDirection
G4VParticleChange * DecayIt(const G4Track &theTrack, const G4Step &theStep)
G4DecayProducts * DoDecay(const G4ParticleDefinition &theParticleDef)
std::mutex G4Mutex
Definition: G4Threading.hh:84