Geant4 Cross Reference |
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 // Author: Mathieu Karamitros 28 // 29 // WARNING : This class is released as a prototype. 30 // It might strongly evolve or even disappear in the next releases. 31 // 32 // History: 33 // ----------- 34 // 13 Nov 2016 M.Karamitros created 35 // 36 // ------------------------------------------------------------------- 37 38 #include "G4DNAChemistryManager.hh" 39 #include "G4DNAMolecularMaterial.hh" 40 #include "G4DNAWaterExcitationStructure.hh" 41 #include "G4ITNavigator.hh" 42 #include "G4Navigator.hh" 43 #include "G4NistManager.hh" 44 #include "G4ParticleChangeForGamma.hh" 45 #include "G4PhysicalConstants.hh" 46 #include "G4SystemOfUnits.hh" 47 #include "G4TransportationManager.hh" 48 49 #include <memory> 50 51 //#define MODEL_VERBOSE 52 53 //------------------------------------------------------------------------------ 54 55 template<typename MODEL> 56 G4TDNAOneStepThermalizationModel<MODEL>:: 57 G4TDNAOneStepThermalizationModel(const G4ParticleDefinition*, 58 const G4String& nam) : 59 G4VEmModel(nam) 60 { 61 fVerboseLevel = 0; 62 SetLowEnergyLimit(0.); 63 G4DNAWaterExcitationStructure exStructure; 64 SetHighEnergyLimit(exStructure.ExcitationEnergy(0)); 65 fpParticleChangeForGamma = nullptr; 66 fpWaterDensity = nullptr; 67 } 68 69 //------------------------------------------------------------------------------ 70 71 template<typename MODEL> 72 G4TDNAOneStepThermalizationModel<MODEL>::~G4TDNAOneStepThermalizationModel() 73 = default; 74 75 //------------------------------------------------------------------------------ 76 template<typename MODEL> 77 void G4TDNAOneStepThermalizationModel<MODEL>:: 78 Initialise(const G4ParticleDefinition* particleDefinition, 79 const G4DataVector&) 80 { 81 #ifdef MODEL_VERBOSE 82 if(fVerboseLevel) 83 G4cout << "Calling G4DNAOneStepThermalizationModel::Initialise()" 84 << G4endl; 85 #endif 86 if (particleDefinition->GetParticleName() != "e-") 87 { 88 G4ExceptionDescription errMsg; 89 errMsg << "G4DNAOneStepThermalizationModel can only be applied " 90 "to electrons"; 91 G4Exception("G4DNAOneStepThermalizationModel::CrossSectionPerVolume", 92 "G4DNAOneStepThermalizationModel001", 93 FatalErrorInArgument,errMsg); 94 return; 95 } 96 97 if(!fIsInitialised) 98 { 99 fIsInitialised = true; 100 fpParticleChangeForGamma = GetParticleChangeForGamma(); 101 } 102 103 G4Navigator* navigator = 104 G4TransportationManager::GetTransportationManager()-> 105 GetNavigatorForTracking(); 106 107 fpNavigator = std::make_unique<G4Navigator>(); 108 109 if(navigator != nullptr){ // add these checks for testing mode 110 auto world=navigator->GetWorldVolume(); 111 if(world != nullptr){ 112 fpNavigator->SetWorldVolume(world); 113 //fNavigator->NewNavigatorState(); 114 } 115 } 116 117 fpWaterDensity = 118 G4DNAMolecularMaterial::Instance()-> 119 GetNumMolPerVolTableFor(G4Material::GetMaterial("G4_WATER")); 120 } 121 122 //------------------------------------------------------------------------------ 123 template<typename MODEL> 124 G4double G4TDNAOneStepThermalizationModel<MODEL>:: 125 CrossSectionPerVolume(const G4Material* material, 126 const G4ParticleDefinition*, 127 G4double ekin, 128 G4double, 129 G4double) 130 { 131 #ifdef MODEL_VERBOSE 132 if(fVerboseLevel > 1) 133 G4cout << "Calling CrossSectionPerVolume() of G4DNAOneStepThermalizationModel" 134 << G4endl; 135 #endif 136 137 if(ekin > HighEnergyLimit()){ 138 return 0.0; 139 } 140 141 G4double waterDensity = (*fpWaterDensity)[material->GetIndex()]; 142 143 if(waterDensity!= 0.0){ 144 return DBL_MAX; 145 } 146 return 0.; 147 } 148 149 //------------------------------------------------------------------------------ 150 template<typename MODEL> 151 double G4TDNAOneStepThermalizationModel<MODEL>::GetRmean(double k){ 152 return MODEL::GetRmean(k); 153 } 154 155 156 //------------------------------------------------------------------------------ 157 158 template<typename MODEL> 159 void G4TDNAOneStepThermalizationModel<MODEL>:: 160 GetPenetration(G4double k, G4ThreeVector& displacement) 161 { 162 return MODEL::GetPenetration(k, displacement); 163 } 164 165 //------------------------------------------------------------------------------ 166 template<typename MODEL> 167 void G4TDNAOneStepThermalizationModel<MODEL>:: 168 SampleSecondaries(std::vector<G4DynamicParticle*>*, 169 const G4MaterialCutsCouple*, 170 const G4DynamicParticle* particle, 171 G4double, 172 G4double) 173 { 174 #ifdef MODEL_VERBOSE 175 if(fVerboseLevel) 176 G4cout << "Calling SampleSecondaries() of G4DNAOneStepThermalizationModel" 177 << G4endl; 178 #endif 179 180 G4double k = particle->GetKineticEnergy(); 181 182 if (k <= HighEnergyLimit()) 183 { 184 fpParticleChangeForGamma->ProposeTrackStatus(fStopAndKill); 185 fpParticleChangeForGamma->ProposeLocalEnergyDeposit(k); 186 187 if(G4DNAChemistryManager::IsActivated()) 188 { 189 G4ThreeVector displacement(0,0,0); 190 GetPenetration(k, displacement); 191 192 //______________________________________________________________ 193 const G4Track * theIncomingTrack = 194 fpParticleChangeForGamma->GetCurrentTrack(); 195 G4ThreeVector finalPosition(theIncomingTrack->GetPosition()+displacement); 196 197 fpNavigator->SetWorldVolume(theIncomingTrack->GetTouchable()-> 198 GetVolume(theIncomingTrack->GetTouchable()-> 199 GetHistoryDepth())); 200 201 double displacementMag = displacement.mag(); 202 double safety = DBL_MAX; 203 G4ThreeVector direction = displacement/displacementMag; 204 205 //-- 206 // 6/09/16 - recupere de molecular dissocation 207 double mag_displacement = displacement.mag(); 208 G4ThreeVector displacement_direction = displacement/mag_displacement; 209 210 // double step = DBL_MAX; 211 // step = fNavigator->CheckNextStep(theIncomingTrack->GetPosition(), 212 // displacement_direction, 213 // mag_displacement, 214 // safety); 215 // 216 // 217 // if(safety < mag_displacement) 218 // { 219 //// mag_displacement = prNewSafety; 220 // finalPosition = theIncomingTrack->GetPosition() 221 // + (displacement/displacementMag)*safety*0.80; 222 // } 223 //-- 224 225 fpNavigator->ResetHierarchyAndLocate(theIncomingTrack->GetPosition(), 226 direction, 227 *((G4TouchableHistory*) 228 theIncomingTrack->GetTouchable())); 229 230 fpNavigator->ComputeStep(theIncomingTrack->GetPosition(), 231 displacement/displacementMag, 232 displacementMag, 233 safety); 234 235 if(safety <= displacementMag) 236 { 237 finalPosition = theIncomingTrack->GetPosition() 238 + (displacement/displacementMag)*safety*0.80; 239 } 240 241 G4DNAChemistryManager::Instance()->CreateSolvatedElectron(theIncomingTrack, 242 &finalPosition); 243 244 fpParticleChangeForGamma->SetProposedKineticEnergy(25.e-3*eV); 245 } 246 } 247 }