Geant4_10
G4DNAMolecularDissociation.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4DNAMolecularDissociation.cc 74551 2013-10-14 12:59:14Z gcosmo $
27 //
28 // Author: Mathieu Karamitros (kara (AT) cenbg . in2p3 . fr)
29 //
30 // WARNING : This class is released as a prototype.
31 // It might strongly evolve or even disapear in the next releases.
32 //
33 // History:
34 // -----------
35 // 10 Oct 2011 M.Karamitros created
36 //
37 // -------------------------------------------------------------------
38 
40 #include "G4SystemOfUnits.hh"
41 #include "G4Track.hh"
42 #include "G4Molecule.hh"
43 #include "G4ITManager.hh"
44 #include "G4ParticleChange.hh"
45 
46 using namespace std;
47 
49  G4ProcessType type) : G4VITRestProcess(processName, type)
50 {
51  // set Process Sub Type
52  SetProcessSubType(59); // DNA sub-type
53  enableAlongStepDoIt = false;
54  enablePostStepDoIt = false;
55  enableAtRestDoIt=true;
56 
57  fVerbose = 0 ;
58 
59 #ifdef G4VERBOSE
60  if (verboseLevel>1)
61  {
62  G4cout << "G4MolecularDecayProcess constructor " << " Name:" << processName << G4endl;
63  }
64 #endif
65 
67 
68  fDecayAtFixedTime = true ;
69 }
70 
72 {
73  DecayDisplacementMap::iterator it = fDecayDisplacementMap.begin();
74 
75  for( ; it != fDecayDisplacementMap.end() ; it++)
76  {
77  if(it->second)
78  {
79  delete it->second;
80  it->second = 0;
81  }
82  }
83  fDecayDisplacementMap.clear();
84 }
85 
86 G4DNAMolecularDissociation::G4DNAMolecularDissociation(const G4DNAMolecularDissociation &right) :
87  G4VITRestProcess(right)
88 {
89  fDecayAtFixedTime = right . fDecayAtFixedTime;
90  fDecayDisplacementMap = right.fDecayDisplacementMap;
91  fVerbose = right.fVerbose ;
92 }
93 
95 {
96  if(aParticleType.GetParticleType()=="Molecule")
97  {
98 #ifdef G4VERBOSE
99  if(fVerbose>1)
100  {
101  G4cout<<"G4MolecularDecay::IsApplicable(";
102  G4cout<<aParticleType.GetParticleName()<<",";
103  G4cout<<aParticleType.GetParticleType()<<")"<<G4endl;
104  }
105 #endif
106  return(true);
107  }
108  else
109  {
110  return false;
111  }
112 }
113 
116 {
117  G4double output = GetMolecule(track)-> GetDecayTime() - track.GetProperTime() ;
118  return (output > 0 ? output : 0 );
119 }
120 
122  const G4Track& track,
123  const G4Step& )
124 {
125  // DEBUG
126  // G4cout << "Is calling G4MolecularDecayProcess::DecayIt" << G4endl;
127 
129  const G4Molecule * theMotherMolecule = GetMolecule(track);
130  const G4MoleculeDefinition* moleculeDefinition = theMotherMolecule->GetDefinition();
131 
132  // DEBUG
133  // G4cout <<"Calling G4MolecularDecayProcess::DecayIt"<<G4endl;
134  // G4cout << "The mother molecule state : " << G4endl;
135  // theMotherMolecule -> PrintState();
136 
137  if(moleculeDefinition-> GetDecayTable())
138  {
139  const vector<const G4MolecularDecayChannel*>* DecayVector =
140  (theMotherMolecule -> GetDecayChannel());
141 
142  if(DecayVector == 0)
143  {
144  G4ExceptionDescription exceptionDescription;
145  theMotherMolecule->GetElectronOccupancy()->DumpInfo();
146  exceptionDescription << "No decay channel was found for the molecule : " << theMotherMolecule-> GetName() << G4endl;
147  G4Exception("G4DNAMolecularDissociation::DecayIt", "G4DNAMolecularDissociation::NoDecayChannel",FatalException,exceptionDescription);
148  return &aParticleChange;
149  }
150 
151  G4int DecayVectorSize = DecayVector-> size();
152  // DEBUG
153  // G4cout<< "Number of decay channels : " << DecayVectorSize<<G4endl;
154  G4double RdmValue = G4UniformRand();
155 
156  const G4MolecularDecayChannel* decayChannel(0);
157  G4int i=0;
158  do
159  {
160  decayChannel = (*DecayVector)[i];
161  if(RdmValue < decayChannel->GetProbability()) break;
162  RdmValue -= decayChannel->GetProbability();
163  i++;
164  }
165  while(i< DecayVectorSize);
166 
167  // DEBUG
168  // G4cout<< "Selected Decay channel : " << decayChannel->GetName()<<G4endl;
169 
170  G4double decayEnergy = decayChannel->GetEnergy();
171  G4int nbProducts = decayChannel->GetNbProducts();
172 
173  if(decayEnergy)
174  {
175  // DEBUG
176  // G4cout<<"Deposit energy :" <<decayChannel->GetEnergy()/eV << " eV"<<G4endl;
177 
179  }
180 
181  if(nbProducts)
182  {
183 
184  // DEBUG
185  // G4cout<<"Number of products :" <<nbProducts<<G4endl;
186 
187  vector<G4ThreeVector> ProductsDisplacement(nbProducts);
188  G4ThreeVector theMotherMoleculeDisplacement;
189 
190  DecayDisplacementMap::iterator it = fDecayDisplacementMap.find(moleculeDefinition);
191 
192  if(it!=fDecayDisplacementMap.end())
193  {
194  G4VMolecularDecayDisplacer* displacer = it->second;
195  ProductsDisplacement = displacer->GetProductsDisplacement(decayChannel);
196  theMotherMoleculeDisplacement = displacer-> GetMotherMoleculeDisplacement(decayChannel);
197  }
198  else
199  {
200  G4ExceptionDescription errMsg;
201  errMsg << "No G4MolecularDecayProcess::theDecayDisplacementMap["
202  << theMotherMolecule->GetName() +"]" ;
203  G4Exception("G4MolecularDecayProcess::DecayIt","DNAMolecularDecay001",FatalErrorInArgument, errMsg);
204  }
205 
207 
208 #ifdef G4VERBOSE
209  if(fVerbose)
210  {
211  G4cout<<"Decay Process : "
212  << theMotherMolecule->GetName()
213  << " (trackID :" << track.GetTrackID() << ") "
214  << decayChannel->GetName()
215  << G4endl;
216  }
217 #endif
218 
219  for (G4int j=0; j<nbProducts ; j++)
220  {
221  G4Molecule* product = new G4Molecule(*decayChannel->GetProduct(j));
222 
223  // create a new track object
224  // Be carefull as this processes is dedicated to be used in the DNA module
225  // The molecular decay will happen one picosecond after the start of the simulation.
226  // This may be seen as a bug and will be hopefully improve in the next releases
227  G4Track* secondary = product->BuildTrack(picosecond,track.GetPosition()
228  + theMotherMoleculeDisplacement + ProductsDisplacement[j]);
229 
230  secondary-> SetTrackStatus(fAlive);
231 #ifdef G4VERBOSE
232  if(fVerbose)
233  {
234  G4cout<<"Product : "<< product->GetName()<<G4endl;
235  }
236 #endif
237  // add the secondary track in the List
238  aParticleChange.G4VParticleChange::AddSecondary(secondary);
239  G4ITManager<G4Molecule>::Instance()->Push(secondary);
240  }
241 #ifdef G4VERBOSE
242  if(fVerbose)
243  G4cout<<"-------------"<<G4endl;
244 #endif
245  }
246  // DEBUG
247  // else if(decayEnergy)
248  // {
249  // G4cout << "No products for this channel" << G4endl ;
250  // G4cout<<"-------------"<<G4endl;
251  // }
252  else if(!decayEnergy && !nbProducts)
253  {
254  G4ExceptionDescription errMsg;
255  errMsg << "There is no products and no energy specified in the molecular decay channel";
256  G4Exception("G4MolecularDecayProcess::DecayIt","DNAMolecularDecay002",FatalErrorInArgument, errMsg);
257  }
258  }
259 
261 
262  return &aParticleChange;
263 }
264 
266 {
267  fDecayDisplacementMap[molDef] = aDisplacer;
268 }
269 
271 {
272  return fDecayDisplacementMap[molDef] ;
273 }
virtual G4double GetMeanLifeTime(const G4Track &, G4ForceCondition *)
G4int verboseLevel
Definition: G4VProcess.hh:368
G4VMolecularDecayDisplacer * GetDecayDisplacer(const G4ParticleDefinition *)
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4double GetProperTime() const
const G4ThreeVector & GetPosition() const
G4bool enableAtRestDoIt
Definition: G4VProcess.hh:350
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
Identical to G4VRestProcess with dependency from G4VITProcess.
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
const G4String & GetName() const
Definition: G4Molecule.cc:243
#define G4UniformRand()
Definition: Randomize.hh:87
G4GLOB_DLL std::ostream G4cout
static G4ITManager< T > * Instance()
bool G4bool
Definition: G4Types.hh:79
virtual G4VParticleChange * DecayIt(const G4Track &, const G4Step &)
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:432
const G4String & GetParticleType() const
virtual std::vector< G4ThreeVector > GetProductsDisplacement(const G4MolecularDecayChannel *) const =0
Definition: G4Step.hh:76
G4bool enablePostStepDoIt
Definition: G4VProcess.hh:352
G4int GetTrackID() const
G4Molecule * GetMolecule(const G4Track &track)
Definition: G4Molecule.cc:67
const G4MoleculeDefinition * GetDefinition() const
Definition: G4Molecule.cc:379
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4Track * BuildTrack(G4double globalTime, const G4ThreeVector &Position)
Definition: G4Molecule.cc:263
virtual void Initialize(const G4Track &)
int picosecond
Definition: hepunit.py:86
void SetNumberOfSecondaries(G4int totSecondaries)
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
G4ParticleChange aParticleChange
Definition: G4VProcess.hh:289
#define G4endl
Definition: G4ios.hh:61
G4bool enableAlongStepDoIt
Definition: G4VProcess.hh:351
void SetDecayDisplacer(const G4ParticleDefinition *, G4VMolecularDecayDisplacer *)
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
const G4ElectronOccupancy * GetElectronOccupancy() const
Definition: G4Molecule.cc:374
G4ForceCondition
const G4Molecule * GetProduct(int) const
virtual G4bool IsApplicable(const G4ParticleDefinition &)
G4DNAMolecularDissociation(const G4String &processName="DNAMolecularDecay", G4ProcessType type=fDecay)
G4ProcessType