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 /// \file radiobiology/src/SD.cc 28 /// \brief Implementation of the RadioBio::SD class 29 30 #include "SD.hh" 31 32 #include "G4AccumulableManager.hh" 33 #include "G4HCofThisEvent.hh" 34 #include "G4SDManager.hh" 35 #include "G4Step.hh" 36 #include "G4ThreeVector.hh" 37 38 #include "Hit.hh" 39 #include "Manager.hh" 40 #include "VRadiobiologicalAccumulable.hh" 41 42 namespace RadioBio 43 { 44 45 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 46 47 SD::SD(const G4String& name, const G4String& hitsCollectionName) : G4VSensitiveDetector(name) 48 { 49 collectionName.insert(hitsCollectionName); 50 } 51 52 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 53 54 void SD::Initialize(G4HCofThisEvent* hce) 55 { 56 // Create hits collection 57 fHitsCollection = new RadioBioHitsCollection(SensitiveDetectorName, collectionName[0]); 58 59 // Add this collection in hce 60 G4int hcID = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]); 61 hce->AddHitsCollection(hcID, fHitsCollection); 62 } 63 64 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 65 66 G4bool SD::ProcessHits(G4Step* aStep, G4TouchableHistory*) 67 { 68 // FIXME why the namespace specification? Should not be necessary... 69 RadioBio::Hit* newHit = new RadioBio::Hit(); 70 71 // Get the pre-step kinetic energy 72 G4double eKinPre = aStep->GetPreStepPoint()->GetKineticEnergy(); 73 // Get the post-step kinetic energy 74 G4double eKinPost = aStep->GetPostStepPoint()->GetKineticEnergy(); 75 // Get the step average kinetic energy 76 G4double eKinMean = (eKinPre + eKinPost) * 0.5; 77 78 const std::vector<const G4Track*>* secondary = aStep->GetSecondaryInCurrentStep(); 79 80 size_t SecondarySize = (*secondary).size(); 81 G4double EnergySecondary = 0.; 82 83 // Get secondary electrons energy deposited 84 if (SecondarySize) // Calculate only secondary particles 85 { 86 for (size_t numsec = 0; numsec < SecondarySize; numsec++) { 87 // Get the PDG code of every secondaty particles in current step 88 G4int PDGSecondary = (*secondary)[numsec]->GetDefinition()->GetPDGEncoding(); 89 90 if (PDGSecondary == 11) // Calculate only secondary electrons 91 { 92 // Calculate the energy deposit of secondary electrons in current step 93 EnergySecondary += (*secondary)[numsec]->GetKineticEnergy(); 94 } 95 } 96 } 97 98 // Update the hit with the necessary quantities 99 newHit->SetTrackID(aStep->GetTrack()->GetTrackID()); 100 newHit->SetPartType(aStep->GetTrack()->GetParticleDefinition()); 101 newHit->SetEkinMean(eKinMean); 102 newHit->SetDeltaE(aStep->GetTotalEnergyDeposit()); 103 newHit->SetEinit(aStep->GetPreStepPoint()->GetKineticEnergy()); 104 newHit->SetTrackLength(aStep->GetStepLength()); 105 newHit->SetElectronEnergy(EnergySecondary); 106 newHit->SetPhysicalVolume(aStep->GetPreStepPoint()->GetPhysicalVolume()); 107 newHit->SetVoxelIndexes(aStep->GetPreStepPoint()->GetTouchableHandle()); 108 109 // Insert the hit in the hitcollection 110 fHitsCollection->insert(newHit); 111 112 // Accumulables are accumulated only if calculation is enabled 113 for (G4int i = 0; i < G4AccumulableManager::Instance()->GetNofAccumulables(); ++i) { 114 // Get the accumulable from the proper manager 115 G4VAccumulable* GenAcc = G4AccumulableManager::Instance()->GetAccumulable(i); 116 117 // Get the quantity from the proper manager using the name 118 auto q = Manager::GetInstance()->GetQuantity(GenAcc->GetName()); 119 120 VRadiobiologicalAccumulable* radioAcc = dynamic_cast<VRadiobiologicalAccumulable*>(GenAcc); 121 122 // If the dynamic_cast did not work, this means that accumulable 123 // was not a VRadiobiologicalAccumulable 124 if (radioAcc == nullptr) continue; 125 126 // If calculation is enabled, accumulate 127 if (q->IsCalculationEnabled()) radioAcc->Accumulate(newHit); 128 } 129 130 return true; 131 } 132 133 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 134 135 void SD::EndOfEvent(G4HCofThisEvent*) 136 { 137 if (verboseLevel > 1) { 138 G4int nofHits = fHitsCollection->entries(); 139 G4cout << G4endl << "-------->Hits Collection: in this event they are " << nofHits 140 << " hits in the detector slices: " << G4endl; 141 for (G4int i = 0; i < nofHits; i++) 142 (*fHitsCollection)[i]->Print(); 143 } 144 } 145 146 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 147 148 } // namespace RadioBio