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 // This example is provided by the Geant4-DNA collaboration 27 // chem6 example is derived from chem4 and chem5 examples 28 // 29 // Any report or published results obtained using the Geant4-DNA software 30 // shall cite the following Geant4-DNA collaboration publication: 31 // J. Appl. Phys. 125 (2019) 104301 32 // Med. Phys. 45 (2018) e722-e739 33 // J. Comput. Phys. 274 (2014) 841-882 34 // Med. Phys. 37 (2010) 4692-4708 35 // Int. J. Model. Simul. Sci. Comput. 1 (2010) 157-178 36 // The Geant4-DNA web site is available at http://geant4-dna.org 37 // 38 // Authors: W. G. Shin and S. Incerti (CENBG, France) 39 // 40 // $Id$ 41 // 42 /// \file PrimaryKiller.cc 43 /// \brief Implementation of the PrimaryKiller class 44 45 #include "PrimaryKiller.hh" 46 47 #include <G4Event.hh> 48 #include <G4RunManager.hh> 49 #include <G4SystemOfUnits.hh> 50 #include <G4UIcmdWith3VectorAndUnit.hh> 51 #include <G4UIcmdWithADoubleAndUnit.hh> 52 #include <G4UnitsTable.hh> 53 54 /** \file PrimaryKiller.cc 55 \class PrimaryKiller 56 57 Kill the primary particle: 58 - either after a given energy loss 59 - or after the primary particle has reached a given energy 60 */ 61 62 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 63 64 PrimaryKiller::PrimaryKiller(G4String name, G4int depth) 65 : G4VPrimitiveScorer(name, depth), G4UImessenger() 66 { 67 fELoss = 0.; // cumulated energy for current event 68 69 fELossRange_Min = DBL_MAX; // fELoss from which the primary is killed 70 fELossRange_Max = DBL_MAX; // fELoss from which the event is aborted 71 fKineticE_Min = 0; // kinetic energy below which the primary is killed 72 fPhantomSize = G4ThreeVector(1 * km, 1 * km, 1 * km); 73 74 fpELossUI = new G4UIcmdWithADoubleAndUnit("/primaryKiller/eLossMin", this); 75 fpAbortEventIfELossUpperThan = new G4UIcmdWithADoubleAndUnit("/primaryKiller/eLossMax", this); 76 fpMinKineticE = new G4UIcmdWithADoubleAndUnit("/primaryKiller/minKineticE", this); 77 fpSizeUI = new G4UIcmdWith3VectorAndUnit("/primaryKiller/setSize", this); 78 fpSizeUI->SetDefaultUnit("um"); 79 } 80 81 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 82 83 PrimaryKiller::~PrimaryKiller() 84 { 85 delete fpELossUI; 86 delete fpAbortEventIfELossUpperThan; 87 delete fpSizeUI; 88 } 89 90 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 91 92 void PrimaryKiller::SetNewValue(G4UIcommand* command, G4String newValue) 93 { 94 if (command == fpELossUI) { 95 fELossRange_Min = fpELossUI->GetNewDoubleValue(newValue); 96 } 97 else if (command == fpAbortEventIfELossUpperThan) { 98 fELossRange_Max = fpAbortEventIfELossUpperThan->GetNewDoubleValue(newValue); 99 } 100 else if (command == fpSizeUI) { 101 fPhantomSize = fpSizeUI->GetNew3VectorValue(newValue); 102 } 103 } 104 105 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 106 107 G4bool PrimaryKiller::ProcessHits(G4Step* aStep, G4TouchableHistory*) 108 { 109 const G4Track* track = aStep->GetTrack(); 110 G4ThreeVector pos = aStep->GetPostStepPoint()->GetPosition(); 111 112 if (std::abs(pos.x()) > fPhantomSize.getX() / 2 || std::abs(pos.y()) > fPhantomSize.getY() / 2 113 || std::abs(pos.z()) > fPhantomSize.getZ() / 2) 114 { 115 ((G4Track*)track)->SetTrackStatus(fStopAndKill); 116 return false; 117 } 118 119 if (track->GetTrackID() != 1 || track->GetParticleDefinition()->GetPDGEncoding() != 11) 120 return FALSE; 121 122 //------------------- 123 124 double kineticE = aStep->GetPostStepPoint()->GetKineticEnergy(); 125 126 G4double eLoss = aStep->GetPreStepPoint()->GetKineticEnergy() - kineticE; 127 128 if (eLoss == 0.) return FALSE; 129 130 //------------------- 131 132 fELoss += eLoss; 133 134 if (fELoss > fELossRange_Max) { 135 G4RunManager::GetRunManager()->AbortEvent(); 136 /* int eventID = 137 G4EventManager::GetEventManager()->GetConstCurrentEvent()->GetEventID(); 138 139 G4cout << " * PrimaryKiller: aborts event " << eventID <<" energy loss " 140 "is too large. \n" 141 << " * Energy loss by primary is: " 142 << G4BestUnit(fELoss, "Energy") 143 << ". Event is aborted if the Eloss is > " 144 << G4BestUnit(fELossRange_Max, "Energy") 145 << G4endl; 146 */ 147 } 148 149 if (fELoss >= fELossRange_Min || kineticE <= fKineticE_Min) { 150 ((G4Track*)track)->SetTrackStatus(fStopAndKill); 151 // G4cout << "kill track at : "<<'\n'; 152 // << G4BestUnit(kineticE, "Energy") 153 // << ", E loss is: " 154 // << G4BestUnit(fELoss, "Energy") 155 // << " /fELossMax: " 156 // << G4BestUnit(fELossMax, "Energy") 157 // << ", EThreshold is: " 158 // << G4BestUnit(fEThreshold, "Energy") 159 // << G4endl; 160 } 161 162 return TRUE; 163 } 164 165 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 166 167 void PrimaryKiller::Initialize(G4HCofThisEvent* /*HCE*/) 168 { 169 fELoss = 0.; 170 } 171 172 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 173 174 void PrimaryKiller::EndOfEvent(G4HCofThisEvent*) {} 175 176 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 177 178 void PrimaryKiller::DrawAll() 179 { 180 ; 181 } 182 183 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 184 185 void PrimaryKiller::PrintAll() {} 186 187 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 188