Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4RadioactiveDecay Class Reference

#include <G4RadioactiveDecay.hh>

Inheritance diagram for G4RadioactiveDecay:
Collaboration diagram for G4RadioactiveDecay:

Public Member Functions

 G4RadioactiveDecay (const G4String &processName="RadioactiveDecay")
 
 ~G4RadioactiveDecay ()
 
G4bool IsApplicable (const G4ParticleDefinition &)
 
G4DecayTableGetDecayTable (const G4ParticleDefinition *)
 
void SelectAVolume (const G4String aVolume)
 
void DeselectAVolume (const G4String aVolume)
 
void SelectAllVolumes ()
 
void DeselectAllVolumes ()
 
void SetDecayBias (G4String filename)
 
void SetHLThreshold (G4double hl)
 
void SetICM (G4bool icm)
 
void SetARM (G4bool arm)
 
void SetSourceTimeProfile (G4String filename)
 
G4bool IsRateTableReady (const G4ParticleDefinition &)
 
void AddDecayRateTable (const G4ParticleDefinition &)
 
void GetDecayRateTable (const G4ParticleDefinition &)
 
void SetDecayRate (G4int, G4int, G4double, G4int, std::vector< G4double >, std::vector< G4double >)
 
std::vector
< G4RadioactivityTable * > 
GetTheRadioactivityTables ()
 
G4DecayTableLoadDecayTable (const G4ParticleDefinition &theParentNucleus)
 
void AddUserDecayDataFile (G4int Z, G4int A, G4String filename)
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
void SetNucleusLimits (G4NucleusLimits theNucleusLimits1)
 
G4NucleusLimits GetNucleusLimits () const
 
void SetAnalogueMonteCarlo (G4bool r)
 
void SetFBeta (G4bool r)
 
G4bool IsAnalogueMonteCarlo ()
 
void SetBRBias (G4bool r)
 
void SetSplitNuclei (G4int r)
 
G4int GetSplitNuclei ()
 
void SetDecayDirection (const G4ThreeVector &theDir)
 
const G4ThreeVectorGetDecayDirection () const
 
void SetDecayHalfAngle (G4double halfAngle=0.*CLHEP::deg)
 
G4double GetDecayHalfAngle () const
 
void SetDecayCollimation (const G4ThreeVector &theDir, G4double halfAngle=0.*CLHEP::deg)
 
void BuildPhysicsTable (const G4ParticleDefinition &)
 
G4VParticleChangeDecayIt (const G4Track &theTrack, const G4Step &theStep)
 
- Public Member Functions inherited from G4VRestDiscreteProcess
 G4VRestDiscreteProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VRestDiscreteProcess (G4VRestDiscreteProcess &)
 
virtual ~G4VRestDiscreteProcess ()
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4int operator== (const G4VProcess &right) const
 
G4int operator!= (const G4VProcess &right) const
 
G4double GetCurrentInteractionLength () const
 
void SetPILfactor (G4double value)
 
G4double GetPILfactor () const
 
G4double AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)
 
G4double AtRestGPIL (const G4Track &track, G4ForceCondition *condition)
 
G4double PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
const G4StringGetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
 
const G4StringGetProcessName () const
 
G4ProcessType GetProcessType () const
 
void SetProcessType (G4ProcessType)
 
G4int GetProcessSubType () const
 
void SetProcessSubType (G4int)
 
virtual void StartTracking (G4Track *)
 
virtual void EndTracking ()
 
virtual void SetProcessManager (const G4ProcessManager *)
 
virtual const G4ProcessManagerGetProcessManager ()
 
virtual void ResetNumberOfInteractionLengthLeft ()
 
G4double GetNumberOfInteractionLengthLeft () const
 
G4double GetTotalNumberOfInteractionLengthTraversed () const
 
G4bool isAtRestDoItIsEnabled () const
 
G4bool isAlongStepDoItIsEnabled () const
 
G4bool isPostStepDoItIsEnabled () const
 
virtual void DumpInfo () const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
virtual void SetMasterProcess (G4VProcess *masterP)
 
const G4VProcessGetMasterProcess () const
 
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
 
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)
 

Protected Member Functions

G4DecayProductsDoDecay (const G4ParticleDefinition &theParticleDef)
 
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)
 
G4double ConvolveSourceTimeProfile (const G4double, const G4double)
 
G4double GetDecayTime ()
 
G4int GetDecayTimeBin (const G4double aDecayTime)
 
void AddDeexcitationSpectrumForBiasMode (G4ParticleDefinition *apartDef, G4double weight, G4double currenTime, std::vector< double > &weights_v, std::vector< double > &times_v, std::vector< G4DynamicParticle * > &secondaries_v)
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChangepParticleChange
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft
 
G4double currentInteractionLength
 
G4double theInitialNumberOfInteractionLength
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType
 
G4int theProcessSubType
 
G4double thePILfactor
 
G4bool enableAtRestDoIt
 
G4bool enableAlongStepDoIt
 
G4bool enablePostStepDoIt
 
G4int verboseLevel
 

Detailed Description

Definition at line 81 of file G4RadioactiveDecay.hh.

Constructor & Destructor Documentation

G4RadioactiveDecay::G4RadioactiveDecay ( const G4String processName = "RadioactiveDecay")

Definition at line 162 of file G4RadioactiveDecay.cc.

163  : G4VRestDiscreteProcess(processName, fDecay), isInitialised(false),
164  forceDecayDirection(0.,0.,0.), forceDecayHalfAngle(0.*deg), verboseLevel(0)
165 {
166 #ifdef G4VERBOSE
167  if (GetVerboseLevel() > 1) {
168  G4cout << "G4RadioactiveDecay constructor: processName = " << processName
169  << G4endl;
170  }
171 #endif
172 
174 
175  theRadioactiveDecaymessenger = new G4RadioactiveDecaymessenger(this);
176  pParticleChange = &fParticleChangeForRadDecay;
177 
178  // Reset the list of user defined data files
179  theUserRadioactiveDataFiles.clear();
180 
181  // Instantiate the map of decay tables
182 #ifdef G4MULTITHREADED
183  G4AutoLock lk(&G4RadioactiveDecay::radioactiveDecayMutex);
184  if(!master_dkmap) master_dkmap = new DecayTableMap;
185 #endif
186  dkmap = new DecayTableMap;
187 
188  // Apply default values.
189  NSourceBin = 1;
190  SBin[0] = 0.* s;
191  SBin[1] = 1.* s;
192  SProfile[0] = 1.;
193  SProfile[1] = 0.;
194  NDecayBin = 1;
195  DBin[0] = 0. * s ;
196  DBin[1] = 1. * s;
197  DProfile[0] = 1.;
198  DProfile[1] = 0.;
199  decayWindows[0] = 0;
201  theRadioactivityTables.push_back(rTable);
202  NSplit = 1;
203  AnalogueMC = true ;
204  FBeta = false ;
205  BRBias = true ;
206  applyICM = true ;
207  applyARM = true ;
208  halflifethreshold = nanosecond;
209 
210  // RDM applies to all logical volumes by default
211  isAllVolumesMode = true;
213 }
std::map< G4String, G4DecayTable * > DecayTableMap
G4int GetVerboseLevel() const
static constexpr double nanosecond
Definition: G4SIunits.hh:158
const XML_Char * s
Definition: expat.h:262
G4GLOB_DLL std::ostream G4cout
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:432
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
#define G4endl
Definition: G4ios.hh:61
static constexpr double deg
Definition: G4SIunits.hh:152

Here is the call graph for this function:

G4RadioactiveDecay::~G4RadioactiveDecay ( )

Definition at line 216 of file G4RadioactiveDecay.cc.

217 {
218  delete theRadioactiveDecaymessenger;
219  for (DecayTableMap::iterator i = dkmap->begin(); i != dkmap->end(); i++) {
220  delete i->second;
221  }
222  dkmap->clear();
223  delete dkmap;
224 }

Member Function Documentation

void G4RadioactiveDecay::AddDecayRateTable ( const G4ParticleDefinition theParentNucleus)

Definition at line 1146 of file G4RadioactiveDecay.cc.

1147 {
1148  // Use extended Bateman equation to calculate the radioactivities of all
1149  // progeny of theParentNucleus. The coefficients required to do this are
1150  // calculated using the method of P. Truscott (Ph.D. thesis and
1151  // DERA Technical Note DERA/CIS/CIS2/7/36/4/10) 11 January 2000.
1152  // Coefficients are then added to the decay rate table vector
1153 
1154  // Create and initialise variables used in the method.
1155  theDecayRateVector.clear();
1156 
1157  G4int nGeneration = 0;
1158 
1159  std::vector<G4double> taos;
1160 
1161  // Dimensionless A coefficients of Eqs. 4.24 and 4.25 of the TN
1162  std::vector<G4double> Acoeffs;
1163 
1164  // According to Eq. 4.26 the first coefficient (A_1:1) is -1
1165  Acoeffs.push_back(-1.);
1166 
1167  G4int A = ((const G4Ions*)(&theParentNucleus))->GetAtomicMass();
1168  G4int Z = ((const G4Ions*)(&theParentNucleus))->GetAtomicNumber();
1169  G4double E = ((const G4Ions*)(&theParentNucleus))->GetExcitationEnergy();
1170  G4double tao = theParentNucleus.GetPDGLifeTime();
1171  if (tao < 0.) tao = 1e-100;
1172  taos.push_back(tao);
1173  G4int nEntry = 0;
1174 
1175  // Fill the decay rate container (G4RadioactiveDecayRate) with the parent
1176  // isotope data
1177  SetDecayRate(Z,A,E,nGeneration,Acoeffs,taos); // Fill TP with parent lifetime
1178 
1179  // store the decay rate in decay rate vector
1180  theDecayRateVector.push_back(theDecayRate);
1181  nEntry++;
1182 
1183  // Now start treating the secondary generations.
1184  G4bool stable = false;
1185  G4int i;
1186  G4int j;
1187  G4VDecayChannel* theChannel = 0;
1188  G4NuclearDecay* theNuclearDecayChannel = 0;
1189 
1190  G4ITDecay* theITChannel = 0;
1191  G4BetaMinusDecay* theBetaMinusChannel = 0;
1192  G4BetaPlusDecay* theBetaPlusChannel = 0;
1193  G4AlphaDecay* theAlphaChannel = 0;
1194  G4ProtonDecay* theProtonChannel = 0;
1195  G4NeutronDecay* theNeutronChannel = 0;
1196  G4RadioactiveDecayMode theDecayMode;
1197  G4double theBR = 0.0;
1198  G4int AP = 0;
1199  G4int ZP = 0;
1200  G4int AD = 0;
1201  G4int ZD = 0;
1202  G4double EP = 0.;
1203  std::vector<G4double> TP;
1204  std::vector<G4double> RP; // A coefficients of the previous generation
1205  G4ParticleDefinition *theDaughterNucleus;
1206  G4double daughterExcitation;
1207  G4double nearestEnergy = 0.0;
1208  G4int nearestLevelIndex = 0;
1209  G4ParticleDefinition *aParentNucleus;
1210  G4IonTable* theIonTable;
1211  G4DecayTable* parentDecayTable;
1212  G4double theRate;
1213  G4double TaoPlus;
1214  G4int nS = 0; // Running index of first decay in a given generation
1215  G4int nT = nEntry; // Total number of decays accumulated over entire history
1216  const G4int nMode = 9;
1217  G4double brs[nMode];
1218  //
1219  theIonTable =
1221 
1222  G4int loop = 0;
1224  ed << " While count exceeded " << G4endl;
1225 
1226  while (!stable) { /* Loop checking, 01.09.2015, D.Wright */
1227  loop++;
1228  if (loop > 10000) {
1229  G4Exception("G4RadioactiveDecay::AddDecayRateTable()", "HAD_RDM_100", JustWarning, ed);
1230  break;
1231  }
1232  nGeneration++;
1233  for (j = nS; j < nT; j++) {
1234  // First time through, get data for parent nuclide
1235  ZP = theDecayRateVector[j].GetZ();
1236  AP = theDecayRateVector[j].GetA();
1237  EP = theDecayRateVector[j].GetE();
1238  RP = theDecayRateVector[j].GetDecayRateC();
1239  TP = theDecayRateVector[j].GetTaos();
1240  if (GetVerboseLevel() > 0) {
1241  G4cout << "G4RadioactiveDecay::AddDecayRateTable : daughters of ("
1242  << ZP << ", " << AP << ", " << EP
1243  << ") are being calculated, generation = " << nGeneration
1244  << G4endl;
1245  }
1246 // G4cout << " Taus = " << G4endl;
1247 // for (G4int ii = 0; ii < TP.size(); ii++) G4cout << TP[ii] << ", " ;
1248 // G4cout << G4endl;
1249 
1250  aParentNucleus = theIonTable->GetIon(ZP,AP,EP);
1251  parentDecayTable = GetDecayTable(aParentNucleus);
1252 
1253  G4DecayTable* summedDecayTable = new G4DecayTable();
1254  // This instance of G4DecayTable is for accumulating BRs and decay
1255  // channels. It will contain one decay channel per type of decay
1256  // (alpha, beta, etc.); its branching ratio will be the sum of all
1257  // branching ratios for that type of decay of the parent. If the
1258  // halflife of a particular channel is longer than some threshold,
1259  // that channel will be inserted specifically and its branching
1260  // ratio will not be included in the above sums.
1261  // This instance is not used to perform actual decays.
1262 
1263  for (G4int k = 0; k < nMode; k++) brs[k] = 0.0;
1264 
1265  // Go through the decay table and sum all channels having the same decay mode
1266  for (i = 0; i < parentDecayTable->entries(); i++) {
1267  theChannel = parentDecayTable->GetDecayChannel(i);
1268  theNuclearDecayChannel = static_cast<G4NuclearDecay*>(theChannel);
1269  theDecayMode = theNuclearDecayChannel->GetDecayMode();
1270  daughterExcitation = theNuclearDecayChannel->GetDaughterExcitation();
1271  theDaughterNucleus = theNuclearDecayChannel->GetDaughterNucleus() ;
1272 
1273  AD = ((const G4Ions*)(theDaughterNucleus))->GetAtomicMass();
1274  ZD = ((const G4Ions*)(theDaughterNucleus))->GetAtomicNumber();
1275  const G4LevelManager* levelManager =
1277 
1278  if (levelManager->NumberOfTransitions() ) {
1279  nearestEnergy = levelManager->NearestLevelEnergy(daughterExcitation);
1280  if (std::abs(daughterExcitation - nearestEnergy) < levelTolerance) {
1281  // Level half-life is in ns and the threshold is set to 1 micros
1282  // by default, user can set it via the UI command
1283  nearestLevelIndex = levelManager->NearestLevelIndex(daughterExcitation);
1284  if (levelManager->LifeTime(nearestLevelIndex)*ns >= halflifethreshold){
1285  // save the metastable nucleus
1286  summedDecayTable->Insert(theChannel);
1287  } else {
1288  brs[theDecayMode] += theChannel->GetBR();
1289  }
1290  } else {
1291  brs[theDecayMode] += theChannel->GetBR();
1292  }
1293  } else {
1294  brs[theDecayMode] += theChannel->GetBR();
1295  }
1296  } // Combine decay channels (loop i)
1297 
1298  brs[2] = brs[2]+brs[3]+brs[4]+brs[5]; // Combine beta+ and EC
1299  brs[3] = brs[4] =brs[5] = 0.0;
1300  for (i= 0; i<nMode; i++){ // loop over decay modes
1301  if (brs[i] > 0.) {
1302  switch ( i ) {
1303  case 0:
1304  // Decay mode is isomeric transition
1305  theITChannel = new G4ITDecay(aParentNucleus, brs[0], 0.0, 0.0);
1306 
1307  summedDecayTable->Insert(theITChannel);
1308  break;
1309 
1310  case 1:
1311  // Decay mode is beta-
1312  theBetaMinusChannel = new G4BetaMinusDecay(aParentNucleus, brs[1],
1313  0.*MeV, 0.*MeV,
1314  noFloat, allowed);
1315  summedDecayTable->Insert(theBetaMinusChannel);
1316  break;
1317 
1318  case 2:
1319  // Decay mode is beta+ + EC.
1320  theBetaPlusChannel = new G4BetaPlusDecay(aParentNucleus, brs[2], // DHW: April 2015
1321  0.*MeV, 0.*MeV,
1322  noFloat, allowed);
1323  summedDecayTable->Insert(theBetaPlusChannel);
1324  break;
1325 
1326  case 6:
1327  // Decay mode is alpha.
1328  theAlphaChannel = new G4AlphaDecay(aParentNucleus, brs[6], 0.*MeV,
1329  0.*MeV, noFloat);
1330  summedDecayTable->Insert(theAlphaChannel);
1331  break;
1332 
1333  case 7:
1334  // Decay mode is proton.
1335  theProtonChannel = new G4ProtonDecay(aParentNucleus, brs[7], 0.*MeV,
1336  0.*MeV, noFloat);
1337  summedDecayTable->Insert(theProtonChannel);
1338  break;
1339  case 8:
1340  // Decay mode is neutron.
1341  theNeutronChannel = new G4NeutronDecay(aParentNucleus, brs[8], 0.*MeV,
1342  0.*MeV, noFloat);
1343  summedDecayTable->Insert(theNeutronChannel);
1344  break;
1345 
1346  default:
1347  break;
1348  }
1349  }
1350  }
1351  // loop over all branches in summedDecayTable
1352  //
1353  for (i = 0; i < summedDecayTable->entries(); i++){
1354  theChannel = summedDecayTable->GetDecayChannel(i);
1355  theNuclearDecayChannel = static_cast<G4NuclearDecay*>(theChannel);
1356  theBR = theChannel->GetBR();
1357  theDaughterNucleus = theNuclearDecayChannel->GetDaughterNucleus();
1358 
1359  // First check if the decay of the original nucleus is an IT channel,
1360  // if true create a new ground-state nucleus
1361  if (theNuclearDecayChannel->GetDecayMode() == IT && nGeneration == 1) {
1362 
1363  A = ((const G4Ions*)(theDaughterNucleus))->GetAtomicMass();
1364  Z = ((const G4Ions*)(theDaughterNucleus))->GetAtomicNumber();
1365  theDaughterNucleus=theIonTable->GetIon(Z,A,0.);
1366  }
1367  if (IsApplicable(*theDaughterNucleus) && theBR &&
1368  aParentNucleus != theDaughterNucleus) {
1369  // need to make sure daughter has decay table
1370  parentDecayTable = GetDecayTable(theDaughterNucleus);
1371 
1372  if (parentDecayTable->entries() ) {
1373  A = ((const G4Ions*)(theDaughterNucleus))->GetAtomicMass();
1374  Z = ((const G4Ions*)(theDaughterNucleus))->GetAtomicNumber();
1375  E = ((const G4Ions*)(theDaughterNucleus))->GetExcitationEnergy();
1376 
1377  TaoPlus = theDaughterNucleus->GetPDGLifeTime();
1378  if (TaoPlus <= 0.) TaoPlus = 1e-100;
1379 
1380  // first set the taos, one simply need to add to the parent ones
1381  taos.clear();
1382  taos = TP; // load lifetimes of all previous generations
1383  size_t k;
1384  //check that TaoPlus differs from other taos from at least 1.e5 relative difference
1385  //for (k = 0; k < TP.size(); k++){
1386  //if (std::abs((TaoPlus-TP[k])/TP[k])<1.e-5 ) TaoPlus=1.00001*TP[k];
1387  //}
1388  taos.push_back(TaoPlus); // add daughter lifetime to list
1389  // now calculate the coefficiencies
1390  //
1391  // they are in two parts, first the less than n ones
1392  // Eq 4.24 of the TN
1393  Acoeffs.clear();
1394  long double ta1,ta2;
1395  ta2 = (long double)TaoPlus;
1396  for (k = 0; k < RP.size(); k++){
1397  ta1 = (long double)TP[k]; // loop over lifetimes of all previous generations
1398  if (ta1 == ta2) {
1399  theRate = 1.e100;
1400  } else {
1401  theRate = ta1/(ta1-ta2);
1402  }
1403  theRate = theRate * theBR * RP[k];
1404  Acoeffs.push_back(theRate);
1405  }
1406 
1407  // the second part: the n:n coefficiency
1408  // Eq 4.25 of the TN. Note Yn+1 is zero apart from Y1 which is -1
1409  // as treated at line 1013
1410  theRate = 0.;
1411  long double aRate, aRate1;
1412  aRate1 = 0.L;
1413  for (k = 0; k < RP.size(); k++){
1414  ta1 = (long double)TP[k];
1415  if (ta1 == ta2 ) {
1416  aRate = 1.e100;
1417  } else {
1418  aRate = ta2/(ta1-ta2);
1419  }
1420  aRate = aRate * (long double)(theBR * RP[k]);
1421  aRate1 += aRate;
1422  }
1423  theRate = -aRate1;
1424  Acoeffs.push_back(theRate);
1425  SetDecayRate (Z,A,E,nGeneration,Acoeffs,taos);
1426  theDecayRateVector.push_back(theDecayRate);
1427  nEntry++;
1428  } // there are entries in the table
1429  } // nuclide is OK to decay
1430  } // end of loop (i) over decay table branches
1431  // delete summedDecayTable;
1432 
1433  } // Getting contents of decay rate vector (end loop on j)
1434  nS = nT;
1435  nT = nEntry;
1436  if (nS == nT) stable = true;
1437  } // while nuclide is not stable
1438 
1439  // end of while loop
1440  // the calculation completed here
1441 
1442 
1443  // fill the first part of the decay rate table
1444  // which is the name of the original particle (isotope)
1445  theDecayRateTable.SetIonName(theParentNucleus.GetParticleName());
1446 
1447  // now fill the decay table with the newly completed decay rate vector
1448  theDecayRateTable.SetItsRates(theDecayRateVector);
1449 
1450  // finally add the decayratetable to the tablevector
1451  theDecayRateTableVector.push_back(theDecayRateTable);
1452 }
G4double GetBR() const
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
const G4LevelManager * GetLevelManager(G4int Z, G4int A)
G4float LifeTime(size_t i) const
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:503
G4RadioactiveDecayMode GetDecayMode()
G4VDecayChannel * GetDecayChannel(G4int index) const
G4int GetVerboseLevel() const
G4bool IsApplicable(const G4ParticleDefinition &)
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
void SetItsRates(G4RadioactiveDecayRates arate)
G4int entries() const
G4IonTable * GetIonTable() const
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
Definition: G4Ions.hh:51
void SetDecayRate(G4int, G4int, G4double, G4int, std::vector< G4double >, std::vector< G4double >)
bool G4bool
Definition: G4Types.hh:79
G4RadioactiveDecayMode
size_t NumberOfTransitions() const
G4float NearestLevelEnergy(G4double energy, size_t index=0) const
G4ParticleDefinition * GetDaughterNucleus()
size_t NearestLevelIndex(G4double energy, size_t index=0) const
G4DecayTable * GetDecayTable(const G4ParticleDefinition *)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
void Insert(G4VDecayChannel *aChannel)
Definition: G4DecayTable.cc:60
G4double GetDaughterExcitation()
static G4ParticleTable * GetParticleTable()
G4double GetPDGLifeTime() const
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
double G4double
Definition: G4Types.hh:76
#define ns
Definition: xmlparse.cc:614
static G4NuclearLevelData * GetInstance()
#define noFloat
Definition: G4Ions.hh:118

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RadioactiveDecay::AddDeexcitationSpectrumForBiasMode ( G4ParticleDefinition apartDef,
G4double  weight,
G4double  currenTime,
std::vector< double > &  weights_v,
std::vector< double > &  times_v,
std::vector< G4DynamicParticle * > &  secondaries_v 
)
protected

Definition at line 2027 of file G4RadioactiveDecay.cc.

2032 {
2033  G4double elevel=((const G4Ions*)(apartDef))->GetExcitationEnergy();
2034  G4double life_time=apartDef->GetPDGLifeTime();
2035  while (life_time <halflifethreshold && elevel>0.) {
2036  G4ITDecay* anITChannel = new G4ITDecay(apartDef, 100., elevel,elevel);
2037  G4DecayProducts* pevap_products = anITChannel->DecayIt(0.);
2038  G4int nb_pevapSecondaries = pevap_products->entries();
2039  for (G4int ind = 0; ind < nb_pevapSecondaries; ind++) {
2040  G4DynamicParticle* a_pevap_secondary= pevap_products->PopProducts();
2041  //Gammas,electrons, alphas coming from excited state
2042  if (a_pevap_secondary->GetDefinition()->GetBaryonNumber() < 5) {
2043  weights_v.push_back(weight);
2044  times_v.push_back(currentTime);
2045  secondaries_v.push_back(a_pevap_secondary);
2046  }
2047  //New excited or ground state
2048  else {
2049  apartDef =a_pevap_secondary->GetDefinition();
2050  elevel=((const G4Ions*)(apartDef))->GetExcitationEnergy();
2051  life_time=apartDef->GetPDGLifeTime();
2052  }
2053  }
2054  delete anITChannel;
2055  }
2056 }
virtual G4DecayProducts * DecayIt(G4double)
Definition: G4ITDecay.cc:81
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
Definition: G4Ions.hh:51
G4DynamicParticle * PopProducts()
G4double GetPDGLifeTime() const
G4int entries() const
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RadioactiveDecay::AddUserDecayDataFile ( G4int  Z,
G4int  A,
G4String  filename 
)

Definition at line 1115 of file G4RadioactiveDecay.cc.

1116 {
1117  if (Z < 1 || A < 2) G4cout << "Z and A not valid!" << G4endl;
1118 
1119  std::ifstream DecaySchemeFile(filename);
1120  if (DecaySchemeFile) {
1121  G4int ID_ion = A*1000 + Z;
1122  theUserRadioactiveDataFiles[ID_ion] = filename;
1123  } else {
1124  G4cout << "The file " << filename << " does not exist!" << G4endl;
1125  }
1126 }
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
#define G4endl
Definition: G4ios.hh:61

Here is the caller graph for this function:

void G4RadioactiveDecay::BuildPhysicsTable ( const G4ParticleDefinition )
virtual

Reimplemented from G4VProcess.

Definition at line 722 of file G4RadioactiveDecay.cc.

723 {
724  if (!isInitialised) {
725  isInitialised = true;
727  G4VAtomDeexcitation* p = theManager->AtomDeexcitation();
728  if (!p) {
730  ed << " Atomic deexcitation is not defined.";
731  G4Exception("G4RadioactiveDecay::BuildPhysicsTable", "HAD_RDM_001",
732  FatalException, ed);
733  /*
734  p = new G4UAtomicDeexcitation();
735  p->SetFluo(true);
736  p->SetAuger(true);
737  p->InitialiseAtomicDeexcitation();
738  theManager->SetAtomDeexcitation(p);
739  */
740  }
741 
743  param->SetUseFilesNEW(true);
744  //param->SetCorrelatedGamma(true); //AR-20Feb2017: Temporary, to fix non-reproducibility problems
745  }
746 }
static G4LossTableManager * Instance()
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
const char * p
Definition: xmltok.h:285
G4DeexPrecoParameters * GetParameters()
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4VAtomDeexcitation * AtomDeexcitation()
static G4NuclearLevelData * GetInstance()

Here is the call graph for this function:

G4ThreeVector G4RadioactiveDecay::ChooseCollimationDirection ( ) const
protected

Definition at line 2001 of file G4RadioactiveDecay.cc.

2001  {
2002  if (origin == forceDecayDirection) return origin; // Don't do collimation
2003  if (forceDecayHalfAngle == 180.*deg) return origin;
2004 
2005  G4ThreeVector dir = forceDecayDirection;
2006 
2007  // Return direction offset by random throw
2008  if (forceDecayHalfAngle > 0.) {
2009  // Generate uniform direction around central axis
2010  G4double phi = 2.*pi*G4UniformRand();
2011  G4double cosMin = std::cos(forceDecayHalfAngle);
2012  G4double cosTheta = (1.-cosMin)*G4UniformRand() + cosMin; // [cosMin,1.)
2013 
2014  dir.setPhi(dir.phi()+phi);
2015  dir.setTheta(dir.theta()+std::acos(cosTheta));
2016  }
2017 
2018 #ifdef G4VERBOSE
2019  if (GetVerboseLevel()>1)
2020  G4cout << " ChooseCollimationDirection returns " << dir << G4endl;
2021 #endif
2022 
2023  return dir;
2024 }
void setPhi(double)
G4int GetVerboseLevel() const
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
double phi() const
void setTheta(double)
double theta() const
#define G4endl
Definition: G4ios.hh:61
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76
static constexpr double deg
Definition: G4SIunits.hh:152

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RadioactiveDecay::CollimateDecay ( G4DecayProducts products)
protected

Definition at line 1958 of file G4RadioactiveDecay.cc.

1958  {
1959  if (origin == forceDecayDirection) return; // No collimation requested
1960  if (180.*deg == forceDecayHalfAngle) return;
1961  if (0 == products || 0 == products->entries()) return;
1962 
1963 #ifdef G4VERBOSE
1964  if (GetVerboseLevel() > 0) G4cout << "Begin of CollimateDecay..." << G4endl;
1965 #endif
1966 
1967  // Particles suitable for directional biasing (for if-blocks below)
1971  static const G4ParticleDefinition* gamma = G4Gamma::Definition();
1974 
1975  G4ThreeVector newDirection; // Re-use to avoid memory churn
1976  for (G4int i=0; i<products->entries(); i++) {
1977  G4DynamicParticle* daughter = (*products)[i];
1978  const G4ParticleDefinition* daughterType =
1979  daughter->GetParticleDefinition();
1980  if (daughterType == electron || daughterType == positron ||
1981  daughterType == neutron || daughterType == gamma ||
1982  daughterType == alpha || daughterType == proton) CollimateDecayProduct(daughter);
1983  }
1984 }
static G4Alpha * Definition()
Definition: G4Alpha.cc:49
static G4Electron * Definition()
Definition: G4Electron.cc:49
G4int GetVerboseLevel() const
static G4Positron * Definition()
Definition: G4Positron.cc:49
int G4int
Definition: G4Types.hh:78
static G4Proton * Definition()
Definition: G4Proton.cc:49
G4GLOB_DLL std::ostream G4cout
void CollimateDecayProduct(G4DynamicParticle *product)
const G4ParticleDefinition * GetParticleDefinition() const
static G4Neutron * Definition()
Definition: G4Neutron.cc:54
#define G4endl
Definition: G4ios.hh:61
G4int entries() const
static constexpr double deg
Definition: G4SIunits.hh:152
static const G4double alpha
static G4Gamma * Definition()
Definition: G4Gamma.cc:49

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RadioactiveDecay::CollimateDecayProduct ( G4DynamicParticle product)
protected

Definition at line 1986 of file G4RadioactiveDecay.cc.

1986  {
1987 #ifdef G4VERBOSE
1988  if (GetVerboseLevel() > 1) {
1989  G4cout << "CollimateDecayProduct for daughter "
1990  << daughter->GetParticleDefinition()->GetParticleName() << G4endl;
1991  }
1992 #endif
1993 
1995  if (origin != collimate) daughter->SetMomentumDirection(collimate);
1996 }
G4int GetVerboseLevel() const
G4GLOB_DLL std::ostream G4cout
G4ThreeVector ChooseCollimationDirection() const
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4RadioactiveDecay::ConvolveSourceTimeProfile ( const G4double  t,
const G4double  tau 
)
protected

Definition at line 394 of file G4RadioactiveDecay.cc.

395 {
396  long double convolvedTime = 0.L;
397  G4int nbin;
398  if ( t > SBin[NSourceBin]) {
399  nbin = NSourceBin;
400  } else {
401  nbin = 0;
402 
403  G4int loop = 0;
405  ed << " While count exceeded " << G4endl;
406  while (t > SBin[nbin]) { /* Loop checking, 01.09.2015, D.Wright */
407  loop++;
408  if (loop > 1000) {
409  G4Exception("G4RadioactiveDecay::ConvolveSourceTimeProfile()",
410  "HAD_RDM_100", JustWarning, ed);
411  break;
412  }
413 
414  nbin++;
415  }
416  nbin--;
417  }
418  long double lt = t ;
419  long double ltau = tau;
420  // G4cout << " Convolve: tau = " << tau << G4endl;
421  if (nbin > 0) {
422  for (G4int i = 0; i < nbin; i++) {
423  convolvedTime += (long double)SProfile[i] *
424  (std::exp(-(lt-(long double)SBin[i+1])/ltau)-std::exp(-(lt-(long double)SBin[i])/ltau));
425  }
426  }
427  convolvedTime += (long double)SProfile[nbin] * (1.L-std::exp(-(lt-(long double)SBin[nbin])/ltau));
428  // Is the above line necessary? If so, the 1.L looks incorrect - should be an exp
429  // Also, it looks like the final integral should be multiplied by ltau
430 
431  if (convolvedTime < 0.) {
432  G4cout << " Convolved time =: " << convolvedTime << " reset to zero! " << G4endl;
433  G4cout << " t = " << t << " tau = " << tau << G4endl;
434  G4cout << SBin[nbin] << " " << SBin[0] << G4endl;
435  convolvedTime = 0.;
436  }
437 #ifdef G4VERBOSE
438  if (GetVerboseLevel() > 1)
439  G4cout << " Convolved time: " << convolvedTime << G4endl;
440 #endif
441  return (G4double)convolvedTime ;
442 }
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4int GetVerboseLevel() const
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
static constexpr double L
Definition: G4SIunits.hh:124

Here is the call graph for this function:

Here is the caller graph for this function:

G4VParticleChange * G4RadioactiveDecay::DecayIt ( const G4Track theTrack,
const G4Step theStep 
)

Definition at line 1572 of file G4RadioactiveDecay.cc.

1573 {
1574  // Initialize G4ParticleChange object, get particle details and decay table
1575 
1576  fParticleChangeForRadDecay.Initialize(theTrack);
1577  fParticleChangeForRadDecay.ProposeWeight(theTrack.GetWeight());
1578  const G4DynamicParticle* theParticle = theTrack.GetDynamicParticle();
1579  const G4ParticleDefinition* theParticleDef = theParticle->GetDefinition();
1580 
1581  // First check whether RDM applies to the current logical volume
1582  if (!isAllVolumesMode) {
1583  if (!std::binary_search(ValidVolumes.begin(), ValidVolumes.end(),
1584  theTrack.GetVolume()->GetLogicalVolume()->GetName())) {
1585 #ifdef G4VERBOSE
1586  if (GetVerboseLevel()>0) {
1587  G4cout <<"G4RadioactiveDecay::DecayIt : "
1588  << theTrack.GetVolume()->GetLogicalVolume()->GetName()
1589  << " is not selected for the RDM"<< G4endl;
1590  G4cout << " There are " << ValidVolumes.size() << " volumes" << G4endl;
1591  G4cout << " The Valid volumes are " << G4endl;
1592  for (size_t i = 0; i< ValidVolumes.size(); i++)
1593  G4cout << ValidVolumes[i] << G4endl;
1594  }
1595 #endif
1596  fParticleChangeForRadDecay.SetNumberOfSecondaries(0);
1597 
1598  // Kill the parent particle.
1599  fParticleChangeForRadDecay.ProposeTrackStatus(fStopAndKill) ;
1600  fParticleChangeForRadDecay.ProposeLocalEnergyDeposit(0.0);
1602  return &fParticleChangeForRadDecay;
1603  }
1604  }
1605 
1606  // Now check if particle is valid for RDM
1607  if (!(IsApplicable(*theParticleDef) ) ) {
1608  // Particle is not an ion or is outside the nucleuslimits for decay
1609 
1610 #ifdef G4VERBOSE
1611  if (GetVerboseLevel()>0) {
1612  G4cerr << "G4RadioactiveDecay::DecayIt : "
1613  << theParticleDef->GetParticleName()
1614  << " is not a valid nucleus for the RDM"<< G4endl;
1615  }
1616 #endif
1617  fParticleChangeForRadDecay.SetNumberOfSecondaries(0);
1618 
1619  // Kill the parent particle
1620  fParticleChangeForRadDecay.ProposeTrackStatus(fStopAndKill) ;
1621  fParticleChangeForRadDecay.ProposeLocalEnergyDeposit(0.0);
1623  return &fParticleChangeForRadDecay;
1624  }
1625  G4DecayTable* theDecayTable = GetDecayTable(theParticleDef);
1626 
1627  if (theDecayTable == 0 || theDecayTable->entries() == 0) {
1628  // No data in the decay table. Set particle change parameters
1629  // to indicate this.
1630 #ifdef G4VERBOSE
1631  if (GetVerboseLevel() > 0) {
1632  G4cerr <<"G4RadioactiveDecay::DecayIt : decay table not defined for ";
1633  G4cerr <<theParticleDef->GetParticleName() <<G4endl;
1634  }
1635 #endif
1636  fParticleChangeForRadDecay.SetNumberOfSecondaries(0);
1637 
1638  // Kill the parent particle.
1639  fParticleChangeForRadDecay.ProposeTrackStatus(fStopAndKill) ;
1640  fParticleChangeForRadDecay.ProposeLocalEnergyDeposit(0.0);
1642  return &fParticleChangeForRadDecay;
1643 
1644  } else {
1645  // Data found. Try to decay nucleus
1646  G4double energyDeposit = 0.0;
1647  G4double finalGlobalTime = theTrack.GetGlobalTime();
1648  G4double finalLocalTime = theTrack.GetLocalTime();
1649  G4int index;
1650  G4ThreeVector currentPosition;
1651  currentPosition = theTrack.GetPosition();
1652 
1653  // Check whether use Analogue or VR implementation
1654  if (AnalogueMC) {
1655 #ifdef G4VERBOSE
1656  if (GetVerboseLevel() > 0)
1657  G4cout <<"DecayIt: Analogue MC version " << G4endl;
1658 # endif
1659 
1660  G4DecayProducts* products = DoDecay(*theParticleDef);
1661 
1662  // Check if the product is the same as input and kill the track if
1663  // necessary to prevent infinite loop (11/05/10, F.Lei)
1664  if ( products->entries() == 1) {
1665  fParticleChangeForRadDecay.SetNumberOfSecondaries(0);
1666  fParticleChangeForRadDecay.ProposeTrackStatus(fStopAndKill);
1667  fParticleChangeForRadDecay.ProposeLocalEnergyDeposit(0.0);
1669  return &fParticleChangeForRadDecay;
1670  }
1671 
1672  // Get parent particle information and boost the decay products to the
1673  // laboratory frame based on this information.
1674 
1675  //The Parent Energy used for the boost should be the total energy of
1676  // the nucleus of the parent ion without the energy of the shell electrons
1677  // (correction for bug 1359 by L. Desorgher)
1678  G4double ParentEnergy = theParticle->GetKineticEnergy()
1679  + theParticle->GetParticleDefinition()->GetPDGMass();
1680  G4ThreeVector ParentDirection(theParticle->GetMomentumDirection());
1681 
1682  if (theTrack.GetTrackStatus() == fStopButAlive) {
1683  //this condition seems to be always True, further investigation is needed (L.Desorgher)
1684 
1685  // The particle is decayed at rest.
1686  // since the time is still for rest particle in G4 we need to add the
1687  // additional time lapsed between the particle come to rest and the
1688  // actual decay. This time is simply sampled with the mean-life of
1689  // the particle. But we need to protect the case PDGTime < 0.
1690  // (F.Lei 11/05/10)
1691  G4double temptime = -std::log( G4UniformRand())
1692  *theParticleDef->GetPDGLifeTime();
1693  if (temptime < 0.) temptime = 0.;
1694  finalGlobalTime += temptime;
1695  finalLocalTime += temptime;
1696  energyDeposit += theParticle->GetKineticEnergy();
1697  }
1698  products->Boost(ParentEnergy, ParentDirection);
1699 
1700  // Add products in theParticleChangeForRadDecay.
1701  G4int numberOfSecondaries = products->entries();
1702  fParticleChangeForRadDecay.SetNumberOfSecondaries(numberOfSecondaries);
1703 #ifdef G4VERBOSE
1704  if (GetVerboseLevel()>1) {
1705  G4cout <<"G4RadioactiveDecay::DecayIt : Decay vertex :";
1706  G4cout <<" Time: " <<finalGlobalTime/ns <<"[ns]";
1707  G4cout <<" X:" <<(theTrack.GetPosition()).x() /cm <<"[cm]";
1708  G4cout <<" Y:" <<(theTrack.GetPosition()).y() /cm <<"[cm]";
1709  G4cout <<" Z:" <<(theTrack.GetPosition()).z() /cm <<"[cm]";
1710  G4cout << G4endl;
1711  G4cout <<"G4Decay::DecayIt : decay products in Lab. Frame" <<G4endl;
1712  products->DumpInfo();
1713  products->IsChecked();
1714  }
1715 #endif
1716  for (index=0; index < numberOfSecondaries; index++) {
1717  G4Track* secondary = new G4Track(products->PopProducts(),
1718  finalGlobalTime, currentPosition);
1719  secondary->SetGoodForTrackingFlag();
1720  secondary->SetTouchableHandle(theTrack.GetTouchableHandle());
1721  fParticleChangeForRadDecay.AddSecondary(secondary);
1722  }
1723  delete products;
1724  // end of analogue MC algorithm
1725 
1726  } else {
1727  // Variance Reduction Method
1728 #ifdef G4VERBOSE
1729  if (GetVerboseLevel()>0)
1730  G4cout << "DecayIt: Variance Reduction version " << G4endl;
1731 #endif
1732  if (!IsRateTableReady(*theParticleDef)) {
1733  // if the decayrates are not ready, calculate them and
1734  // add to the rate table vector
1735  AddDecayRateTable(*theParticleDef);
1736  }
1737  //retrieve the rates
1738  GetDecayRateTable(*theParticleDef);
1739 
1740  // declare some of the variables required in the implementation
1741  G4ParticleDefinition* parentNucleus;
1742  G4IonTable* theIonTable;
1743  G4int PZ;
1744  G4int PA;
1745  G4double PE;
1746  G4String keyName;
1747  std::vector<G4double> PT;
1748  std::vector<G4double> PR;
1749  G4double taotime;
1750  long double decayRate;
1751 
1752  size_t i;
1753  size_t j;
1754  G4int numberOfSecondaries;
1755  G4int totalNumberOfSecondaries = 0;
1756  G4double currentTime = 0.;
1757  G4int ndecaych;
1758  G4DynamicParticle* asecondaryparticle;
1759  std::vector<G4DynamicParticle*> secondaryparticles;
1760  std::vector<G4double> pw;
1761  std::vector<G4double> ptime;
1762  pw.clear();
1763  ptime.clear();
1764 
1765  //now apply the nucleus splitting
1766  for (G4int n = 0; n < NSplit; n++) {
1767  // Get the decay time following the decay probability function
1768  // suppllied by user
1769  G4double theDecayTime = GetDecayTime();
1770  G4int nbin = GetDecayTimeBin(theDecayTime);
1771 
1772  // calculate the first part of the weight function
1773  G4double weight1 = 1.;
1774  if (nbin == 1) {
1775  weight1 = 1./DProfile[nbin-1]
1776  *(DBin[nbin]-DBin[nbin-1])/NSplit;
1777  } else if (nbin > 1) {
1778  weight1 = 1./(DProfile[nbin]-DProfile[nbin-2])
1779  *(DBin[nbin]-DBin[nbin-1])/NSplit;
1780  }
1781 
1782  // it should be calculated in seconds
1783  weight1 /= s ;
1784 
1785  // loop over all the possible secondaries of the nucleus
1786  // the first one is itself.
1787  for (i = 0; i < theDecayRateVector.size(); i++) {
1788  PZ = theDecayRateVector[i].GetZ();
1789  PA = theDecayRateVector[i].GetA();
1790  PE = theDecayRateVector[i].GetE();
1791  PT = theDecayRateVector[i].GetTaos();
1792  PR = theDecayRateVector[i].GetDecayRateC();
1793 
1794  // Calculate the decay rate of the isotope
1795  // decayRate is the radioactivity of isotope (PZ,PA,PE) at the
1796  // time 'theDecayTime'
1797  // it will be used to calculate the statistical weight of the
1798  // decay products of this isotope
1799 
1800 // G4cout <<"PA= "<< PA << " PZ= " << PZ << " PE= "<< PE <<G4endl;
1801  decayRate = 0.L;
1802  for (j = 0; j < PT.size(); j++) {
1803 // G4cout << " RDM::DecayIt: tau input to Convolve: " << PT[j] << G4endl;
1804  taotime = ConvolveSourceTimeProfile(theDecayTime,PT[j]);
1805 // taotime = GetTaoTime(theDecayTime,PT[j]);
1806  decayRate -= PR[j] * (long double)taotime;
1807  // Eq.4.23 of of the TN
1808  // note the negative here is required as the rate in the
1809  // equation is defined to be negative,
1810  // i.e. decay away, but we need positive value here.
1811 
1812  // G4cout << j << "\t"<< PT[j]/s <<"\t"<<PR[j]<< "\t"
1813  // << decayRate << G4endl;
1814  }
1815 
1816  // add the isotope to the radioactivity tables
1817  // G4cout <<theDecayTime/s <<"\t"<<nbin<<G4endl;
1818  // G4cout << theTrack.GetWeight() <<"\t"<<weight1<<"\t"<<decayRate<< G4endl;
1819  theRadioactivityTables[decayWindows[nbin-1]]->AddIsotope(PZ,PA,PE,weight1*decayRate,theTrack.GetWeight());
1820 
1821  // Now calculate the statistical weight
1822  // One needs to fold the source bias function with the decaytime
1823  // also need to include the track weight! (F.Lei, 28/10/10)
1824  G4double weight = weight1*decayRate*theTrack.GetWeight();
1825 
1826  // decay the isotope
1828  parentNucleus = theIonTable->GetIon(PZ,PA,PE);
1829 
1830  // Create a temprary products buffer.
1831  // Its contents to be transfered to the products at the end of the loop
1832  G4DecayProducts* tempprods = 0;
1833 
1834  // Decide whether to apply branching ratio bias or not
1835  if (BRBias) {
1836  G4DecayTable* decayTable = GetDecayTable(parentNucleus);
1837 
1838  ndecaych = G4int(decayTable->entries()*G4UniformRand());
1839  G4VDecayChannel* theDecayChannel = decayTable->GetDecayChannel(ndecaych);
1840  if (theDecayChannel == 0) {
1841  // Decay channel not found.
1842 #ifdef G4VERBOSE
1843  if (GetVerboseLevel()>0) {
1844  G4cerr << " G4RadioactiveDecay::DoIt : cannot determine decay channel ";
1845  G4cerr << " for this nucleus; decay as if no biasing active ";
1846  G4cerr << G4endl;
1847  decayTable ->DumpInfo();
1848  }
1849 #endif
1850  tempprods = DoDecay(*parentNucleus); // DHW 6 Dec 2010 - do decay as if no biasing
1851  // to avoid deref of temppprods = 0
1852  } else {
1853  // A decay channel has been identified, so execute the DecayIt.
1854  G4double tempmass = parentNucleus->GetPDGMass();
1855  tempprods = theDecayChannel->DecayIt(tempmass);
1856  weight *= (theDecayChannel->GetBR())*(decayTable->entries());
1857  }
1858  } else {
1859  tempprods = DoDecay(*parentNucleus);
1860  }
1861 
1862 
1863  // save the secondaries for buffers
1864  numberOfSecondaries = tempprods->entries();
1865  currentTime = finalGlobalTime + theDecayTime;
1866  for (index = 0; index < numberOfSecondaries; index++) {
1867  asecondaryparticle = tempprods->PopProducts();
1868  if (asecondaryparticle->GetDefinition()->GetBaryonNumber() < 5) {
1869  pw.push_back(weight);
1870  ptime.push_back(currentTime);
1871  secondaryparticles.push_back(asecondaryparticle);
1872  }
1873  //Generate gammas and XRays from excited nucleus, added by L.Desorgher
1874  else if (((const G4Ions*)(asecondaryparticle->GetDefinition()))->GetExcitationEnergy()>0. && weight>0.){//Compute the gamma
1875  G4ParticleDefinition* apartDef =asecondaryparticle->GetDefinition();
1876  AddDeexcitationSpectrumForBiasMode(apartDef,weight,currentTime,pw,ptime,secondaryparticles);
1877  }
1878  }
1879  delete tempprods;
1880 
1881  } // end of i loop
1882  } // end of n loop
1883 
1884  // now deal with the secondaries in the two stl containers
1885  // and submmit them back to the tracking manager
1886  totalNumberOfSecondaries = pw.size();
1887  fParticleChangeForRadDecay.SetNumberOfSecondaries(totalNumberOfSecondaries);
1888  for (index=0; index < totalNumberOfSecondaries; index++) {
1889  G4Track* secondary = new G4Track(secondaryparticles[index],
1890  ptime[index], currentPosition);
1891  secondary->SetGoodForTrackingFlag();
1892  secondary->SetTouchableHandle(theTrack.GetTouchableHandle());
1893  secondary->SetWeight(pw[index]);
1894  fParticleChangeForRadDecay.AddSecondary(secondary);
1895  }
1896  // make sure the original track is set to stop and its kinematic energy collected
1897  //
1898  //theTrack.SetTrackStatus(fStopButAlive);
1899  //energyDeposit += theParticle->GetKineticEnergy();
1900 
1901  } // End of Variance Reduction
1902 
1903  // Kill the parent particle
1904  fParticleChangeForRadDecay.ProposeTrackStatus(fStopAndKill) ;
1905  fParticleChangeForRadDecay.ProposeLocalEnergyDeposit(energyDeposit);
1906  fParticleChangeForRadDecay.ProposeLocalTime(finalLocalTime);
1907  // Reset NumberOfInteractionLengthLeft.
1909 
1910  return &fParticleChangeForRadDecay ;
1911  }
1912 }
G4DecayProducts * DoDecay(const G4ParticleDefinition &theParticleDef)
G4double GetLocalTime() const
G4double GetKineticEnergy() const
G4bool IsRateTableReady(const G4ParticleDefinition &)
const G4DynamicParticle * GetDynamicParticle() const
G4double ConvolveSourceTimeProfile(const G4double, const G4double)
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:503
const G4ThreeVector & GetPosition() const
void AddSecondary(G4Track *aSecondary)
G4TrackStatus GetTrackStatus() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
G4VDecayChannel * GetDecayChannel(G4int index) const
G4int GetVerboseLevel() const
G4ParticleDefinition * GetDefinition() const
tuple x
Definition: test.py:50
G4bool IsApplicable(const G4ParticleDefinition &)
void Boost(G4double totalEnergy, const G4ThreeVector &momentumDirection)
void ClearNumberOfInteractionLengthLeft()
Definition: G4VProcess.hh:447
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
void SetWeight(G4double aValue)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
G4int entries() const
const XML_Char * s
Definition: expat.h:262
void ProposeWeight(G4double finalWeight)
virtual void Initialize(const G4Track &)
G4IonTable * GetIonTable() const
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
Definition: G4Ions.hh:51
void DumpInfo() const
const G4ThreeVector & GetMomentumDirection() const
void DumpInfo() const
static constexpr double cm
Definition: G4SIunits.hh:119
G4int GetDecayTimeBin(const G4double aDecayTime)
const G4ParticleDefinition * GetParticleDefinition() const
const G4int n
G4double GetGlobalTime() const
const G4TouchableHandle & GetTouchableHandle() const
G4DecayTable * GetDecayTable(const G4ParticleDefinition *)
G4LogicalVolume * GetLogicalVolume() const
G4double GetPDGMass() const
static G4ParticleTable * GetParticleTable()
void GetDecayRateTable(const G4ParticleDefinition &)
void SetNumberOfSecondaries(G4int totSecondaries)
G4DynamicParticle * PopProducts()
void AddDecayRateTable(const G4ParticleDefinition &)
G4VPhysicalVolume * GetVolume() const
tuple z
Definition: test.py:28
G4double GetWeight() const
G4double GetPDGLifeTime() const
#define G4endl
Definition: G4ios.hh:61
G4bool IsChecked() const
G4int entries() const
const G4String & GetName() const
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
void AddDeexcitationSpectrumForBiasMode(G4ParticleDefinition *apartDef, G4double weight, G4double currenTime, std::vector< double > &weights_v, std::vector< double > &times_v, std::vector< G4DynamicParticle * > &secondaries_v)
#define ns
Definition: xmlparse.cc:614
void SetGoodForTrackingFlag(G4bool value=true)
G4GLOB_DLL std::ostream G4cerr
void G4RadioactiveDecay::DeselectAllVolumes ( )

Definition at line 344 of file G4RadioactiveDecay.cc.

345 {
346  ValidVolumes.clear();
347  isAllVolumesMode=false;
348 #ifdef G4VERBOSE
349  if (GetVerboseLevel() > 0) G4cout << "RDM removed from all volumes" << G4endl;
350 #endif
351 }
G4int GetVerboseLevel() const
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

void G4RadioactiveDecay::DeselectAVolume ( const G4String  aVolume)

Definition at line 289 of file G4RadioactiveDecay.cc.

290 {
291  G4LogicalVolumeStore* theLogicalVolumes;
292  G4LogicalVolume* volume;
293  theLogicalVolumes = G4LogicalVolumeStore::GetInstance();
294  for (size_t i = 0; i < theLogicalVolumes->size(); i++){
295  volume=(*theLogicalVolumes)[i];
296  if (volume->GetName() == aVolume) {
297  std::vector<G4String>::iterator location;
298  location = std::find(ValidVolumes.begin(),ValidVolumes.end(),aVolume);
299  if (location != ValidVolumes.end()) {
300  ValidVolumes.erase(location);
301  std::sort(ValidVolumes.begin(), ValidVolumes.end());
302  isAllVolumesMode =false;
303  } else {
304  G4cerr << " DeselectVolume:" << aVolume << " is not in the list "
305  << G4endl;
306  }
307 #ifdef G4VERBOSE
308  if (GetVerboseLevel() > 0)
309  G4cout << " DeselectVolume: " << aVolume << " is removed from list "
310  << G4endl;
311 #endif
312  } else if (i == theLogicalVolumes->size()) {
313  G4cerr << " DeselectVolume:" << aVolume
314  << "is not a valid logical volume name" << G4endl;
315  }
316  }
317 }
G4int GetVerboseLevel() const
G4GLOB_DLL std::ostream G4cout
static G4LogicalVolumeStore * GetInstance()
#define G4endl
Definition: G4ios.hh:61
const G4String & GetName() const
G4GLOB_DLL std::ostream G4cerr

Here is the call graph for this function:

G4DecayProducts * G4RadioactiveDecay::DoDecay ( const G4ParticleDefinition theParticleDef)
protected

Definition at line 1916 of file G4RadioactiveDecay.cc.

1917 {
1918  G4DecayProducts* products = 0;
1919  G4DecayTable* theDecayTable = GetDecayTable(&theParticleDef);
1920 
1921  // Choose a decay channel.
1922 #ifdef G4VERBOSE
1923  if (GetVerboseLevel() > 0) G4cout << "Select a channel..." << G4endl;
1924 #endif
1925 
1926  // G4DecayTable::SelectADecayChannel checks to see if sum of daughter masses
1927  // exceeds parent mass. Pass it the parent mass + maximum Q value to account
1928  // for difference in mass defect.
1929  G4double parentPlusQ = theParticleDef.GetPDGMass() + 30.*MeV;
1930  G4VDecayChannel* theDecayChannel = theDecayTable->SelectADecayChannel(parentPlusQ);
1931 
1932  if (theDecayChannel == 0) {
1933  // Decay channel not found.
1935  ed << " Cannot determine decay channel for " << theParticleDef.GetParticleName() << G4endl;
1936  G4Exception("G4RadioactiveDecay::DoDecay", "HAD_RDM_013",
1937  FatalException, ed);
1938  } else {
1939  // A decay channel has been identified, so execute the DecayIt.
1940 #ifdef G4VERBOSE
1941  if (GetVerboseLevel() > 1) {
1942  G4cerr << "G4RadioactiveDecay::DoIt : selected decay channel addr:";
1943  G4cerr << theDecayChannel << G4endl;
1944  }
1945 #endif
1946  products = theDecayChannel->DecayIt(theParticleDef.GetPDGMass() );
1947 
1948  // Apply directional bias if requested by user
1949  CollimateDecay(products);
1950  }
1951 
1952  return products;
1953 }
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4int GetVerboseLevel() const
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
G4VDecayChannel * SelectADecayChannel(G4double parentMass=-1.)
Definition: G4DecayTable.cc:81
G4DecayTable * GetDecayTable(const G4ParticleDefinition *)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4double GetPDGMass() const
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
virtual G4DecayProducts * DecayIt(G4double parentMass=-1.0)=0
double G4double
Definition: G4Types.hh:76
void CollimateDecay(G4DecayProducts *products)
G4GLOB_DLL std::ostream G4cerr

Here is the call graph for this function:

Here is the caller graph for this function:

const G4ThreeVector& G4RadioactiveDecay::GetDecayDirection ( ) const
inline

Definition at line 216 of file G4RadioactiveDecay.hh.

216  {
217  return forceDecayDirection;
218  }
G4double G4RadioactiveDecay::GetDecayHalfAngle ( ) const
inline

Definition at line 224 of file G4RadioactiveDecay.hh.

224 {return forceDecayHalfAngle;}
void G4RadioactiveDecay::GetDecayRateTable ( const G4ParticleDefinition aParticle)

Definition at line 370 of file G4RadioactiveDecay.cc.

371 {
372  G4String aParticleName = aParticle.GetParticleName();
373 
374  for (size_t i = 0; i < theDecayRateTableVector.size(); i++) {
375  if (theDecayRateTableVector[i].GetIonName() == aParticleName) {
376  theDecayRateVector = theDecayRateTableVector[i].GetItsRates();
377  }
378  }
379 #ifdef G4VERBOSE
380  if (GetVerboseLevel() > 0) {
381  G4cout << "The DecayRate Table for " << aParticleName << " is selected."
382  << G4endl;
383  }
384 #endif
385 }
G4int GetVerboseLevel() const
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

G4DecayTable * G4RadioactiveDecay::GetDecayTable ( const G4ParticleDefinition aNucleus)

Definition at line 249 of file G4RadioactiveDecay.cc.

250 {
251  G4String key = aNucleus->GetParticleName();
252  DecayTableMap::iterator table_ptr = dkmap->find(key);
253 
254  G4DecayTable* theDecayTable = 0;
255  if (table_ptr == dkmap->end() ) { // If table not there,
256  theDecayTable = LoadDecayTable(*aNucleus); // load from file and
257  if(theDecayTable) (*dkmap)[key] = theDecayTable; // store in library
258  } else {
259  theDecayTable = table_ptr->second;
260  }
261 
262  return theDecayTable;
263 }
const G4String & GetParticleName() const
G4DecayTable * LoadDecayTable(const G4ParticleDefinition &theParentNucleus)

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4RadioactiveDecay::GetDecayTime ( )
protected

Definition at line 562 of file G4RadioactiveDecay.cc.

563 {
564  G4double decaytime = 0.;
565  G4double rand = G4UniformRand();
566  G4int i = 0;
567 
568  G4int loop = 0;
570  ed << " While count exceeded " << G4endl;
571  while ( DProfile[i] < rand) { /* Loop checking, 01.09.2015, D.Wright */
572  i++;
573  loop++;
574  if (loop > 100000) {
575  G4Exception("G4RadioactiveDecay::GetDecayTime()", "HAD_RDM_100", JustWarning, ed);
576  break;
577  }
578  }
579 
580  rand = G4UniformRand();
581  decaytime = DBin[i] + rand*(DBin[i+1]-DBin[i]);
582 #ifdef G4VERBOSE
583  if (GetVerboseLevel() > 1)
584  G4cout <<" Decay time: " <<decaytime/s <<"[s]" <<G4endl;
585 #endif
586  return decaytime;
587 }
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4int GetVerboseLevel() const
int G4int
Definition: G4Types.hh:78
const XML_Char * s
Definition: expat.h:262
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4int G4RadioactiveDecay::GetDecayTimeBin ( const G4double  aDecayTime)
protected

Definition at line 590 of file G4RadioactiveDecay.cc.

591 {
592  G4int i = 0;
593 
594  G4int loop = 0;
596  ed << " While count exceeded " << G4endl;
597  while ( aDecayTime > DBin[i] ) { /* Loop checking, 01.09.2015, D.Wright */
598  i++;
599  loop++;
600  if (loop > 100000) {
601  G4Exception("G4RadioactiveDecay::GetDecayTimeBin()", "HAD_RDM_100", JustWarning, ed);
602  break;
603  }
604  }
605 
606  return i;
607 }
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
int G4int
Definition: G4Types.hh:78
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4RadioactiveDecay::GetMeanFreePath ( const G4Track theTrack,
G4double  previousStepSize,
G4ForceCondition condition 
)
protectedvirtual

Implements G4VRestDiscreteProcess.

Definition at line 657 of file G4RadioactiveDecay.cc.

659 {
660  const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
661  const G4ParticleDefinition* aParticleDef = aParticle->GetDefinition();
662  G4double tau = aParticleDef->GetPDGLifeTime();
663  G4double aMass = aParticle->GetMass();
664 
665 #ifdef G4VERBOSE
666  if (GetVerboseLevel() > 2) {
667  G4cout << "G4RadioactiveDecay::GetMeanFreePath() " << G4endl;
668  G4cout << " KineticEnergy: " << aParticle->GetKineticEnergy()/GeV
669  << " GeV, Mass: " << aMass/GeV << " GeV, tau: " << tau << " ns "
670  << G4endl;
671  }
672 #endif
673  G4double pathlength = DBL_MAX;
674  if (tau != -1) {
675  // Ion can decay
676 
677  if (tau < -1000.0) {
678  pathlength = DBL_MIN; // nuclide had very short lifetime or wasn't in table
679 
680  } else if (tau < 0.0) {
681  G4cout << aParticleDef->GetParticleName() << " has lifetime " << tau << G4endl;
683  ed << "Ion has negative lifetime " << tau
684  << " but is not stable. Setting mean free path to DBL_MAX" << G4endl;
685  G4Exception("G4RadioactiveDecay::GetMeanFreePath()", "HAD_RDM_011",
686  JustWarning, ed);
687  pathlength = DBL_MAX;
688 
689  } else {
690  // Calculate mean free path
691  G4double betaGamma = aParticle->GetTotalMomentum()/aMass;
692  pathlength = c_light*tau*betaGamma;
693 
694  if (pathlength < DBL_MIN) {
695  pathlength = DBL_MIN;
696 #ifdef G4VERBOSE
697  if (GetVerboseLevel() > 2) {
698  G4cout << "G4Decay::GetMeanFreePath: "
699  << aParticleDef->GetParticleName()
700  << " stops, kinetic energy = "
701  << aParticle->GetKineticEnergy()/keV <<" keV " << G4endl;
702  }
703 #endif
704  }
705  }
706  }
707 
708 #ifdef G4VERBOSE
709  if (GetVerboseLevel() > 1) {
710  G4cout << "mean free path: "<< pathlength/m << " m" << G4endl;
711  }
712 #endif
713  return pathlength;
714 }
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4double GetKineticEnergy() const
G4int GetVerboseLevel() const
G4ParticleDefinition * GetDefinition() const
const G4String & GetParticleName() const
G4double GetTotalMomentum() const
G4GLOB_DLL std::ostream G4cout
static constexpr double m
Definition: G4SIunits.hh:129
G4double GetMass() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define DBL_MIN
Definition: templates.hh:75
static constexpr double GeV
Definition: G4SIunits.hh:217
G4double GetPDGLifeTime() const
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83
static constexpr double keV
Definition: G4SIunits.hh:216
float c_light
Definition: hepunit.py:257

Here is the call graph for this function:

G4double G4RadioactiveDecay::GetMeanLifeTime ( const G4Track theTrack,
G4ForceCondition condition 
)
protectedvirtual

Implements G4VRestDiscreteProcess.

Definition at line 615 of file G4RadioactiveDecay.cc.

617 {
618  // For variance reduction the time is set to 0 so as to force the particle
619  // to decay immediately.
620  // In analogueMC mode it returns the particle's mean-life.
621 
622  G4double meanlife = 0.;
623  if (AnalogueMC) {
624  const G4DynamicParticle* theParticle = theTrack.GetDynamicParticle();
625  const G4ParticleDefinition* theParticleDef = theParticle->GetDefinition();
626  G4double theLife = theParticleDef->GetPDGLifeTime();
627 #ifdef G4VERBOSE
628  if (GetVerboseLevel() > 2) {
629  G4cout << "G4RadioactiveDecay::GetMeanLifeTime() " << G4endl;
630  G4cout << "KineticEnergy: " << theParticle->GetKineticEnergy()/GeV
631  << " GeV, Mass: " << theParticle->GetMass()/GeV
632  << " GeV, Life time: " << theLife/ns << " ns " << G4endl;
633  }
634 #endif
635  if (theParticleDef->GetPDGStable()) {meanlife = DBL_MAX;}
636  else if (theLife < 0.0) {meanlife = DBL_MAX;}
637  else {meanlife = theLife;}
638  // Set meanlife to zero for excited istopes which are not in the
639  // RDM database
640  if (((const G4Ions*)(theParticleDef))->GetExcitationEnergy() > 0. &&
641  meanlife == DBL_MAX) {meanlife = 0.;}
642  }
643 #ifdef G4VERBOSE
644  if (GetVerboseLevel() > 1)
645  G4cout << " mean life time: " << meanlife/s << " s " << G4endl;
646 #endif
647 
648  return meanlife;
649 }
G4double GetKineticEnergy() const
const G4DynamicParticle * GetDynamicParticle() const
G4bool GetPDGStable() const
G4int GetVerboseLevel() const
G4ParticleDefinition * GetDefinition() const
const XML_Char * s
Definition: expat.h:262
G4GLOB_DLL std::ostream G4cout
Definition: G4Ions.hh:51
G4double GetMass() const
static constexpr double GeV
Definition: G4SIunits.hh:217
G4double GetPDGLifeTime() const
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83
#define ns
Definition: xmlparse.cc:614

Here is the call graph for this function:

G4NucleusLimits G4RadioactiveDecay::GetNucleusLimits ( ) const
inline

Definition at line 178 of file G4RadioactiveDecay.hh.

179  {return theNucleusLimits;}
G4int G4RadioactiveDecay::GetSplitNuclei ( )
inline

Definition at line 209 of file G4RadioactiveDecay.hh.

209 {return NSplit;}
std::vector<G4RadioactivityTable*> G4RadioactiveDecay::GetTheRadioactivityTables ( )
inline

Definition at line 155 of file G4RadioactiveDecay.hh.

156  {return theRadioactivityTables;}

Here is the caller graph for this function:

G4int G4RadioactiveDecay::GetVerboseLevel ( ) const
inline

Definition at line 170 of file G4RadioactiveDecay.hh.

170 {return verboseLevel;}

Here is the caller graph for this function:

G4bool G4RadioactiveDecay::IsAnalogueMonteCarlo ( )
inline

Definition at line 193 of file G4RadioactiveDecay.hh.

193 {return AnalogueMC;}

Here is the caller graph for this function:

G4bool G4RadioactiveDecay::IsApplicable ( const G4ParticleDefinition aParticle)
virtual

Reimplemented from G4VProcess.

Definition at line 227 of file G4RadioactiveDecay.cc.

228 {
229  // All particles other than G4Ions, are rejected by default
230  if (((const G4Ions*)(&aParticle))->GetExcitationEnergy() > 0.) {return true;}
231  if (aParticle.GetParticleName() == "GenericIon") {
232  return true;
233  } else if (!(aParticle.GetParticleType() == "nucleus")
234  || aParticle.GetPDGLifeTime() < 0. ) {
235  return false;
236  }
237 
238  // Determine whether the nuclide falls into the correct A and Z range
239  G4int A = ((const G4Ions*) (&aParticle))->GetAtomicMass();
240  G4int Z = ((const G4Ions*) (&aParticle))->GetAtomicNumber();
241 
242  if (A > theNucleusLimits.GetAMax() || A < theNucleusLimits.GetAMin())
243  {return false;}
244  else if (Z > theNucleusLimits.GetZMax() || Z < theNucleusLimits.GetZMin())
245  {return false;}
246  return true;
247 }
G4int GetZMax() const
int G4int
Definition: G4Types.hh:78
G4int GetAMin() const
const G4String & GetParticleName() const
double A(double temperature)
Definition: G4Ions.hh:51
G4int GetAMax() const
const G4String & GetParticleType() const
G4double GetPDGLifeTime() const
G4int GetZMin() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4RadioactiveDecay::IsRateTableReady ( const G4ParticleDefinition aParticle)

Definition at line 355 of file G4RadioactiveDecay.cc.

356 {
357  // Check whether the radioactive decay rates table for the ion has already
358  // been calculated.
359  G4String aParticleName = aParticle.GetParticleName();
360  for (size_t i = 0; i < theDecayRateTableVector.size(); i++) {
361  if (theDecayRateTableVector[i].GetIonName() == aParticleName) return true;
362  }
363  return false;
364 }
const G4String & GetParticleName() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4DecayTable * G4RadioactiveDecay::LoadDecayTable ( const G4ParticleDefinition theParentNucleus)

Definition at line 756 of file G4RadioactiveDecay.cc.

757 {
758  // Generate input data file name using Z and A of the parent nucleus
759  // file containing radioactive decay data.
760  G4int A = ((const G4Ions*)(&theParentNucleus))->GetAtomicMass();
761  G4int Z = ((const G4Ions*)(&theParentNucleus))->GetAtomicNumber();
762 
763  G4double levelEnergy = ((const G4Ions*)(&theParentNucleus))->GetExcitationEnergy();
764  G4Ions::G4FloatLevelBase floatingLevel =
765  ((const G4Ions*)(&theParentNucleus))->GetFloatLevelBase();
766 
767 #ifdef G4MULTITHREADED
768  G4AutoLock lk(&G4RadioactiveDecay::radioactiveDecayMutex);
769 
770  G4String key = theParentNucleus.GetParticleName();
771  DecayTableMap::iterator master_table_ptr = master_dkmap->find(key);
772 
773  if (master_table_ptr != master_dkmap->end() ) { // If table is there
774  return master_table_ptr->second;
775  }
776 #endif
777 
778  //Check if data have been provided by the user
779  G4String file = theUserRadioactiveDataFiles[1000*A+Z];
780 
781  if (file == "") {
782  if (!getenv("G4RADIOACTIVEDATA") ) {
783  G4cout << "Please setenv G4RADIOACTIVEDATA to point to the radioactive decay data files."
784  << G4endl;
785  throw G4HadronicException(__FILE__, __LINE__, " Please setenv G4RADIOACTIVEDATA to point to the radioactive decay data files.");
786  }
787  G4String dirName = getenv("G4RADIOACTIVEDATA");
788 
789  std::ostringstream os;
790  os << dirName << "/z" << Z << ".a" << A << '\0';
791  file = os.str();
792  }
793 
794  G4DecayTable* theDecayTable = new G4DecayTable();
795  G4bool found(false); // True if energy level matches one in table
796 
797  std::ifstream DecaySchemeFile;
798  DecaySchemeFile.open(file);
799 
800  if (DecaySchemeFile.good()) {
801  // Initialize variables used for reading in radioactive decay data
802  G4bool floatMatch(false);
803  const G4int nMode = 9;
804  G4double modeTotalBR[nMode] = {0.0};
805  G4double modeSumBR[nMode];
806  for (G4int i = 0; i < nMode; i++) {
807  modeSumBR[i] = 0.0;
808  }
809 
810  char inputChars[120]={' '};
811  G4String inputLine;
812  G4String recordType("");
813  G4String floatingFlag("");
814  G4String daughterFloatFlag("");
815  G4Ions::G4FloatLevelBase daughterFloatLevel;
816  G4RadioactiveDecayMode theDecayMode;
817  G4double decayModeTotal(0.0);
818  G4double parentExcitation(0.0);
819  G4double a(0.0);
820  G4double b(0.0);
821  G4double c(0.0);
822  G4double dummy(0.0);
823  G4BetaDecayType betaType(allowed);
824 
825  // Loop through each data file record until you identify the decay
826  // data relating to the nuclide of concern.
827 
828  G4bool complete(false); // bool insures only one set of values read for any
829  // given parent energy level
830  G4int loop = 0;
832  ed << " While count exceeded " << G4endl;
833 
834  while (!complete && !DecaySchemeFile.getline(inputChars, 120).eof()) { /* Loop checking, 01.09.2015, D.Wright */
835  loop++;
836  if (loop > 100000) {
837  G4Exception("G4RadioactiveDecay::LoadDecayTable()", "HAD_RDM_100", JustWarning, ed);
838  break;
839  }
840 
841  inputLine = inputChars;
842  inputLine = inputLine.strip(1);
843  if (inputChars[0] != '#' && inputLine.length() != 0) {
844  std::istringstream tmpStream(inputLine);
845 
846  if (inputChars[0] == 'P') {
847  // Nucleus is a parent type. Check excitation level to see if it
848  // matches that of theParentNucleus
849  tmpStream >> recordType >> parentExcitation >> floatingFlag >> dummy;
850  // "dummy" takes the place of half-life
851  // Now read in from ENSDFSTATE in particle category
852 
853  if (found) {
854  complete = true;
855  } else {
856  // Take first level which matches excitation energy regardless of floating level
857  found = (std::abs(parentExcitation*keV - levelEnergy) < levelTolerance);
858  if (floatingLevel != noFloat) {
859  // If floating level specificed, require match of both energy and floating level
860  floatMatch = (floatingLevel == G4Ions::FloatLevelBase(floatingFlag.back()) );
861  if (!floatMatch) found = false;
862  }
863  }
864 
865  } else if (found) {
866  // The right part of the radioactive decay data file has been found. Search
867  // through it to determine the mode of decay of the subsequent records.
868 
869  // Store for later the total decay probability for each decay mode
870  if (inputLine.length() < 72) {
871  tmpStream >> theDecayMode >> dummy >> decayModeTotal;
872  switch (theDecayMode) {
873  case IT:
874  {
875  G4ITDecay* anITChannel = new G4ITDecay(&theParentNucleus, decayModeTotal,
876  0.0, 0.0);
877  anITChannel->SetHLThreshold(halflifethreshold);
878  anITChannel->SetARM(applyARM);
879  theDecayTable->Insert(anITChannel);
880 // anITChannel->DumpNuclearInfo();
881  }
882  break;
883  case BetaMinus:
884  modeTotalBR[1] = decayModeTotal; break;
885  case BetaPlus:
886  modeTotalBR[2] = decayModeTotal; break;
887  case KshellEC:
888  modeTotalBR[3] = decayModeTotal; break;
889  case LshellEC:
890  modeTotalBR[4] = decayModeTotal; break;
891  case MshellEC:
892  modeTotalBR[5] = decayModeTotal; break;
893  case Alpha:
894  modeTotalBR[6] = decayModeTotal; break;
895  case Proton:
896  modeTotalBR[7] = decayModeTotal; break;
897  case Neutron:
898  modeTotalBR[8] = decayModeTotal; break;
899  case BDProton:
900  break;
901  case BDNeutron:
902  break;
903  case Beta2Minus:
904  break;
905  case Beta2Plus:
906  break;
907  case Proton2:
908  break;
909  case Neutron2:
910  break;
911  case SpFission:
912  break;
913  case RDM_ERROR:
914 
915  default:
916  G4Exception("G4RadioactiveDecay::LoadDecayTable()", "HAD_RDM_000",
917  FatalException, "Selected decay mode does not exist");
918  } // switch
919 
920  } else {
921  if (inputLine.length() < 84) {
922  tmpStream >> theDecayMode >> a >> daughterFloatFlag >> b >> c;
923  betaType = allowed;
924  } else {
925  tmpStream >> theDecayMode >> a >> daughterFloatFlag >> b >> c >> betaType;
926  }
927 
928  // Allowed transitions are the default. Forbidden transitions are
929  // indicated in the last column.
930  a /= 1000.;
931  c /= 1000.;
932  daughterFloatLevel = G4Ions::FloatLevelBase(daughterFloatFlag.back());
933 
934  switch (theDecayMode) {
935  case BetaMinus:
936  {
937  G4BetaMinusDecay* aBetaMinusChannel =
938  new G4BetaMinusDecay(&theParentNucleus, b, c*MeV, a*MeV,
939  daughterFloatLevel, betaType);
940 // aBetaMinusChannel->DumpNuclearInfo();
941  aBetaMinusChannel->SetHLThreshold(halflifethreshold);
942  theDecayTable->Insert(aBetaMinusChannel);
943  modeSumBR[1] += b;
944  }
945  break;
946 
947  case BetaPlus:
948  {
949  G4BetaPlusDecay* aBetaPlusChannel =
950  new G4BetaPlusDecay(&theParentNucleus, b, c*MeV, a*MeV,
951  daughterFloatLevel, betaType);
952 // aBetaPlusChannel->DumpNuclearInfo();
953  aBetaPlusChannel->SetHLThreshold(halflifethreshold);
954  theDecayTable->Insert(aBetaPlusChannel);
955  modeSumBR[2] += b;
956  }
957  break;
958 
959  case KshellEC: // K-shell electron capture
960  {
961  G4ECDecay* aKECChannel =
962  new G4ECDecay(&theParentNucleus, b, c*MeV, a*MeV,
963  daughterFloatLevel, KshellEC);
964 // aKECChannel->DumpNuclearInfo();
965  aKECChannel->SetHLThreshold(halflifethreshold);
966  aKECChannel->SetARM(applyARM);
967  theDecayTable->Insert(aKECChannel);
968  modeSumBR[3] += b;
969  }
970  break;
971 
972  case LshellEC: // L-shell electron capture
973  {
974  G4ECDecay* aLECChannel =
975  new G4ECDecay(&theParentNucleus, b, c*MeV, a*MeV,
976  daughterFloatLevel, LshellEC);
977 // aLECChannel->DumpNuclearInfo();
978  aLECChannel->SetHLThreshold(halflifethreshold);
979  aLECChannel->SetARM(applyARM);
980  theDecayTable->Insert(aLECChannel);
981  modeSumBR[4] += b;
982  }
983  break;
984 
985  case MshellEC: // M-shell electron capture
986  {
987  G4ECDecay* aMECChannel =
988  new G4ECDecay(&theParentNucleus, b, c*MeV, a*MeV,
989  daughterFloatLevel, MshellEC);
990 // aMECChannel->DumpNuclearInfo();
991  aMECChannel->SetHLThreshold(halflifethreshold);
992  aMECChannel->SetARM(applyARM);
993  theDecayTable->Insert(aMECChannel);
994  modeSumBR[5] += b;
995  }
996  break;
997 
998  case Alpha:
999  {
1000  G4AlphaDecay* anAlphaChannel =
1001  new G4AlphaDecay(&theParentNucleus, b, c*MeV, a*MeV,
1002  daughterFloatLevel);
1003 // anAlphaChannel->DumpNuclearInfo();
1004  anAlphaChannel->SetHLThreshold(halflifethreshold);
1005  theDecayTable->Insert(anAlphaChannel);
1006  modeSumBR[6] += b;
1007  }
1008  break;
1009 
1010  case Proton:
1011  {
1012  G4ProtonDecay* aProtonChannel =
1013  new G4ProtonDecay(&theParentNucleus, b, c*MeV, a*MeV,
1014  daughterFloatLevel);
1015 // aProtonChannel->DumpNuclearInfo();
1016  aProtonChannel->SetHLThreshold(halflifethreshold);
1017  theDecayTable->Insert(aProtonChannel);
1018  modeSumBR[7] += b;
1019  }
1020  break;
1021 
1022  case Neutron:
1023  {
1024  G4NeutronDecay* aNeutronChannel =
1025  new G4NeutronDecay(&theParentNucleus, b, c*MeV, a*MeV,
1026  daughterFloatLevel);
1027 // aNeutronChannel->DumpNuclearInfo();
1028  aNeutronChannel->SetHLThreshold(halflifethreshold);
1029  theDecayTable->Insert(aNeutronChannel);
1030  modeSumBR[8] += b;
1031  }
1032  break;
1033 
1034  case BDProton:
1035  // Not yet implemented
1036  // G4cout << " beta-delayed proton decay, a = " << a << ", b = " << b << ", c = " << c << G4endl;
1037  break;
1038  case BDNeutron:
1039  // Not yet implemented
1040  // G4cout << " beta-delayed neutron decay, a = " << a << ", b = " << b << ", c = " << c << G4endl;
1041  break;
1042  case Beta2Minus:
1043  // Not yet implemented
1044  // G4cout << " Double beta- decay, a = " << a << ", b = " << b << ", c = " << c << G4endl;
1045  break;
1046  case Beta2Plus:
1047  // Not yet implemented
1048  // G4cout << " Double beta+ decay, a = " << a << ", b = " << b << ", c = " << c << G4endl;
1049  break;
1050  case Proton2:
1051  // Not yet implemented
1052  // G4cout << " Double proton decay, a = " << a << ", b = " << b << ", c = " << c << G4endl;
1053  break;
1054  case Neutron2:
1055  // Not yet implemented
1056  // G4cout << " Double beta- decay, a = " << a << ", b = " << b << ", c = " << c << G4endl;
1057  break;
1058  case SpFission:
1059  // Not yet implemented
1060  //G4cout<<"Sp fission channel"<<a<<'\t'<<b<<'\t'<<c<<std::endl;
1061  break;
1062  case RDM_ERROR:
1063 
1064  default:
1065  G4Exception("G4RadioactiveDecay::LoadDecayTable()", "HAD_RDM_000",
1066  FatalException, "Selected decay mode does not exist");
1067  } // switch
1068  } // line < 72
1069  } // if char == P
1070  } // if char != #
1071  } // While
1072 
1073  // Go through the decay table and make sure that the branching ratios are
1074  // correctly normalised.
1075 
1076  G4VDecayChannel* theChannel = 0;
1077  G4NuclearDecay* theNuclearDecayChannel = 0;
1078  G4String mode = "";
1079 
1080  G4double theBR = 0.0;
1081  for (G4int i = 0; i < theDecayTable->entries(); i++) {
1082  theChannel = theDecayTable->GetDecayChannel(i);
1083  theNuclearDecayChannel = static_cast<G4NuclearDecay*>(theChannel);
1084  theDecayMode = theNuclearDecayChannel->GetDecayMode();
1085 
1086  if (theDecayMode != IT) {
1087  theBR = theChannel->GetBR();
1088  theChannel->SetBR(theBR*modeTotalBR[theDecayMode]/modeSumBR[theDecayMode]);
1089  }
1090  }
1091  } // decay file exists
1092 
1093  DecaySchemeFile.close();
1094 
1095  if (!found && levelEnergy > 0) {
1096  // Case where IT cascade for excited isotopes has no entries in RDM database
1097  // Decay mode is isomeric transition.
1098  G4ITDecay* anITChannel = new G4ITDecay(&theParentNucleus, 1.0, 0.0, 0.0);
1099  anITChannel->SetHLThreshold(halflifethreshold);
1100  anITChannel->SetARM(applyARM);
1101  theDecayTable->Insert(anITChannel);
1102  }
1103 
1104  if (theDecayTable && GetVerboseLevel() > 1) {
1105  theDecayTable->DumpInfo();
1106  }
1107 
1108 #ifdef G4MULTITHREADED
1109  //(*master_dkmap)[key] = theDecayTable; // store in master library
1110 #endif
1111  return theDecayTable;
1112 }
G4double GetBR() const
void SetBR(G4double value)
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
void SetHLThreshold(G4double HLT)
void SetARM(G4bool onoff)
Definition: G4ITDecay.hh:59
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
G4String strip(G4int strip_Type=trailing, char c=' ')
G4RadioactiveDecayMode GetDecayMode()
G4VDecayChannel * GetDecayChannel(G4int index) const
G4int GetVerboseLevel() const
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
G4int entries() const
tuple b
Definition: test.py:12
static G4Ions::G4FloatLevelBase FloatLevelBase(char flbChar)
Definition: G4Ions.hh:189
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
Definition: G4Ions.hh:51
bool G4bool
Definition: G4Types.hh:79
void DumpInfo() const
G4RadioactiveDecayMode
G4BetaDecayType
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
void Insert(G4VDecayChannel *aChannel)
Definition: G4DecayTable.cc:60
void SetARM(G4bool onoff)
Definition: G4ECDecay.hh:56
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
double G4double
Definition: G4Types.hh:76
tuple c
Definition: test.py:13
static constexpr double keV
Definition: G4SIunits.hh:216
#define noFloat
Definition: G4Ions.hh:118
G4FloatLevelBase
Definition: G4Ions.hh:95

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RadioactiveDecay::SelectAllVolumes ( )

Definition at line 320 of file G4RadioactiveDecay.cc.

321 {
322  G4LogicalVolumeStore* theLogicalVolumes;
323  G4LogicalVolume* volume;
324  theLogicalVolumes = G4LogicalVolumeStore::GetInstance();
325  ValidVolumes.clear();
326 #ifdef G4VERBOSE
327  if (GetVerboseLevel()>0)
328  G4cout << " RDM Applies to all Volumes" << G4endl;
329 #endif
330  for (size_t i = 0; i < theLogicalVolumes->size(); i++){
331  volume = (*theLogicalVolumes)[i];
332  ValidVolumes.push_back(volume->GetName());
333 #ifdef G4VERBOSE
334  if (GetVerboseLevel()>0)
335  G4cout << " RDM Applies to Volume " << volume->GetName() << G4endl;
336 #endif
337  }
338  std::sort(ValidVolumes.begin(), ValidVolumes.end());
339  // sort needed in order to allow binary_search
340  isAllVolumesMode=true;
341 }
G4int GetVerboseLevel() const
G4GLOB_DLL std::ostream G4cout
static G4LogicalVolumeStore * GetInstance()
#define G4endl
Definition: G4ios.hh:61
const G4String & GetName() const

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RadioactiveDecay::SelectAVolume ( const G4String  aVolume)

Definition at line 266 of file G4RadioactiveDecay.cc.

267 {
268  G4LogicalVolumeStore* theLogicalVolumes;
269  G4LogicalVolume* volume;
270  theLogicalVolumes = G4LogicalVolumeStore::GetInstance();
271  for (size_t i = 0; i < theLogicalVolumes->size(); i++) {
272  volume=(*theLogicalVolumes)[i];
273  if (volume->GetName() == aVolume) {
274  ValidVolumes.push_back(aVolume);
275  std::sort(ValidVolumes.begin(), ValidVolumes.end());
276  // sort need for performing binary_search
277 #ifdef G4VERBOSE
278  if (GetVerboseLevel()>0)
279  G4cout << " RDM Applies to : " << aVolume << G4endl;
280 #endif
281  } else if(i == theLogicalVolumes->size()) {
282  G4cerr << "SelectAVolume: "<< aVolume
283  << " is not a valid logical volume name" << G4endl;
284  }
285  }
286 }
G4int GetVerboseLevel() const
G4GLOB_DLL std::ostream G4cout
static G4LogicalVolumeStore * GetInstance()
#define G4endl
Definition: G4ios.hh:61
const G4String & GetName() const
G4GLOB_DLL std::ostream G4cerr

Here is the call graph for this function:

void G4RadioactiveDecay::SetAnalogueMonteCarlo ( G4bool  r)
inline

Definition at line 183 of file G4RadioactiveDecay.hh.

183  {
184  AnalogueMC = r;
185  if (!AnalogueMC) halflifethreshold = 1e-6*CLHEP::s;
186  }
static constexpr double s

Here is the caller graph for this function:

void G4RadioactiveDecay::SetARM ( G4bool  arm)
inline

Definition at line 128 of file G4RadioactiveDecay.hh.

128 {applyARM = arm;}

Here is the caller graph for this function:

void G4RadioactiveDecay::SetBRBias ( G4bool  r)
inline

Definition at line 197 of file G4RadioactiveDecay.hh.

197  {
198  BRBias = r;
200  }
void SetAnalogueMonteCarlo(G4bool r)

Here is the call graph for this function:

void G4RadioactiveDecay::SetDecayBias ( G4String  filename)

Definition at line 1511 of file G4RadioactiveDecay.cc.

1512 {
1513 
1514  std::ifstream infile(filename, std::ios::in);
1515  if (!infile) G4Exception("G4RadioactiveDecay::SetDecayBias()", "HAD_RDM_003",
1516  FatalException, "Unable to open bias data file" );
1517 
1518  G4double bin, flux;
1519  G4int dWindows = 0;
1520  G4int i ;
1521 
1522  theRadioactivityTables.clear();
1523 
1524  NDecayBin = -1;
1525 
1526  G4int loop = 0;
1528  ed << " While count exceeded " << G4endl;
1529 
1530  while (infile >> bin >> flux ) { /* Loop checking, 01.09.2015, D.Wright */
1531  NDecayBin++;
1532  loop++;
1533  if (loop > 10000) {
1534  G4Exception("G4RadioactiveDecay::SetDecayBias()", "HAD_RDM_100", JustWarning, ed);
1535  break;
1536  }
1537 
1538  if (NDecayBin > 99) {
1539  G4Exception("G4RadioactiveDecay::SetDecayBias()", "HAD_RDM_004",
1540  FatalException, "Input bias file too big (>100 rows)" );
1541  } else {
1542  DBin[NDecayBin] = bin * s;
1543  DProfile[NDecayBin] = flux;
1544  if (flux > 0.) {
1545  decayWindows[NDecayBin] = dWindows;
1546  dWindows++;
1547  G4RadioactivityTable *rTable = new G4RadioactivityTable() ;
1548  theRadioactivityTables.push_back(rTable);
1549  }
1550  }
1551  }
1552  for ( i = 1; i<= NDecayBin; i++) DProfile[i] += DProfile[i-1];
1553  for ( i = 0; i<= NDecayBin; i++) DProfile[i] /= DProfile[NDecayBin];
1554  // converted to accumulated probabilities
1555 
1557  infile.close();
1558 
1559 #ifdef G4VERBOSE
1560  if (GetVerboseLevel() > 1)
1561  G4cout <<" Decay Bias Profile Nbin = " << NDecayBin <<G4endl;
1562 #endif
1563 }
tuple bin
Definition: plottest35.py:22
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4int GetVerboseLevel() const
int G4int
Definition: G4Types.hh:78
const XML_Char * s
Definition: expat.h:262
G4GLOB_DLL std::ostream G4cout
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
void SetAnalogueMonteCarlo(G4bool r)

Here is the call graph for this function:

void G4RadioactiveDecay::SetDecayCollimation ( const G4ThreeVector theDir,
G4double  halfAngle = 0.*CLHEP::deg 
)
inline

Definition at line 226 of file G4RadioactiveDecay.hh.

227  {
228  SetDecayDirection(theDir);
229  SetDecayHalfAngle(halfAngle);
230  }
void SetDecayDirection(const G4ThreeVector &theDir)
void SetDecayHalfAngle(G4double halfAngle=0.*CLHEP::deg)

Here is the call graph for this function:

void G4RadioactiveDecay::SetDecayDirection ( const G4ThreeVector theDir)
inline

Definition at line 212 of file G4RadioactiveDecay.hh.

212  {
213  forceDecayDirection = theDir.unit();
214  }
Hep3Vector unit() const

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RadioactiveDecay::SetDecayHalfAngle ( G4double  halfAngle = 0.*CLHEP::deg)
inline

Definition at line 220 of file G4RadioactiveDecay.hh.

220  {
221  forceDecayHalfAngle = std::min(std::max(0.*CLHEP::deg,halfAngle),180.*CLHEP::deg);
222  }
static constexpr double deg
T max(const T t1, const T t2)
brief Return the largest of the two arguments
T min(const T t1, const T t2)
brief Return the smallest of the two arguments

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RadioactiveDecay::SetDecayRate ( G4int  theZ,
G4int  theA,
G4double  theE,
G4int  theG,
std::vector< G4double theCoefficients,
std::vector< G4double theTaos 
)

Definition at line 1130 of file G4RadioactiveDecay.cc.

1134 {
1135  //fill the decay rate vector
1136  theDecayRate.SetZ(theZ);
1137  theDecayRate.SetA(theA);
1138  theDecayRate.SetE(theE);
1139  theDecayRate.SetGeneration(theG);
1140  theDecayRate.SetDecayRateC(theCoefficients);
1141  theDecayRate.SetTaos(theTaos);
1142 }
void SetTaos(std::vector< G4double > value)
void SetDecayRateC(std::vector< G4double > value)

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RadioactiveDecay::SetFBeta ( G4bool  r)
inline

Definition at line 190 of file G4RadioactiveDecay.hh.

190 { FBeta = r; }
void G4RadioactiveDecay::SetHLThreshold ( G4double  hl)
inline

Definition at line 122 of file G4RadioactiveDecay.hh.

122 {halflifethreshold = hl;}
void G4RadioactiveDecay::SetICM ( G4bool  icm)
inline

Definition at line 125 of file G4RadioactiveDecay.hh.

125 {applyICM = icm;}
void G4RadioactiveDecay::SetNucleusLimits ( G4NucleusLimits  theNucleusLimits1)
inline

Definition at line 173 of file G4RadioactiveDecay.hh.

174  {theNucleusLimits = theNucleusLimits1 ;}
void G4RadioactiveDecay::SetSourceTimeProfile ( G4String  filename)

Definition at line 1461 of file G4RadioactiveDecay.cc.

1462 {
1463  std::ifstream infile ( filename, std::ios::in );
1464  if (!infile) {
1466  ed << " Could not open file " << filename << G4endl;
1467  G4Exception("G4RadioactiveDecay::SetSourceTimeProfile()", "HAD_RDM_001",
1468  FatalException, ed);
1469  }
1470 
1471  G4double bin, flux;
1472  NSourceBin = -1;
1473 
1474  G4int loop = 0;
1476  ed << " While count exceeded " << G4endl;
1477 
1478  while (infile >> bin >> flux ) { /* Loop checking, 01.09.2015, D.Wright */
1479  loop++;
1480  if (loop > 10000) {
1481  G4Exception("G4RadioactiveDecay::SetSourceTimeProfile()", "HAD_RDM_100", JustWarning, ed);
1482  break;
1483  }
1484 
1485  NSourceBin++;
1486  if (NSourceBin > 99) {
1487  G4Exception("G4RadioactiveDecay::SetSourceTimeProfile()", "HAD_RDM_002",
1488  FatalException, "Input source time file too big (>100 rows)");
1489 
1490  } else {
1491  SBin[NSourceBin] = bin * s;
1492  SProfile[NSourceBin] = flux;
1493  }
1494  }
1496  infile.close();
1497 
1498 #ifdef G4VERBOSE
1499  if (GetVerboseLevel() > 1)
1500  G4cout <<" Source Timeprofile Nbin = " << NSourceBin <<G4endl;
1501 #endif
1502 }
tuple bin
Definition: plottest35.py:22
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4int GetVerboseLevel() const
int G4int
Definition: G4Types.hh:78
const XML_Char * s
Definition: expat.h:262
G4GLOB_DLL std::ostream G4cout
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
void SetAnalogueMonteCarlo(G4bool r)

Here is the call graph for this function:

void G4RadioactiveDecay::SetSplitNuclei ( G4int  r)
inline

Definition at line 203 of file G4RadioactiveDecay.hh.

203  {
204  NSplit = r;
206  }
void SetAnalogueMonteCarlo(G4bool r)

Here is the call graph for this function:

void G4RadioactiveDecay::SetVerboseLevel ( G4int  value)
inline

Definition at line 167 of file G4RadioactiveDecay.hh.

167 {verboseLevel = value;}
const XML_Char int const XML_Char * value
Definition: expat.h:331

The documentation for this class was generated from the following files: