Skip to content

Commit

Permalink
radioactive decays to 11.2.2
Browse files Browse the repository at this point in the history
  • Loading branch information
martinunland committed Aug 29, 2024
1 parent b7d6075 commit 6445079
Show file tree
Hide file tree
Showing 5 changed files with 643 additions and 941 deletions.
4 changes: 2 additions & 2 deletions simulations/CMakeLists.txt
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
add_subdirectory(effective_area)
#add_subdirectory(radioactive_decays)
#add_subdirectory(effective_area)
add_subdirectory(radioactive_decays)
#add_subdirectory(supernova)
#add_subdirectory(efficiency_calibration)
172 changes: 48 additions & 124 deletions simulations/radioactive_decays/include/OMSimG4RadioactiveDecay.hh
Original file line number Diff line number Diff line change
Expand Up @@ -14,48 +14,10 @@
* @ingroup radioactive
*/



// ********************************************************************
// * License and Disclaimer *
// * *
// * The Geant4 software is copyright of the Copyright Holders of *
// * the Geant4 Collaboration. It is provided under the terms and *
// * conditions of the Geant4 Software License, included in the file *
// * LICENSE and available at http://cern.ch/geant4/license . These *
// * include a list of copyright holders. *
// * *
// * Neither the authors of this software system, nor their employing *
// * institutes,nor the agencies providing financial support for this *
// * work make any representation or warranty, express or implied, *
// * regarding this software system or assume any liability for its *
// * use. Please see the license in the file LICENSE and URL above *
// * for the full disclaimer and the limitation of liability. *
// * *
// * This code implementation is the result of the scientific and *
// * technical work of the GEANT4 collaboration. *
// * By using, copying, modifying or distributing the software (or *
// * any work based on the software) you agree to acknowledge its *
// * use in resulting scientific publications, and indicate your *
// * acceptance of all terms of the Geant4 Software license. *
// ********************************************************************
//
////////////////////////////////////////////////////////////////////////////////
// //
// File: G4RadioactiveDecay.hh //
// Author: D.H. Wright (SLAC) //
// Date: 9 August 2017 //
// Description: version the G4RadioactiveDecay process by F. Lei and //
// P.R. Truscott with biasing and activation calculations //
// removed to a derived class. It performs alpha, beta, //
// electron capture and isomeric transition decays of //
// radioactive nuclei. //
// //
////////////////////////////////////////////////////////////////////////////////

#ifndef G4RadioactiveDecay_h
#define G4RadioactiveDecay_h 1


#include <vector>
#include <map>
#include <CLHEP/Units/SystemOfUnits.h>
Expand All @@ -67,12 +29,14 @@

#include "G4NucleusLimits.hh"
#include "G4ThreeVector.hh"
#include "G4Threading.hh"
#include "G4RadioactiveDecayMode.hh"

class G4Fragment;
class G4RadioactiveDecayMessenger;
class G4PhotonEvaporation;
class G4Ions;
class G4DecayTable;
class G4ITDecay;

typedef std::map<G4String, G4DecayTable*> DecayTableMap;

Expand All @@ -91,16 +55,28 @@ class G4RadioactiveDecay : public G4VRestDiscreteProcess

public: // with description

G4RadioactiveDecay(const G4String& processName="RadioactiveDecay");
~G4RadioactiveDecay();

virtual void ProcessDescription(std::ostream& outFile) const;
G4RadioactiveDecay(const G4String& processName="RadioactiveDecay",
const G4double timeThreshold=-1.0);
~G4RadioactiveDecay() override;

G4bool IsApplicable(const G4ParticleDefinition&);
G4bool IsApplicable(const G4ParticleDefinition&) override;
// Return true if the specified isotope is
// 1) defined as "nucleus" and
// 2) it is within theNucleusLimit

G4VParticleChange* AtRestDoIt(const G4Track& theTrack,
const G4Step& theStep) override;

G4VParticleChange* PostStepDoIt(const G4Track& theTrack,
const G4Step& theStep) override;

void BuildPhysicsTable(const G4ParticleDefinition &) override;

void ProcessDescription(std::ostream& outFile) const override;

virtual G4VParticleChange* DecayIt(const G4Track& theTrack,
const G4Step& theStep);

// Return decay table if it exists, if not, load it from file
G4DecayTable* GetDecayTable(const G4ParticleDefinition*);

Expand All @@ -119,19 +95,13 @@ class G4RadioactiveDecay : public G4VRestDiscreteProcess
// Enable/disable ARM
void SetARM(G4bool arm) {applyARM = arm;}

G4DecayTable* LoadDecayTable(const G4ParticleDefinition& theParentNucleus);
G4DecayTable* LoadDecayTable(const G4Ions*);
// Load the decay data of isotope theParentNucleus

void AddUserDecayDataFile(G4int Z, G4int A,G4String filename);
void AddUserDecayDataFile(G4int Z, G4int A, const G4String& filename);
// Allow the user to replace the radio-active decay data provided in Geant4
// by its own data file for a given isotope

inline void SetVerboseLevel(G4int value) {verboseLevel = value;}
// Sets the VerboseLevel which controls duggering display

inline G4int GetVerboseLevel() const {return verboseLevel;}
// Returns the VerboseLevel which controls level of debugging output

inline void SetNucleusLimits(G4NucleusLimits theNucleusLimits1)
{theNucleusLimits = theNucleusLimits1 ;}
// Sets theNucleusLimits which specifies the range of isotopes
Expand Down Expand Up @@ -169,114 +139,68 @@ class G4RadioactiveDecay : public G4VRestDiscreteProcess
}
inline G4double GetThresholdForVeryLongDecayTime() const {return fThresholdForVeryLongDecayTime;}

void BuildPhysicsTable(const G4ParticleDefinition &);
void StreamInfo(std::ostream& os, const G4String& endline);

G4VParticleChange* DecayIt(const G4Track& theTrack,
const G4Step& theStep);
G4RadioactiveDecay(const G4RadioactiveDecay& right) = delete;
G4RadioactiveDecay& operator=(const G4RadioactiveDecay& right) = delete;

protected:
G4VParticleChange* TerminateDecay();
G4double GetMeanFreePath(const G4Track& theTrack, G4double previousStepSize,
G4ForceCondition* condition) override;

void DecayAnalog(const G4Track& theTrack);
G4double GetMeanLifeTime(const G4Track& theTrack,
G4ForceCondition* condition) override;

G4VParticleChange* TerminateDecay();
// sampling of products
void DecayAnalog(const G4Track& theTrack, G4DecayTable*);

G4DecayProducts* DoDecay(const G4ParticleDefinition& theParticleDef);
// sampling products at rest
G4DecayProducts* DoDecay(const G4ParticleDefinition&, G4DecayTable*);

// Apply directional bias for "visible" daughters (e+-, gamma, n, p, alpha)
void CollimateDecay(G4DecayProducts* products);
void CollimateDecayProduct(G4DynamicParticle* product);
G4ThreeVector ChooseCollimationDirection() const;

G4double GetMeanFreePath(const G4Track& theTrack, G4double previousStepSize,
G4ForceCondition* condition);

G4double GetMeanLifeTime(const G4Track& theTrack,
G4ForceCondition* condition);

// ParticleChange for decay process
G4ParticleChangeForRadDecay fParticleChangeForRadDecay;

G4RadioactiveDecayMessenger* theRadioactiveDecayMessenger;
G4PhotonEvaporation* photonEvaporation;
G4ITDecay* decayIT;

std::vector<G4String> ValidVolumes;
bool isAllVolumesMode;
bool isAllVolumesMode{true};

static const G4double levelTolerance;

// Library of decay tables
DecayTableMap* dkmap;
#ifdef G4MULTITHREADED
static DecayTableMap* master_dkmap;
#endif

private:

void StreamInfo(std::ostream& os, const G4String& endline);

G4RadioactiveDecay(const G4RadioactiveDecay& right);
G4RadioactiveDecay& operator=(const G4RadioactiveDecay& right);

G4NucleusLimits theNucleusLimits;

G4bool isInitialised;

G4bool applyARM;
G4bool isInitialised{false};
G4bool applyARM{true};

// Parameters for pre-collimated (biased) decay products
G4ThreeVector forceDecayDirection;
G4double forceDecayHalfAngle;
G4ThreeVector forceDecayDirection{G4ThreeVector(0., 0., 0.)};
G4double forceDecayHalfAngle{0.0};
static const G4ThreeVector origin; // (0,0,0) for convenience

// Radioactive decay database directory path
G4String dirPath;

//User define radioactive decay data files replacing some files in the G4RADECAY database
std::map<G4int, G4String> theUserRadioactiveDataFiles;

//The last RadDecayMode
G4RadioactiveDecayMode theRadDecayMode;

// // Library of decay tables
// DecayTableMap* dkmap;
// #ifdef G4MULTITHREADED
// static DecayTableMap* master_dkmap;
// #endif
static G4String dirPath;

// Remainder of life time at rest
G4double fRemainderLifeTime;
G4int verboseLevel;
// User define radioactive decay data files replacing some files in the G4RADECAY database
static std::map<G4int, G4String>* theUserRDataFiles;

// The last RadDecayMode
G4RadioactiveDecayMode theRadDecayMode{G4RadioactiveDecayMode::IT};

// Ignore radioactive decays at rest of nuclides happening after this (very long) time threshold
G4double fThresholdForVeryLongDecayTime;

// inline implementations
inline
G4double AtRestGetPhysicalInteractionLength(const G4Track& track,
G4ForceCondition* condition)
{
fRemainderLifeTime =
G4VRestDiscreteProcess::AtRestGetPhysicalInteractionLength(track, condition);
return fRemainderLifeTime;
}

inline
G4VParticleChange* AtRestDoIt(const G4Track& theTrack,
const G4Step& theStep)
{return DecayIt(theTrack, theStep);}

inline
G4VParticleChange* PostStepDoIt(const G4Track& theTrack,
const G4Step& theStep)
{return DecayIt(theTrack, theStep);}

#ifdef G4MULTITHREADED
public:
static G4Mutex radioactiveDecayMutex;
protected:
G4int& NumberOfInstances();
#endif
};

#endif

#endif
3 changes: 2 additions & 1 deletion simulations/radioactive_decays/src/OMSimDecaysGPS.cc
Original file line number Diff line number Diff line change
Expand Up @@ -4,8 +4,9 @@
#include <G4SystemOfUnits.hh>
#include <G4Poisson.hh>

void OMSimDecaysGPS::setProductionRadius(G4double p_productionRadius) : m_productionRadius(p_productionRadius)
void OMSimDecaysGPS::setProductionRadius(G4double p_productionRadius)
{
m_productionRadius = p_productionRadius;
}

G4String OMSimDecaysGPS::getDecayTerminationNuclide()
Expand Down
18 changes: 6 additions & 12 deletions simulations/radioactive_decays/src/OMSimEventAction.cc
Original file line number Diff line number Diff line change
Expand Up @@ -5,31 +5,25 @@

#include <G4RunManager.hh>

/**
* @brief Custom actions at the beginning of the event.
* This function sets the current event ID in the EventInfoManager.
* @param evt Pointer to the current event.
*/
void OMSimEventAction::BeginOfEventAction(const G4Event *evt)
{
EventInfoManager::getInstance().setCurrentEventID(evt->GetEventID());
}

/**
* @brief Custom actions at the end of the event.
*
* Depending on the arguments set, this function will write hit and decay
* information with the analysis manager and reset hit and analysis data for the next event.
* @param evt Pointer to the current event.
* @param p_event Pointer to the current event.
*/
void OMSimEventAction::EndOfEventAction(const G4Event *evt)
void OMSimEventAction::EndOfEventAction(const G4Event *p_event)
{
if (!OMSimCommandArgsTable::getInstance().get<bool>("multiplicity_study"))
{
log_debug("End of event, saving information and reseting (thread {})", G4Threading::G4GetThreadId());
OMSimDecaysAnalysis &lAnalysisManager = OMSimDecaysAnalysis::getInstance();
lAnalysisManager.writeThreadHitInformation();
lAnalysisManager.writeThreadDecayInformation();
lAnalysisManager.reset();
OMSimDecaysAnalysis &analysisManager = OMSimDecaysAnalysis::getInstance();
analysisManager.writeThreadHitInformation();
analysisManager.writeThreadDecayInformation();
analysisManager.reset();
}
}
Loading

0 comments on commit 6445079

Please sign in to comment.