88 using namespace CLHEP;
93 #ifdef G4MULTITHREADED
101 forceDecayDirection(0.,0.,0.), forceDecayHalfAngle(0.*
deg), dirPath(
""),
106 G4cout <<
"G4RadioactiveDecayBase constructor: processName = " << processName
125 char* path_var = getenv(
"G4RADIOACTIVEDATA");
128 "Environment variable G4RADIOACTIVEDATA is not set");
131 std::ostringstream os;
133 std::ifstream testFile;
134 testFile.open(os.str() );
135 if (!testFile.is_open() )
137 "Environment variable G4RADIOACTIVEDATA is set, but does not point to correct directory");
144 #ifdef G4MULTITHREADED
145 G4AutoLock lk(&G4RadioactiveDecayBase::radioactiveDecayMutex);
162 outFile <<
"The radioactive decay process (G4RadioactiveDecay) handles the\n"
163 <<
"alpha, beta+, beta-, electron capture and isomeric transition\n"
164 <<
"decays of nuclei (G4GenericIon) with masses A > 4.\n"
165 <<
"The required half-lives and decay schemes are retrieved from\n"
166 <<
"the RadioactiveDecay database which was derived from ENSDF.\n";
174 for (DecayTableMap::iterator i =
dkmap->begin(); i !=
dkmap->end(); i++) {
185 if (((
const G4Ions*)(&aParticle))->GetExcitationEnergy() > 0.) {
return true;}
194 G4int A = ((
const G4Ions*) (&aParticle))->GetAtomicMass();
195 G4int Z = ((
const G4Ions*) (&aParticle))->GetAtomicNumber();
207 DecayTableMap::iterator table_ptr =
dkmap->find(key);
210 if (table_ptr ==
dkmap->end() ) {
212 if(theDecayTable) (*dkmap)[key] = theDecayTable;
214 theDecayTable = table_ptr->second;
217 return theDecayTable;
226 for (
size_t i = 0; i < theLogicalVolumes->size(); i++) {
227 volume=(*theLogicalVolumes)[i];
228 if (volume->
GetName() == aVolume) {
236 }
else if(i == theLogicalVolumes->size()) {
237 G4cerr <<
"SelectAVolume: "<< aVolume
238 <<
" is not a valid logical volume name" <<
G4endl;
249 for (
size_t i = 0; i < theLogicalVolumes->size(); i++){
250 volume=(*theLogicalVolumes)[i];
251 if (volume->
GetName() == aVolume) {
252 std::vector<G4String>::iterator location;
259 G4cerr <<
" DeselectVolume:" << aVolume <<
" is not in the list "
264 G4cout <<
" DeselectVolume: " << aVolume <<
" is removed from list "
267 }
else if (i == theLogicalVolumes->size()) {
268 G4cerr <<
" DeselectVolume:" << aVolume
269 <<
"is not a valid logical volume name" <<
G4endl;
285 for (
size_t i = 0; i < theLogicalVolumes->size(); i++){
286 volume = (*theLogicalVolumes)[i];
329 G4cout <<
"G4RadioactiveDecay::GetMeanLifeTime() " <<
G4endl;
331 <<
" GeV, Mass: " << theParticle->
GetMass()/
GeV
332 <<
" GeV, Life time: " << theLife/
ns <<
" ns " <<
G4endl;
336 else if (theLife < 0.0) {meanlife =
DBL_MAX;}
337 else {meanlife = theLife;}
340 if (((
const G4Ions*)(theParticleDef))->GetExcitationEnergy() > 0. &&
341 meanlife ==
DBL_MAX) {meanlife = 0.;}
345 G4cout <<
" mean life time: " << meanlife/
s <<
" s " <<
G4endl;
367 G4cout <<
"G4RadioactiveDecay::GetMeanFreePath() " <<
G4endl;
369 <<
" GeV, Mass: " << aMass/
GeV <<
" GeV, tau: " << tau <<
" ns "
380 }
else if (tau < 0.0) {
383 ed <<
"Ion has negative lifetime " << tau
384 <<
" but is not stable. Setting mean free path to DBL_MAX" <<
G4endl;
385 G4Exception(
"G4RadioactiveDecay::GetMeanFreePath()",
"HAD_RDM_011",
392 pathlength =
c_light*tau*betaGamma;
398 G4cout <<
"G4Decay::GetMeanFreePath: "
400 <<
" stops, kinetic energy = "
410 G4cout <<
"mean free path: "<< pathlength/
m <<
" m" <<
G4endl;
444 G4int A = ((
const G4Ions*)(&theParentNucleus))->GetAtomicMass();
445 G4int Z = ((
const G4Ions*)(&theParentNucleus))->GetAtomicNumber();
447 G4double levelEnergy = ((
const G4Ions*)(&theParentNucleus))->GetExcitationEnergy();
449 ((
const G4Ions*)(&theParentNucleus))->GetFloatLevelBase();
451 #ifdef G4MULTITHREADED
452 G4AutoLock lk(&G4RadioactiveDecayBase::radioactiveDecayMutex);
455 DecayTableMap::iterator master_table_ptr = master_dkmap->find(key);
457 if (master_table_ptr != master_dkmap->end() ) {
458 return master_table_ptr->second;
474 std::ostringstream os;
475 os <<
dirPath <<
"/z" << Z <<
".a" << A <<
'\0';
482 std::ifstream DecaySchemeFile;
483 DecaySchemeFile.open(file);
485 if (DecaySchemeFile.good()) {
488 const G4int nMode = 9;
489 G4double modeTotalBR[nMode] = {0.0};
491 for (
G4int i = 0; i < nMode; i++) {
495 char inputChars[120]={
' '};
517 ed <<
" While count exceeded " <<
G4endl;
519 while (!complete && !DecaySchemeFile.getline(inputChars, 120).eof()) {
526 inputLine = inputChars;
527 inputLine = inputLine.
strip(1);
528 if (inputChars[0] !=
'#' && inputLine.length() != 0) {
529 std::istringstream tmpStream(inputLine);
531 if (inputChars[0] ==
'P') {
534 tmpStream >> recordType >> parentExcitation >> floatingFlag >> dummy;
543 if (floatingLevel !=
noFloat) {
546 if (!floatMatch) found =
false;
555 if (inputLine.length() < 72) {
556 tmpStream >> theDecayMode >> dummy >> decayModeTotal;
557 switch (theDecayMode) {
564 theDecayTable->
Insert(anITChannel);
569 modeTotalBR[1] = decayModeTotal;
break;
571 modeTotalBR[2] = decayModeTotal;
break;
573 modeTotalBR[3] = decayModeTotal;
break;
575 modeTotalBR[4] = decayModeTotal;
break;
577 modeTotalBR[5] = decayModeTotal;
break;
579 modeTotalBR[6] = decayModeTotal;
break;
581 modeTotalBR[7] = decayModeTotal;
break;
583 modeTotalBR[8] = decayModeTotal;
break;
601 G4Exception(
"G4RadioactiveDecay::LoadDecayTable()",
"HAD_RDM_000",
606 if (inputLine.length() < 84) {
607 tmpStream >> theDecayMode >> a >> daughterFloatFlag >> b >> c;
610 tmpStream >> theDecayMode >> a >> daughterFloatFlag >> b >> c >> betaType;
619 switch (theDecayMode) {
624 daughterFloatLevel, betaType);
627 theDecayTable->
Insert(aBetaMinusChannel);
636 daughterFloatLevel, betaType);
639 theDecayTable->
Insert(aBetaPlusChannel);
652 theDecayTable->
Insert(aKECChannel);
665 theDecayTable->
Insert(aLECChannel);
678 theDecayTable->
Insert(aMECChannel);
690 theDecayTable->
Insert(anAlphaChannel);
702 theDecayTable->
Insert(aProtonChannel);
714 theDecayTable->
Insert(aNeutronChannel);
750 G4Exception(
"G4RadioactiveDecay::LoadDecayTable()",
"HAD_RDM_000",
768 theNuclearDecayChannel =
static_cast<G4NuclearDecay*
>(theChannel);
771 if (theDecayMode !=
IT) {
772 theBR = theChannel->
GetBR();
773 theChannel->
SetBR(theBR*modeTotalBR[theDecayMode]/modeSumBR[theDecayMode]);
778 DecaySchemeFile.close();
780 if (!found && levelEnergy > 0) {
787 theDecayTable->
Insert(anITChannel);
794 #ifdef G4MULTITHREADED
797 return theDecayTable;
803 if (Z < 1 || A < 2)
G4cout <<
"Z and A not valid!" <<
G4endl;
805 std::ifstream DecaySchemeFile(filename);
806 if (DecaySchemeFile) {
807 G4int ID_ion = A*1000 +
Z;
810 G4cout <<
"The file " << filename <<
" does not exist!" <<
G4endl;
836 G4cout <<
"G4RadioactiveDecay::DecayIt : "
838 <<
" is not selected for the RDM"<<
G4endl;
861 G4cerr <<
"G4RadioactiveDecay::DecayIt : "
863 <<
" is not a valid nucleus for the RDM"<<
G4endl;
876 if (theDecayTable == 0 || theDecayTable->
entries() == 0) {
881 G4cerr <<
"G4RadioactiveDecay::DecayIt : decay table not defined for ";
906 if (products->
entries() == 1) {
935 if (temptime < 0.) temptime = 0.;
936 finalGlobalTime += temptime;
937 finalLocalTime += temptime;
940 products->
Boost(ParentEnergy, ParentDirection);
948 G4cout <<
"G4RadioactiveDecay::DecayIt : Decay vertex :";
949 G4cout <<
" Time: " <<finalGlobalTime/
ns <<
"[ns]";
954 G4cout <<
"G4Decay::DecayIt : decay products in Lab. Frame" <<
G4endl;
959 for (index=0; index < numberOfSecondaries; index++) {
961 finalGlobalTime, currentPosition);
997 if (theDecayChannel == 0) {
1001 G4Exception(
"G4RadioactiveDecay::DoDecay",
"HAD_RDM_013",
1007 G4cerr <<
"G4RadioactiveDecay::DoIt : selected decay channel addr:";
1026 if (0 == products || 0 == products->
entries())
return;
1045 if (daughterType == electron || daughterType == positron ||
1046 daughterType == neutron || daughterType == gamma ||
1054 G4cout <<
"CollimateDecayProduct for daughter "
1085 G4cout <<
" ChooseCollimationDirection returns " << dir <<
G4endl;
virtual void ProcessDescription(std::ostream &outFile) const
static G4GenericIon * GenericIon()
G4double GetLocalTime() const
G4DecayTable * GetDecayTable(const G4ParticleDefinition *)
static G4Gamma * Definition()
static G4Proton * Definition()
void BuildPhysicsTable(const G4ParticleDefinition &)
virtual void SetICM(G4bool)
std::ostringstream G4ExceptionDescription
G4NucleusLimits theNucleusLimits
std::vector< ExP01TrackerHit * > a
G4LogicalVolume * GetLogicalVolume() const
static G4Alpha * Definition()
static constexpr double MeV
G4double GetPDGLifeTime() const
static G4Neutron * Definition()
void RegisterExtraProcess(G4VProcess *)
G4DeexPrecoParameters * GetParameters()
static constexpr double keV
void ProposeWeight(G4double finalWeight)
G4ThreeVector ChooseCollimationDirection() const
void AddSecondary(G4Track *aSecondary)
const G4ThreeVector & GetMomentumDirection() const
static G4NuclearLevelData * GetInstance()
G4double forceDecayHalfAngle
G4RadioactiveDecayBaseMessenger * theRadioactiveDecayBaseMessenger
const G4String & GetParticleName() const
const G4String & GetParticleType() const
const G4TouchableHandle & GetTouchableHandle() const
G4PhotonEvaporation * photonEvaporation
G4VParticleChange * DecayIt(const G4Track &theTrack, const G4Step &theStep)
void Insert(G4VDecayChannel *aChannel)
G4double GetMeanLifeTime(const G4Track &theTrack, G4ForceCondition *condition)
G4DynamicParticle * PopProducts()
G4double GetPDGMass() const
void Boost(G4double totalEnergy, const G4ThreeVector &momentumDirection)
G4int GetVerboseLevel() const
void SetMomentumDirection(const G4ThreeVector &aDirection)
G4double GetWeight() const
std::map< G4int, G4String > theUserRadioactiveDataFiles
void SelectAVolume(const G4String aVolume)
static constexpr double m
void SetBR(G4double value)
static const G4ThreeVector origin
#define G4MUTEX_INITIALIZER
G4ParticleDefinition * GetDefinition() const
G4VParticleChange * pParticleChange
static constexpr double deg
void CollimateDecay(G4DecayProducts *products)
void SetARM(G4bool onoff)
G4VDecayChannel * GetDecayChannel(G4int index) const
std::vector< G4String > ValidVolumes
G4DecayProducts * DoDecay(const G4ParticleDefinition &theParticleDef)
static const G4double alpha
G4double GetGlobalTime() const
G4String strip(G4int strip_Type=trailing, char c=' ')
double A(double temperature)
void SetTouchableHandle(const G4TouchableHandle &apValue)
const G4ThreeVector & GetPosition() const
static constexpr double eV
G4ParticleChangeForRadDecay fParticleChangeForRadDecay
void ProposeLocalTime(G4double t)
G4ThreeVector forceDecayDirection
G4GLOB_DLL std::ostream G4cerr
G4bool IsApplicable(const G4ParticleDefinition &)
G4bool GetPDGStable() const
static G4Electron * Definition()
void DeselectAVolume(const G4String aVolume)
static G4HadronicProcessStore * Instance()
void SetGoodForTrackingFlag(G4bool value=true)
static const G4double levelTolerance
void SetARM(G4bool onoff)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
void SetCorrelatedGamma(G4bool)
static constexpr double c_light
void DeselectAllVolumes()
void AddUserDecayDataFile(G4int Z, G4int A, G4String filename)
void SetProcessSubType(G4int)
G4double GetKineticEnergy() const
static constexpr double cm
G4GLOB_DLL std::ostream G4cout
void CollimateDecayProduct(G4DynamicParticle *product)
static constexpr double deg
void ClearNumberOfInteractionLengthLeft()
virtual G4DecayProducts * DecayIt(G4double parentMass=-1.0)=0
virtual void RDMForced(G4bool)
std::map< G4String, G4DecayTable * > DecayTableMap
G4VDecayChannel * SelectADecayChannel(G4double parentMass=-1.)
G4double GetTotalMomentum() const
static constexpr double pi
const G4ParticleDefinition * GetParticleDefinition() const
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
G4TrackStatus GetTrackStatus() const
G4DecayTable * LoadDecayTable(const G4ParticleDefinition &theParentNucleus)
static G4Positron * Definition()
static G4Ions::G4FloatLevelBase FloatLevelBase(char flbChar)
G4RadioactiveDecayMode GetDecayMode()
virtual void Initialize(const G4Track &)
G4RadioactiveDecayBase(const G4String &processName="RadioactiveDecayBase")
void ProposeTrackStatus(G4TrackStatus status)
static constexpr double GeV
const G4DynamicParticle * GetDynamicParticle() const
static G4LogicalVolumeStore * GetInstance()
void RegisterParticleForExtraProcess(G4VProcess *, const G4ParticleDefinition *)
G4double GetMeanFreePath(const G4Track &theTrack, G4double previousStepSize, G4ForceCondition *condition)
const G4String & GetName() const
G4VPhysicalVolume * GetVolume() const
~G4RadioactiveDecayBase()
void SetNumberOfSecondaries(G4int totSecondaries)