#include <G4DNAMolecularDecay.hh>
Inheritance diagram for G4DNAMolecularDecay:
Public Member Functions | |
G4DNAMolecularDecay (const G4String &processName="DNAMolecularDecay", G4ProcessType type=fDecay) | |
virtual | ~G4DNAMolecularDecay () |
virtual G4bool | IsApplicable (const G4ParticleDefinition &) |
G4double | AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition) |
G4VParticleChange * | AtRestDoIt (const G4Track &track, const G4Step &step) |
void | SetVerbose (G4int) |
void | SetDecayDisplacer (const G4ParticleDefinition *, G4VMolecularDecayDisplacer *) |
G4VMolecularDecayDisplacer * | GetDecayDisplacer (const G4ParticleDefinition *) |
Protected Member Functions | |
virtual G4VParticleChange * | DecayIt (const G4Track &, const G4Step &) |
virtual G4double | GetMeanLifeTime (const G4Track &, G4ForceCondition *) |
Definition at line 56 of file G4DNAMolecularDecay.hh.
G4DNAMolecularDecay::G4DNAMolecularDecay | ( | const G4String & | processName = "DNAMolecularDecay" , |
|
G4ProcessType | type = fDecay | |||
) |
Definition at line 48 of file G4DNAMolecularDecay.cc.
References G4VProcess::aParticleChange, G4VProcess::enableAlongStepDoIt, G4VProcess::enableAtRestDoIt, G4VProcess::enablePostStepDoIt, G4cout, G4DNAMolecularDecay(), G4endl, G4VProcess::pParticleChange, G4VProcess::SetProcessSubType(), and G4VProcess::verboseLevel.
Referenced by G4DNAMolecularDecay().
00049 : G4VITRestProcess(processName, type) 00050 { 00051 // set Process Sub Type 00052 SetProcessSubType(59); // DNA sub-type 00053 enableAlongStepDoIt = false; 00054 enablePostStepDoIt = false; 00055 enableAtRestDoIt=true; 00056 00057 fVerbose = 0 ; 00058 00059 #ifdef G4VERBOSE 00060 if (verboseLevel>1) 00061 { 00062 G4cout << "G4MolecularDecayProcess constructor " << " Name:" << processName << G4endl; 00063 } 00064 #endif 00065 00066 pParticleChange = &aParticleChange; 00067 00068 fDecayAtFixedTime = true ; 00069 }
G4DNAMolecularDecay::~G4DNAMolecularDecay | ( | ) | [virtual] |
Definition at line 71 of file G4DNAMolecularDecay.cc.
00072 { 00073 DecayDisplacementMap::iterator it = fDecayDisplacementMap.begin(); 00074 00075 for( ; it != fDecayDisplacementMap.end() ; it++) 00076 { 00077 if(it->second) 00078 { 00079 delete it->second; 00080 it->second = 0; 00081 } 00082 } 00083 fDecayDisplacementMap.clear(); 00084 }
G4VParticleChange * G4DNAMolecularDecay::AtRestDoIt | ( | const G4Track & | track, | |
const G4Step & | step | |||
) | [inline, virtual] |
Reimplemented from G4VITRestProcess.
Definition at line 120 of file G4DNAMolecularDecay.hh.
References G4VITProcess::ClearInteractionTimeLeft(), G4VITProcess::ClearNumberOfInteractionLengthLeft(), and DecayIt().
00124 { 00125 ClearNumberOfInteractionLengthLeft(); 00126 ClearInteractionTimeLeft(); 00127 return DecayIt(track, step); 00128 }
G4double G4DNAMolecularDecay::AtRestGetPhysicalInteractionLength | ( | const G4Track & | track, | |
G4ForceCondition * | condition | |||
) | [inline, virtual] |
Reimplemented from G4VITRestProcess.
Definition at line 108 of file G4DNAMolecularDecay.hh.
References G4VITRestProcess::AtRestGetPhysicalInteractionLength(), and GetMeanLifeTime().
00111 { 00112 if(fDecayAtFixedTime) 00113 { 00114 return GetMeanLifeTime(track, condition); 00115 } 00116 00117 return G4VITRestProcess::AtRestGetPhysicalInteractionLength(track, condition); 00118 }
G4VParticleChange * G4DNAMolecularDecay::DecayIt | ( | const G4Track & | , | |
const G4Step & | ||||
) | [protected, virtual] |
Definition at line 121 of file G4DNAMolecularDecay.cc.
References G4VProcess::aParticleChange, G4ElectronOccupancy::DumpInfo(), fAlive, FatalErrorInArgument, FatalException, fStopAndKill, G4cout, G4endl, G4Exception(), G4UniformRand, G4Molecule::GetDefinition(), G4Molecule::GetElectronOccupancy(), G4MolecularDecayChannel::GetEnergy(), GetMolecule(), G4MolecularDecayChannel::GetName(), G4Molecule::GetName(), G4MolecularDecayChannel::GetNbProducts(), G4Track::GetPosition(), G4MolecularDecayChannel::GetProbability(), G4MolecularDecayChannel::GetProduct(), G4Track::GetTrackID(), G4ParticleChange::Initialize(), G4ITManager< T >::Instance(), G4VParticleChange::ProposeLocalEnergyDeposit(), G4VParticleChange::ProposeTrackStatus(), and G4VParticleChange::SetNumberOfSecondaries().
Referenced by AtRestDoIt().
00124 { 00125 // DEBUG 00126 // G4cout << "Is calling G4MolecularDecayProcess::DecayIt" << G4endl; 00127 00128 aParticleChange.Initialize(track); 00129 const G4Molecule * theMotherMolecule = GetMolecule(track); 00130 const G4MoleculeDefinition* moleculeDefinition = theMotherMolecule->GetDefinition(); 00131 00132 // DEBUG 00133 // G4cout <<"Calling G4MolecularDecayProcess::DecayIt"<<G4endl; 00134 // G4cout << "The mother molecule state : " << G4endl; 00135 // theMotherMolecule -> PrintState(); 00136 00137 if(moleculeDefinition-> GetDecayTable()) 00138 { 00139 const vector<const G4MolecularDecayChannel*>* DecayVector = 00140 (theMotherMolecule -> GetDecayChannel()); 00141 00142 if(DecayVector == 0) 00143 { 00144 G4ExceptionDescription exceptionDescription; 00145 theMotherMolecule->GetElectronOccupancy()->DumpInfo(); 00146 exceptionDescription << "No decay channel was found for the molecule : " << theMotherMolecule-> GetName() << G4endl; 00147 G4Exception("G4DNAMolecularDecay::DecayIt", "G4DNAMolecularDecay::NoDecayChannel",FatalException,exceptionDescription); 00148 return &aParticleChange; 00149 } 00150 00151 G4int DecayVectorSize = DecayVector-> size(); 00152 // DEBUG 00153 // G4cout<< "Number of decay channels : " << DecayVectorSize<<G4endl; 00154 G4double RdmValue = G4UniformRand(); 00155 00156 const G4MolecularDecayChannel* decayChannel(0); 00157 G4int i=0; 00158 do 00159 { 00160 decayChannel = (*DecayVector)[i]; 00161 if(RdmValue < decayChannel->GetProbability()) break; 00162 RdmValue -= decayChannel->GetProbability(); 00163 i++; 00164 } 00165 while(i< DecayVectorSize); 00166 00167 // DEBUG 00168 // G4cout<< "Selected Decay channel : " << decayChannel->GetName()<<G4endl; 00169 00170 G4double decayEnergy = decayChannel->GetEnergy(); 00171 G4int nbProducts = decayChannel->GetNbProducts(); 00172 00173 if(decayEnergy) 00174 { 00175 // DEBUG 00176 // G4cout<<"Deposit energy :" <<decayChannel->GetEnergy()/eV << " eV"<<G4endl; 00177 00178 aParticleChange.ProposeLocalEnergyDeposit(decayChannel->GetEnergy()); 00179 } 00180 00181 if(nbProducts) 00182 { 00183 00184 // DEBUG 00185 // G4cout<<"Number of products :" <<nbProducts<<G4endl; 00186 00187 vector<G4ThreeVector> ProductsDisplacement(nbProducts); 00188 G4ThreeVector theMotherMoleculeDisplacement; 00189 00190 DecayDisplacementMap::iterator it = fDecayDisplacementMap.find(moleculeDefinition); 00191 00192 if(it!=fDecayDisplacementMap.end()) 00193 { 00194 G4VMolecularDecayDisplacer* displacer = it->second; 00195 ProductsDisplacement = displacer->GetProductsDisplacement(decayChannel); 00196 theMotherMoleculeDisplacement = displacer-> GetMotherMoleculeDisplacement(decayChannel); 00197 } 00198 else 00199 { 00200 G4ExceptionDescription errMsg; 00201 errMsg << "No G4MolecularDecayProcess::theDecayDisplacementMap[" 00202 << theMotherMolecule->GetName() +"]" ; 00203 G4Exception("G4MolecularDecayProcess::DecayIt","DNAMolecularDecay001",FatalErrorInArgument, errMsg); 00204 } 00205 00206 aParticleChange.SetNumberOfSecondaries(nbProducts); 00207 00208 #ifdef G4VERBOSE 00209 if(fVerbose) 00210 { 00211 G4cout<<"Decay Process : " 00212 << theMotherMolecule->GetName() 00213 << " (trackID :" << track.GetTrackID() << ") " 00214 << decayChannel->GetName() 00215 << G4endl; 00216 } 00217 #endif 00218 00219 for (G4int j=0; j<nbProducts ; j++) 00220 { 00221 G4Molecule* product = new G4Molecule(*decayChannel->GetProduct(j)); 00222 00223 // create a new track object 00224 // Be carefull as this processes is dedicated to be used in the DNA module 00225 // The molecular decay will happen one picosecond after the start of the simulation. 00226 // This may be seen as a bug and will be hopefully improve in the next releases 00227 G4Track* secondary = product->BuildTrack(picosecond,track.GetPosition() 00228 + theMotherMoleculeDisplacement + ProductsDisplacement[j]); 00229 00230 secondary-> SetTrackStatus(fAlive); 00231 #ifdef G4VERBOSE 00232 if(fVerbose) 00233 { 00234 G4cout<<"Product : "<< product->GetName()<<G4endl; 00235 } 00236 #endif 00237 // add the secondary track in the List 00238 aParticleChange.G4VParticleChange::AddSecondary(secondary); 00239 G4ITManager<G4Molecule>::Instance()->Push(secondary); 00240 } 00241 #ifdef G4VERBOSE 00242 if(fVerbose) 00243 G4cout<<"-------------"<<G4endl; 00244 #endif 00245 } 00246 // DEBUG 00247 // else if(decayEnergy) 00248 // { 00249 // G4cout << "No products for this channel" << G4endl ; 00250 // G4cout<<"-------------"<<G4endl; 00251 // } 00252 else if(!decayEnergy && !nbProducts) 00253 { 00254 G4ExceptionDescription errMsg; 00255 errMsg << "There is no products and no energy specified in the molecular decay channel"; 00256 G4Exception("G4MolecularDecayProcess::DecayIt","DNAMolecularDecay002",FatalErrorInArgument, errMsg); 00257 } 00258 } 00259 00260 aParticleChange.ProposeTrackStatus(fStopAndKill); 00261 00262 return &aParticleChange; 00263 }
G4VMolecularDecayDisplacer * G4DNAMolecularDecay::GetDecayDisplacer | ( | const G4ParticleDefinition * | ) |
G4double G4DNAMolecularDecay::GetMeanLifeTime | ( | const G4Track & | , | |
G4ForceCondition * | ||||
) | [protected, virtual] |
Implements G4VITRestProcess.
Definition at line 114 of file G4DNAMolecularDecay.cc.
References GetMolecule(), and G4Track::GetProperTime().
Referenced by AtRestGetPhysicalInteractionLength().
00116 { 00117 G4double output = GetMolecule(track)-> GetDecayTime() - track.GetProperTime() ; 00118 return (output > 0 ? output : 0 ); 00119 }
G4bool G4DNAMolecularDecay::IsApplicable | ( | const G4ParticleDefinition & | ) | [virtual] |
Reimplemented from G4VProcess.
Definition at line 94 of file G4DNAMolecularDecay.cc.
References G4cout, G4endl, G4ParticleDefinition::GetParticleName(), and G4ParticleDefinition::GetParticleType().
00095 { 00096 if(aParticleType.GetParticleType()=="Molecule") 00097 { 00098 #ifdef G4VERBOSE 00099 if(fVerbose>1) 00100 { 00101 G4cout<<"G4MolecularDecay::IsApplicable("; 00102 G4cout<<aParticleType.GetParticleName()<<","; 00103 G4cout<<aParticleType.GetParticleType()<<")"<<G4endl; 00104 } 00105 #endif 00106 return(true); 00107 } 00108 else 00109 { 00110 return false; 00111 } 00112 }
void G4DNAMolecularDecay::SetDecayDisplacer | ( | const G4ParticleDefinition * | , | |
G4VMolecularDecayDisplacer * | ||||
) |
void G4DNAMolecularDecay::SetVerbose | ( | G4int | ) | [inline] |