Geant4  9.6.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4DNAChemistryManager.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: G4DNAChemistryManager.cc 64057 2012-10-30 15:04:49Z gcosmo $
27 //
28 // Author: Mathieu Karamitros (kara@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 
39 #include "G4DNAChemistryManager.hh"
40 #include "G4SystemOfUnits.hh"
41 #include "G4Molecule.hh"
42 #include "G4ITTrackHolder.hh"
43 #include "G4H2O.hh"
47 #include "G4Electron_aq.hh"
48 #include "G4ITManager.hh"
50 #include "G4MoleculeCounter.hh"
51 
52 using namespace std;
53 
54 auto_ptr<G4DNAChemistryManager> G4DNAChemistryManager::fInstance(0);
55 
56 G4DNAChemistryManager::G4DNAChemistryManager() :
57  fActiveChemistry(false)
58 {
59  fExcitationLevel = 0;
60  fIonisationLevel = 0;
61  fWriteFile = false;
62 }
63 
65 {
66  if(!fInstance.get()) fInstance = auto_ptr<G4DNAChemistryManager>(new G4DNAChemistryManager());
67  return fInstance.get();
68 }
69 
70 G4DNAChemistryManager::~G4DNAChemistryManager()
71 {
72  if (fOutput.is_open())
73  {
74  fOutput.close();
75  }
76  if(fIonisationLevel) delete fIonisationLevel;
77  if(fExcitationLevel) delete fExcitationLevel;
81  fInstance.release();
83 }
84 
86 {
87  if(fInstance.get())
88  fInstance.reset();
89 }
90 
92  ios_base::openmode mode)
93 {
94  fOutput.open(output.data(), mode);
95  fOutput << std::setprecision(6) << std::scientific;
96 
97  fOutput << setw(11) << left << "#Parent ID"
98  << setw(10) << "Molecule"
99  << setw(14) << "Elec Modif"
100  << setw(13) << "Energy (eV)"
101  << setw(22) << "X pos of parent [nm]"
102  << setw(22) << "Y pos of parent [nm]"
103  << setw(22) << "Z pos of parent [nm]"
104  << setw(14) << "X pos [nm]"
105  << setw(14) << "Y pos [nm]"
106  << setw(14) << "Z pos [nm]"
107  << G4endl
108  << setw(21) << "#"
109  << setw(13) << "1)io/ex=0/1"
110  << G4endl
111  << setw(21) << "#"
112  << setw(13) << "2)level=0...5"
113  << G4endl;
114  fWriteFile = true;
115 }
116 
118 {
119  if (fOutput.is_open())
120  {
121  fOutput.close();
122  }
123  fWriteFile = false;
124 }
125 
127 {
128  if(!fExcitationLevel)
129  {
130  fExcitationLevel = new G4DNAWaterExcitationStructure;
131  }
132  return fExcitationLevel;
133 }
134 
136 {
137  if(!fIonisationLevel)
138  {
139  fIonisationLevel = new G4DNAWaterIonisationStructure;
140  }
141  return fIonisationLevel;
142 }
143 
145  G4int electronicLevel,
146  const G4Track* theIncomingTrack)
147 {
148  if(fWriteFile)
149  {
150  G4double energy = -1.;
151 
152  switch (modification)
153  {
155  energy = -1;
156  break;
157  case eExcitedMolecule :
158  energy = GetExcitationLevel()->ExcitationEnergy(electronicLevel);
159  break;
160  case eIonizedMolecule :
161  energy = GetIonisationLevel()->IonisationEnergy(electronicLevel);
162  break;
163  }
164 
165  fOutput << setw(11) << left << theIncomingTrack->GetTrackID()
166  << setw(10) << "H2O"
167  << left << modification
168  << internal <<":"
169  << right <<electronicLevel
170  << left
171  << setw(11) << ""
172  << std::setprecision(2) << std::fixed
173  << setw(13) << energy/eV
174  << std::setprecision(6) << std::scientific
175  << setw(22) << (theIncomingTrack->GetPosition().x())/nanometer
176  << setw(22) << (theIncomingTrack->GetPosition().y())/nanometer
177  << setw(22) << (theIncomingTrack->GetPosition().z())/nanometer
178  << G4endl;
179  }
180 
181  if(fActiveChemistry)
182  {
183  G4Molecule * H2O = new G4Molecule (G4H2O::Definition());
184 
185  switch (modification)
186  {
188  H2O -> AddElectron(5,1);
189  break;
190  case eExcitedMolecule :
191  H2O -> ExciteMolecule(electronicLevel);
192  break;
193  case eIonizedMolecule :
194  H2O -> IonizeMolecule(electronicLevel);
195  break;
196  }
197 
198  G4Track * H2OTrack = H2O->BuildTrack(1*picosecond,
199  theIncomingTrack->GetPosition());
200 
201  H2OTrack -> SetParentID(theIncomingTrack->GetTrackID());
202  H2OTrack -> SetTrackStatus(fStopButAlive);
203  H2OTrack -> SetKineticEnergy(0.);
204 
206  }
207 }
208 
210  G4ThreeVector* finalPosition)
211 // finalPosition is a pointer because this argument is optional
212 {
213  if(fWriteFile)
214  {
215  fOutput << setw(11)<< theIncomingTrack->GetTrackID()
216  << setw(10)<< "e_aq"
217  << setw(14)<< -1
218  << std::setprecision(2) << std::fixed
219  << setw(13)<< theIncomingTrack->GetKineticEnergy()/eV
220  << std::setprecision(6) << std::scientific
221  << setw(22)<< (theIncomingTrack->GetPosition().x())/nanometer
222  << setw(22)<< (theIncomingTrack->GetPosition().y())/nanometer
223  << setw(22)<< (theIncomingTrack->GetPosition().z())/nanometer ;
224 
225  if(finalPosition != 0)
226  {
227  fOutput<< setw(14)<< (finalPosition->x())/nanometer
228  << setw(14)<< (finalPosition->y())/nanometer
229  << setw(14)<< (finalPosition->z())/nanometer ;
230  }
231 
232  fOutput << G4endl;
233  }
234 
235  if(fActiveChemistry)
236  {
238  G4Track * e_aqTrack(0);
239  if(finalPosition)
240  {
241  e_aqTrack = e_aq->BuildTrack(picosecond,*finalPosition);
242  }
243  else
244  {
245  e_aqTrack = e_aq->BuildTrack(picosecond,theIncomingTrack->GetPosition());
246  }
247  e_aqTrack -> SetTrackStatus(fAlive);
248  e_aqTrack -> SetParentID(theIncomingTrack->GetTrackID());
249  G4ITTrackHolder::Instance()->PushTrack(e_aqTrack);
250  G4ITManager<G4Molecule>::Instance()->Push(e_aqTrack);
251  }
252 }
253 
254 
255 void G4DNAChemistryManager::PushMolecule(G4Molecule*& molecule, double time,
256  const G4ThreeVector& position, int parentID)
257 {
258  if(fWriteFile)
259  {
260  fOutput << setw(11)<< parentID
261  << setw(10)<< molecule->GetName()
262  << setw(14)<< -1
263  << std::setprecision(2) << std::fixed
264  << setw(13)<< -1
265  << std::setprecision(6) << std::scientific
266  << setw(22)<< (position.x())/nanometer
267  << setw(22)<< (position.y())/nanometer
268  << setw(22)<< (position.z())/nanometer;
269  fOutput << G4endl;
270  }
271 
272  if(fActiveChemistry)
273  {
274  G4Track* track = molecule->BuildTrack(time,position);
275  track -> SetTrackStatus(fAlive);
276  track -> SetParentID(parentID);
278  G4ITManager<G4Molecule>::Instance()->Push(track);
279  }
280  else
281  {
282  delete molecule;
283  molecule = 0;
284  }
285 }
286 
288  const G4Track* theIncomingTrack)
289 {
290  if(fWriteFile)
291  {
292  fOutput << setw(11)<< theIncomingTrack->GetTrackID()
293  << setw(10)<< molecule->GetName()
294  << setw(14)<< -1
295  << std::setprecision(2) << std::fixed
296  << setw(13)<< theIncomingTrack->GetKineticEnergy()/eV
297  << std::setprecision(6) << std::scientific
298  << setw(22)<< (theIncomingTrack->GetPosition().x())/nanometer
299  << setw(22)<< (theIncomingTrack->GetPosition().y())/nanometer
300  << setw(22)<< (theIncomingTrack->GetPosition().z())/nanometer ;
301  fOutput << G4endl;
302  }
303 
304  if(fActiveChemistry)
305  {
306  G4Track* track = molecule->BuildTrack(theIncomingTrack->GetGlobalTime(),theIncomingTrack->GetPosition());
307  track -> SetTrackStatus(fAlive);
308  track -> SetParentID(theIncomingTrack->GetTrackID());
310  G4ITManager<G4Molecule>::Instance()->Push(track);
311  }
312  else
313  {
314  delete molecule;
315  molecule = 0;
316  }
317 }