91 using namespace CLHEP;
96 #ifdef G4MULTITHREADED
101 G4int& G4RadioactiveDecayBase::NumberOfInstances()
103 static G4int numberOfInstances = 0;
104 return numberOfInstances;
110 forceDecayDirection(0.,0.,0.), forceDecayHalfAngle(0.*
deg), dirPath(
""),
115 G4cout <<
"G4RadioactiveDecayBase constructor: processName = " << processName
134 char* path_var = std::getenv(
"G4RADIOACTIVEDATA");
137 "Environment variable G4RADIOACTIVEDATA is not set");
140 std::ostringstream os;
142 std::ifstream testFile;
143 testFile.open(os.str() );
144 if (!testFile.is_open() )
146 "Environment variable G4RADIOACTIVEDATA is set, but does not point to correct directory");
153 #ifdef G4MULTITHREADED
154 G4AutoLock lk(&G4RadioactiveDecayBase::radioactiveDecayMutex);
155 NumberOfInstances()++;
172 outFile <<
"The radioactive decay process (G4RadioactiveDecay) handles the\n"
173 <<
"alpha, beta+, beta-, electron capture and isomeric transition\n"
174 <<
"decays of nuclei (G4GenericIon) with masses A > 4.\n"
175 <<
"The required half-lives and decay schemes are retrieved from\n"
176 <<
"the RadioactiveDecay database which was derived from ENSDF.\n";
184 for (DecayTableMap::iterator i =
dkmap->begin(); i !=
dkmap->end(); i++) {
189 #ifdef G4MULTITHREADED
190 G4AutoLock lk(&G4RadioactiveDecayBase::radioactiveDecayMutex);
191 --NumberOfInstances();
192 if(NumberOfInstances()==0)
194 for (DecayTableMap::iterator i = master_dkmap->begin(); i != master_dkmap->end(); i++) {
197 master_dkmap->clear();
207 if (((
const G4Ions*)(&aParticle))->GetExcitationEnergy() > 0.) {
return true;}
216 G4int A = ((
const G4Ions*) (&aParticle))->GetAtomicMass();
217 G4int Z = ((
const G4Ions*) (&aParticle))->GetAtomicNumber();
229 DecayTableMap::iterator table_ptr =
dkmap->find(key);
232 if (table_ptr ==
dkmap->end() ) {
234 if(theDecayTable) (*dkmap)[key] = theDecayTable;
236 theDecayTable = table_ptr->second;
238 return theDecayTable;
247 for (
size_t i = 0; i < theLogicalVolumes->size(); i++) {
248 volume = (*theLogicalVolumes)[i];
249 if (volume->
GetName() == aVolume) {
255 G4cout <<
" Radioactive decay applied to " << aVolume <<
G4endl;
257 }
else if (i == theLogicalVolumes->size() ) {
259 ed << aVolume <<
" is not a valid logical volume name. Decay not activated for it."
261 G4Exception(
"G4RadioactiveDecayBase::SelectAVolume()",
"HAD_RDM_300",
273 for (
size_t i = 0; i < theLogicalVolumes->size(); i++) {
274 volume = (*theLogicalVolumes)[i];
275 if (volume->
GetName() == aVolume) {
276 std::vector<G4String>::iterator location;
283 G4cout <<
" G4RadioactiveDecay::DeselectAVolume: " << aVolume
284 <<
" is removed from list " <<
G4endl;
287 ed << aVolume <<
" is not in the list. No action taken." <<
G4endl;
288 G4Exception(
"G4RadioactiveDecayBase::DeselectAVolume()",
"HAD_RDM_300",
291 }
else if (i == theLogicalVolumes->size()) {
293 ed << aVolume <<
" is not a valid logical volume name. No action taken."
295 G4Exception(
"G4RadioactiveDecayBase::DeselectAVolume()",
"HAD_RDM_300",
312 for (
size_t i = 0; i < theLogicalVolumes->size(); i++){
313 volume = (*theLogicalVolumes)[i];
351 G4cout <<
"G4RadioactiveDecay::GetMeanLifeTime() " <<
G4endl;
353 <<
" GeV, Mass: " << theParticle->
GetMass()/
GeV
354 <<
" GeV, Life time: " << theLife/
ns <<
" ns " <<
G4endl;
358 else if (theLife < 0.0) {meanlife =
DBL_MAX;}
359 else {meanlife = theLife;}
362 if (((
const G4Ions*)(theParticleDef))->GetExcitationEnergy() > 0. &&
363 meanlife ==
DBL_MAX) {meanlife = 0.;}
366 G4cout <<
" mean life time: " << meanlife/
s <<
" s " <<
G4endl;
388 G4cout <<
"G4RadioactiveDecay::GetMeanFreePath() " <<
G4endl;
390 <<
" GeV, Mass: " << aMass/
GeV <<
" GeV, tau: " << tau <<
" ns "
401 }
else if (tau < 0.0) {
404 ed <<
"Ion has negative lifetime " << tau
405 <<
" but is not stable. Setting mean free path to DBL_MAX" <<
G4endl;
406 G4Exception(
"G4RadioactiveDecay::GetMeanFreePath()",
"HAD_RDM_011",
413 pathlength =
c_light*tau*betaGamma;
419 G4cout <<
"G4Decay::GetMeanFreePath: "
421 <<
" stops, kinetic energy = "
431 G4cout <<
"mean free path: "<< pathlength/
m <<
" m" <<
G4endl;
469 os <<
"======================================================================"
471 os <<
"====== Radioactive Decay Physics Parameters ======="
473 os <<
"======================================================================"
475 os <<
"Max life time "
477 os <<
"Internal e- conversion flag "
479 os <<
"Stored internal conversion coefficients "
481 os <<
"Enable correlated gamma emission "
483 os <<
"Max 2J for sampling of angular correlations "
485 os <<
"Atomic de-excitation enabled "
486 << emparam->
Fluo() << endline;
487 os <<
"Auger electron emission enabled "
488 << emparam->
Auger() << endline;
489 os <<
"Auger cascade enabled "
491 os <<
"Check EM cuts disabled for atomic de-excitation "
493 os <<
"Use Bearden atomic level energies "
495 os <<
"======================================================================"
513 G4int A = ((
const G4Ions*)(&theParentNucleus))->GetAtomicMass();
514 G4int Z = ((
const G4Ions*)(&theParentNucleus))->GetAtomicNumber();
516 G4double levelEnergy = ((
const G4Ions*)(&theParentNucleus))->GetExcitationEnergy();
518 ((
const G4Ions*)(&theParentNucleus))->GetFloatLevelBase();
520 #ifdef G4MULTITHREADED
521 G4AutoLock lk(&G4RadioactiveDecayBase::radioactiveDecayMutex);
524 DecayTableMap::iterator master_table_ptr = master_dkmap->find(key);
526 if (master_table_ptr != master_dkmap->end() ) {
527 return master_table_ptr->second;
535 std::ostringstream os;
536 os <<
dirPath <<
"/z" << Z <<
".a" << A <<
'\0';
543 std::ifstream DecaySchemeFile;
544 DecaySchemeFile.open(file);
546 if (DecaySchemeFile.good()) {
549 const G4int nMode = 12;
550 G4double modeTotalBR[nMode] = {0.0};
552 for (
G4int i = 0; i < nMode; i++) {
556 char inputChars[120]={
' '};
577 while (!complete && !DecaySchemeFile.getline(inputChars, 120).eof()) {
580 G4Exception(
"G4RadioactiveDecay::LoadDecayTable()",
"HAD_RDM_100",
585 inputLine = inputChars;
586 inputLine = inputLine.
strip(1);
587 if (inputChars[0] !=
'#' && inputLine.length() != 0) {
588 std::istringstream tmpStream(inputLine);
590 if (inputChars[0] ==
'P') {
593 tmpStream >> recordType >> parentExcitation >> floatingFlag >> dummy;
602 if (floatingLevel !=
noFloat) {
605 if (!floatMatch) found =
false;
614 if (inputLine.length() < 72) {
615 tmpStream >> theDecayMode >> dummy >> decayModeTotal;
616 switch (theDecayMode) {
623 theDecayTable->
Insert(anITChannel);
628 modeTotalBR[1] = decayModeTotal;
break;
630 modeTotalBR[2] = decayModeTotal;
break;
632 modeTotalBR[3] = decayModeTotal;
break;
634 modeTotalBR[4] = decayModeTotal;
break;
636 modeTotalBR[5] = decayModeTotal;
break;
638 modeTotalBR[6] = decayModeTotal;
break;
640 modeTotalBR[7] = decayModeTotal;
break;
642 modeTotalBR[8] = decayModeTotal;
break;
644 modeTotalBR[9] = decayModeTotal;
break;
658 modeTotalBR[10] = decayModeTotal;
break;
660 modeTotalBR[11] = decayModeTotal;
break;
664 G4Exception(
"G4RadioactiveDecay::LoadDecayTable()",
"HAD_RDM_000",
669 if (inputLine.length() < 84) {
670 tmpStream >> theDecayMode >> a >> daughterFloatFlag >> b >>
c;
673 tmpStream >> theDecayMode >> a >> daughterFloatFlag >> b >> c >> betaType;
683 switch (theDecayMode) {
688 daughterFloatLevel, betaType);
691 theDecayTable->
Insert(aBetaMinusChannel);
700 daughterFloatLevel, betaType);
703 theDecayTable->
Insert(aBetaPlusChannel);
716 theDecayTable->
Insert(aKECChannel);
729 theDecayTable->
Insert(aLECChannel);
742 theDecayTable->
Insert(aMECChannel);
755 theDecayTable->
Insert(aNECChannel);
767 theDecayTable->
Insert(anAlphaChannel);
779 theDecayTable->
Insert(aProtonChannel);
791 theDecayTable->
Insert(aNeutronChannel);
826 theDecayTable->
Insert(aSpontFissChannel);
837 theDecayTable->
Insert(aTritonChannel);
845 G4Exception(
"G4RadioactiveDecay::LoadDecayTable()",
"HAD_RDM_000",
863 theNuclearDecayChannel =
static_cast<G4NuclearDecay*
>(theChannel);
866 if (theDecayMode !=
IT) {
867 theBR = theChannel->
GetBR();
868 theChannel->
SetBR(theBR*modeTotalBR[theDecayMode]/modeSumBR[theDecayMode]);
873 DecaySchemeFile.close();
875 if (!found && levelEnergy > 0) {
882 theDecayTable->
Insert(anITChannel);
889 #ifdef G4MULTITHREADED
892 return theDecayTable;
898 if (Z < 1 || A < 2)
G4cout <<
"Z and A not valid!" <<
G4endl;
900 std::ifstream DecaySchemeFile(filename);
901 if (DecaySchemeFile) {
902 G4int ID_ion = A*1000 +
Z;
906 ed << filename <<
" does not exist! " <<
G4endl;
907 G4Exception(
"G4RadioactiveDecayBase::AddUserDecayDataFile()",
"HAD_RDM_001",
933 G4cout <<
"G4RadioactiveDecay::DecayIt : "
935 <<
" is not selected for the RDM"<<
G4endl;
957 G4cout <<
"G4RadioactiveDecay::DecayIt : "
959 <<
" is not an ion or is outside (Z,A) limits set for the decay. "
960 <<
" Set particle change accordingly. "
975 if (theDecayTable == 0 || theDecayTable->
entries() == 0) {
980 G4cout <<
"G4RadioactiveDecay::DecayIt : "
981 <<
"decay table not defined for "
983 <<
". Set particle change accordingly. "
1094 if (products->
entries() == 1) {
1126 if (temptime < 0.) temptime = 0.;
1127 finalGlobalTime += temptime;
1128 finalLocalTime += temptime;
1131 products->
Boost(ParentEnergy, ParentDirection);
1138 G4cout <<
"G4RadioactiveDecay::DecayAnalog: Decay vertex :";
1139 G4cout <<
" Time: " << finalGlobalTime/
ns <<
"[ns]";
1144 G4cout <<
"G4Decay::DecayIt : decay products in Lab. Frame" <<
G4endl;
1149 for (
G4int index = 0; index < numberOfSecondaries; index++) {
1159 && index <numberOfSecondaries-1){
1192 if (theDecayChannel == 0) {
1196 G4Exception(
"G4RadioactiveDecay::DoDecay",
"HAD_RDM_013",
1202 G4cout <<
"G4RadioactiveDecay::DoIt : selected decay channel addr: "
1203 << theDecayChannel <<
G4endl;
1222 if (0 == products || 0 == products->
entries())
return;
1242 if (daughterType == electron || daughterType == positron ||
1243 daughterType == neutron || daughterType == gamma ||
1244 daughterType == alpha || daughterType == triton || daughterType == proton)
CollimateDecayProduct(daughter);
1251 G4cout <<
"CollimateDecayProduct for daughter "
1282 G4cout <<
" ChooseCollimationDirection returns " << dir <<
G4endl;