Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/examples/extended/medical/dna/dnadamage2/src/G4EmDNAChemistryForPlasmids.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  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 //
 27 // G4EmDNAChemistryForPlasmids.cc
 28 //
 29 //  Created on: Feb 10, 2021
 30 //      Authors: J. Naoki D. Kondo
 31 //               W. G. Shin, J. Ramos-Mendez and B. Faddegon
 32 //
 33 /// \file G4EmDNAChemistryForPlasmids.cc
 34 /// \brief Implementation of the Chemistry parameters with DNA reactions
 35 
 36 #include "G4EmDNAChemistryForPlasmids.hh"
 37 
 38 #include "G4DNAChemistryManager.hh"
 39 #include "G4DNAGenericIonsManager.hh"
 40 #include "G4DNAWaterDissociationDisplacer.hh"
 41 #include "G4DNAWaterExcitationStructure.hh"
 42 #include "G4PhysicalConstants.hh"
 43 #include "G4ProcessManager.hh"
 44 #include "G4SystemOfUnits.hh"
 45 
 46 // *** Processes and models for Geant4-DNA
 47 
 48 #include "G4DNABrownianTransportation.hh"
 49 #include "G4DNAElectronHoleRecombination.hh"
 50 #include "G4DNAElectronSolvation.hh"
 51 #include "G4DNAIRT.hh"
 52 #include "G4DNAMolecularDissociation.hh"
 53 #include "G4DNAMolecularIRTModel.hh"
 54 #include "G4DNAMolecularReactionTable.hh"
 55 #include "G4DNASancheExcitationModel.hh"
 56 #include "G4DNAVibExcitation.hh"
 57 #include "G4VDNAReactionModel.hh"
 58 
 59 // particles
 60 
 61 #include "PlasmidMolecules.hh"
 62 #include "ScavengerMolecules.hh"
 63 
 64 #include "G4BuilderType.hh"
 65 #include "G4Electron.hh"
 66 #include "G4Electron_aq.hh"
 67 #include "G4FakeMolecule.hh"
 68 #include "G4GenericIon.hh"
 69 #include "G4H2.hh"
 70 #include "G4H2O.hh"
 71 #include "G4H2O2.hh"
 72 #include "G4H3O.hh"
 73 #include "G4HO2.hh"
 74 #include "G4Hydrogen.hh"
 75 #include "G4MoleculeTable.hh"
 76 #include "G4O2.hh"
 77 #include "G4O3.hh"
 78 #include "G4OH.hh"
 79 #include "G4Oxygen.hh"
 80 #include "G4PhysicsListHelper.hh"
 81 #include "G4Proton.hh"
 82 
 83 /****/
 84 #include "G4DNAMoleculeEncounterStepper.hh"
 85 #include "G4MolecularConfiguration.hh"
 86 #include "G4ProcessTable.hh"
 87 /****/
 88 // factory
 89 #include "G4PhysicsConstructorFactory.hh"
 90 #include "G4ChemDissociationChannels_option1.hh"
 91 G4_DECLARE_PHYSCONSTR_FACTORY(G4EmDNAChemistryForPlasmids);
 92 
 93 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
 94 
 95 G4EmDNAChemistryForPlasmids::G4EmDNAChemistryForPlasmids() : G4VUserChemistryList(true)
 96 {
 97   G4DNAChemistryManager::Instance()->SetChemistryList(this);
 98 
 99   fDMSO = 0;
100   fOxygen = 0;
101 }
102 
103 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
104 
105 G4EmDNAChemistryForPlasmids::G4EmDNAChemistryForPlasmids(G4double dmso, G4double oxygen)
106   : G4VUserChemistryList(true)
107 {
108   G4DNAChemistryManager::Instance()->SetChemistryList(this);
109 
110   fDMSO = dmso;
111   fOxygen = oxygen;
112 }
113 
114 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
115 
116 void G4EmDNAChemistryForPlasmids::ConstructMolecule()
117 {
118   G4ChemDissociationChannels_option1::ConstructMolecule();
119   auto molTab = G4MoleculeTable::Instance();
120 
121   molTab->CreateConfiguration("DMSO", G4DMSO::Definition(), 0, 0 * (m2 / s));
122 
123   molTab->CreateConfiguration("O2", G4O2::Definition());
124   molTab->GetConfiguration("O2")->SetVanDerVaalsRadius(0.17 * nm);
125   molTab->CreateConfiguration("Oxygen", G4OxygenB::Definition(), 0,
126                                                    0 * (m2 / s));
127   molTab->CreateConfiguration("Deoxyribose", G4DNA_Deoxyribose::Definition(),
128                                                    0, 1E-150 * (m2 / s));
129   molTab->CreateConfiguration(
130     "Damaged_DeoxyriboseOH", G4DNA_DamagedDeoxyriboseOH::Definition(), 0, 1E-150 * (m2 / s));
131   molTab->CreateConfiguration(
132     "Damaged_DeoxyriboseH", G4DNA_DamagedDeoxyriboseH::Definition(), 0, 1E-150 * (m2 / s));
133   molTab->CreateConfiguration(
134     "Damaged_DeoxyriboseEAQ", G4DNA_DamagedDeoxyriboseEAQ::Definition(), 0, 1E-150 * (m2 / s));
135 }
136 
137 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
138 
139 void G4EmDNAChemistryForPlasmids::ConstructDissociationChannels()
140 {
141   G4ChemDissociationChannels_option1::ConstructDissociationChannels();
142 }
143 
144 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
145 
146 void G4EmDNAChemistryForPlasmids::ConstructReactionTable(
147   G4DNAMolecularReactionTable* theReactionTable)
148 {
149   // fReactionTable = theReactionTable;
150 
151   //-----------------------------------
152   // Get the molecular configuration
153   auto molTab = G4MoleculeTable::Instance();
154 
155   auto OH = molTab->GetConfiguration("°OH");
156   auto OHm = molTab->GetConfiguration("OHm");
157   auto e_aq = molTab->GetConfiguration("e_aq");
158   auto H2 = molTab->GetConfiguration("H2");
159   auto H3Op = molTab->GetConfiguration("H3Op");
160   auto H = molTab->GetConfiguration("H");
161   auto H2O2 = molTab->GetConfiguration("H2O2");
162 
163   //----------------------------------------------------------------//
164   // Type II                                                        //
165   //----------------------------------------------------------------//
166   // e_aq + *OH -> OH-    prob: 0.49
167   auto reactionData =
168     new G4DNAMolecularReactionData(2.95e10 * (1e-3 * m3 / (mole * s)), e_aq, OH);
169   reactionData->AddProduct(OHm);
170   reactionData->SetReactionType(1);
171   theReactionTable->SetReaction(reactionData);
172   //----------------------------------------------------------------//
173   // e_aq + H2O2 -> OH- + *OH        prob: 0.11
174   reactionData = new G4DNAMolecularReactionData(1.10e10 * (1e-3 * m3 / (mole * s)), e_aq, H2O2);
175   reactionData->AddProduct(OHm);
176   reactionData->AddProduct(OH);
177   reactionData->SetReactionType(1);
178   theReactionTable->SetReaction(reactionData);
179   //----------------------------------------------------------------//
180   // *OH + *H -> H2O        prob: 0.33
181   reactionData = new G4DNAMolecularReactionData(1.55e10 * (1e-3 * m3 / (mole * s)), OH, H);
182   reactionData->SetReactionType(1);
183   theReactionTable->SetReaction(reactionData);
184   //----------------------------------------------------------------//
185   // H + H2O2 -> OH                prob: 0.00
186   reactionData = new G4DNAMolecularReactionData(0.009e10 * (1e-3 * m3 / (mole * s)), H, H2O2);
187   reactionData->AddProduct(OH);
188   reactionData->SetReactionType(1);
189   theReactionTable->SetReaction(reactionData);
190   //----------------------------------------------------------------//
191   // *OH + *OH -> H2O2                prob: 0.55
192   reactionData = new G4DNAMolecularReactionData(0.55e10 * (1e-3 * m3 / (mole * s)), OH, OH);
193   reactionData->AddProduct(H2O2);
194   reactionData->SetReactionType(1);
195   theReactionTable->SetReaction(reactionData);
196   //----------------------------------------------------------------//
197 
198   //----------------------------------------------------------------//
199   // Type III                                                       //
200   //----------------------------------------------------------------//
201   // e_aq + e_aq + 2H2O -> H2 + 2OH-
202   reactionData = new G4DNAMolecularReactionData(0.636e10 * (1e-3 * m3 / (mole * s)), e_aq, e_aq);
203   reactionData->AddProduct(OHm);
204   reactionData->AddProduct(OHm);
205   reactionData->AddProduct(H2);
206   reactionData->SetReactionType(0);
207   theReactionTable->SetReaction(reactionData);
208   //------------------------------------------------------------------
209   // H3O+ + OH- -> 2H2O
210   reactionData = new G4DNAMolecularReactionData(11.3e10 * (1e-3 * m3 / (mole * s)), H3Op, OHm);
211   reactionData->SetReactionType(0);
212   theReactionTable->SetReaction(reactionData);
213   //------------------------------------------------------------------
214 
215   //----------------------------------------------------------------//
216   // Type IV                                                        //
217   //----------------------------------------------------------------//
218   // e_aq + H3O+ -> H* + H2O        prob: 0.04
219   reactionData = new G4DNAMolecularReactionData(2.11e10 * (1e-3 * m3 / (mole * s)), e_aq, H3Op);
220   reactionData->AddProduct(H);
221   reactionData->SetReactionType(1);
222   theReactionTable->SetReaction(reactionData);
223   //----------------------------------------------------------------//
224 
225   //----------------------------------------------------------------//
226   // Type V                                                         //
227   // First order reaction                                           //
228   //----------------------------------------------------------------//
229   // e_aq + *H -> OH- + H2
230   reactionData = new G4DNAMolecularReactionData(2.5e10 * (1e-3 * m3 / (mole * s)), e_aq, H3Op);
231   reactionData->AddProduct(OHm);
232   reactionData->AddProduct(H2);
233   reactionData->SetReactionType(0);
234   theReactionTable->SetReaction(reactionData);
235   //----------------------------------------------------------------//
236   // *H + *H -> H2
237   reactionData = new G4DNAMolecularReactionData(0.503e10 * (1e-3 * m3 / (mole * s)), H, H);
238   reactionData->AddProduct(H2);
239   reactionData->SetReactionType(0);
240   theReactionTable->SetReaction(reactionData);
241 
242   if (fDMSO > 0) DeclareDMSOAndDNAReactions(theReactionTable);
243 
244   if (fOxygen > 0) DeclareOxygenReactions(theReactionTable);
245 }
246 
247 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
248 
249 void G4EmDNAChemistryForPlasmids::ConstructProcess()
250 {
251   G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
252 
253   //===============================================================
254   // Extend vibrational to low energy
255   // Anyway, solvation of electrons is taken into account from 7.4 eV
256   // So below this threshold, for now, no accurate modeling is done
257   //
258   G4VProcess* process =
259     G4ProcessTable::GetProcessTable()->FindProcess("e-_G4DNAVibExcitation", "e-");
260 
261   if (process) {
262     G4DNAVibExcitation* vibExcitation = (G4DNAVibExcitation*)process;
263     G4VEmModel* model = vibExcitation->EmModel();
264     G4DNASancheExcitationModel* sancheExcitationMod =
265       dynamic_cast<G4DNASancheExcitationModel*>(model);
266     if (sancheExcitationMod) {
267       sancheExcitationMod->ExtendLowEnergyLimit(0.025 * eV);
268     }
269   }
270 
271   //===============================================================
272   // *** Electron Solvatation ***
273   //
274   process = G4ProcessTable::GetProcessTable()->FindProcess("e-_G4DNAElectronSolvation", "e-");
275 
276   if (process == 0) {
277     ph->RegisterProcess(new G4DNAElectronSolvation("e-_G4DNAElectronSolvation"),
278                         G4Electron::Definition());
279   }
280 
281   //===============================================================
282   // Define processes for molecules
283   //
284   G4MoleculeTable* theMoleculeTable = G4MoleculeTable::Instance();
285   G4MoleculeDefinitionIterator iterator = theMoleculeTable->GetDefintionIterator();
286   iterator.reset();
287   while (iterator()) {
288     G4MoleculeDefinition* moleculeDef = iterator.value();
289 
290     if (moleculeDef != G4H2O::Definition()) {
291       // G4cout << "Brownian motion added for: "<< moleculeDef->GetName() << G4endl;
292       //      G4DNABrownianTransportation* brown = new G4DNABrownianTransportation();
293       //      ph->RegisterProcess(brown, moleculeDef);
294     }
295     else {
296       moleculeDef->GetProcessManager()->AddRestProcess(new G4DNAElectronHoleRecombination(), 2);
297       G4DNAMolecularDissociation* dissociationProcess =
298         new G4DNAMolecularDissociation("H2O_DNAMolecularDecay");
299       dissociationProcess->SetDisplacer(moleculeDef, new G4DNAWaterDissociationDisplacer);
300       dissociationProcess->SetVerboseLevel(3);
301 
302       moleculeDef->GetProcessManager()->AddRestProcess(dissociationProcess, 1);
303     }
304     /*
305      * Warning : end of particles and processes are needed by
306      * EM Physics builders
307      */
308   }
309 
310   G4DNAChemistryManager::Instance()->Initialize();
311 }
312 
313 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
314 
315 void G4EmDNAChemistryForPlasmids::ConstructTimeStepModel(G4DNAMolecularReactionTable*
316                                                          /*reactionTable*/)
317 {
318   RegisterTimeStepModel(new G4DNAMolecularIRTModel(), 0);
319 }
320 
321 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
322 
323 void G4EmDNAChemistryForPlasmids::DeclareDMSOAndDNAReactions(
324   G4DNAMolecularReactionTable* theReactionTable)
325 {
326   auto molTab = G4MoleculeTable::Instance();
327   auto DMSO = molTab->GetConfiguration("DMSO");
328   auto OH = molTab->GetConfiguration("°OH");
329   auto H = molTab->GetConfiguration("H");
330   auto e_aq = molTab->GetConfiguration("e_aq");
331   auto deoxyribose =
332       molTab->GetConfiguration("Deoxyribose");
333   auto damage_deoxyribose_OH =
334       molTab->GetConfiguration("Damaged_DeoxyriboseOH");
335   auto damage_deoxyribose_H =
336       molTab->GetConfiguration("Damaged_DeoxyriboseH");
337   auto damage_deoxyribose_eaq =
338       molTab->GetConfiguration("Damaged_DeoxyriboseEAQ");
339 
340   G4double DNA_OH_Rate = 1.32E7 * std::pow(fDMSO * 7.1E9, 0.29);
341 
342   // DMSO Reactions
343   // OH + DMSO
344   auto reactionData = new G4DNAMolecularReactionData(fDMSO * 7.1E9 / s, OH, DMSO);
345   theReactionTable->SetReaction(reactionData);
346   //----------------------------------------------------------------//
347   // H + DMSO
348   reactionData = new G4DNAMolecularReactionData(fDMSO * 2.7E7 / s, H, DMSO);
349   theReactionTable->SetReaction(reactionData);
350   //----------------------------------------------------------------//
351   // e-_aq + DMSO
352   reactionData = new G4DNAMolecularReactionData(fDMSO * 3.8E6 / s, e_aq, DMSO);
353   theReactionTable->SetReaction(reactionData);
354   //----------------------------------------------------------------//
355 
356   //----------------------------------------------------------------//
357   // DNA Type                                                       //
358   //----------------------------------------------------------------//
359   // DNA + OH -> Damage
360   reactionData =
361     new G4DNAMolecularReactionData(DNA_OH_Rate * (1e-3 * m3 / (mole * s)), deoxyribose, OH);
362   reactionData->AddProduct(damage_deoxyribose_OH);
363   theReactionTable->SetReaction(reactionData);
364   //----------------------------------------------------------------//
365   // DNA + H
366   reactionData = new G4DNAMolecularReactionData(0.03e9 * (1e-3 * m3 / (mole * s)), deoxyribose, H);
367   reactionData->AddProduct(damage_deoxyribose_H);
368   theReactionTable->SetReaction(reactionData);
369   //----------------------------------------------------------------//
370   // DNA + e-_aq
371   reactionData =
372     new G4DNAMolecularReactionData(0.01e9 * (1e-3 * m3 / (mole * s)), deoxyribose, e_aq);
373   reactionData->AddProduct(damage_deoxyribose_eaq);
374   theReactionTable->SetReaction(reactionData);
375   //----------------------------------------------------------------//
376 }
377 
378 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
379 
380 void G4EmDNAChemistryForPlasmids::DeclareOxygenReactions(
381   G4DNAMolecularReactionTable* theReactionTable)
382 {
383   auto molTab = G4MoleculeTable::Instance();
384   auto Oxygen = molTab->GetConfiguration("Oxygen");
385   auto H = molTab->GetConfiguration("H");
386   auto e_aq = molTab->GetConfiguration("e_aq");
387   auto O2m = molTab->GetConfiguration("O2m");
388   auto HO2 = molTab->GetConfiguration("HO2°");
389   auto O2 = molTab->GetConfiguration("O2");
390   auto OH = molTab->GetConfiguration("°OH");
391 
392   // Oxygen Reactions
393   // e_aq + O2B
394   auto reactionData =
395     new G4DNAMolecularReactionData((fOxygen * 1.9E10) / s, e_aq, Oxygen);
396   reactionData->AddProduct(O2m);
397   theReactionTable->SetReaction(reactionData);
398   //------------------------------------------------------------------
399   // H + O2B
400   reactionData = new G4DNAMolecularReactionData((fOxygen * 2.1E10) / s, H, Oxygen);
401   reactionData->AddProduct(HO2);
402   theReactionTable->SetReaction(reactionData);
403   //------------------------------------------------------------------
404   // e_aq + O2
405   reactionData = new G4DNAMolecularReactionData(1.9E10 * (1e-3 * m3 / (mole * s)), e_aq, O2);
406   reactionData->AddProduct(O2m);
407   reactionData->SetReactionType(1);
408   theReactionTable->SetReaction(reactionData);
409   //------------------------------------------------------------------
410   // H + O2
411   reactionData = new G4DNAMolecularReactionData(2.1E10 * (1e-3 * m3 / (mole * s)), H, O2);
412   reactionData->AddProduct(HO2);
413   reactionData->SetReactionType(1);
414   theReactionTable->SetReaction(reactionData);
415   //------------------------------------------------------------------
416   // OH + HO2
417   reactionData = new G4DNAMolecularReactionData(7.9E9 * (1e-3 * m3 / (mole * s)), OH, HO2);
418   reactionData->AddProduct(O2);
419   reactionData->SetReactionType(1);
420   theReactionTable->SetReaction(reactionData);
421   //------------------------------------------------------------------
422 }
423 
424 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
425