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 /// \file electromagnetic/TestEm11/src/doiPETRun.cc 27 /// \brief Implementation of the doiPETRun class 28 // 29 // $Id: doiPETRun.cc 71376 2013-06-14 07:44:50Z maire $ 30 // 31 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 32 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 33 34 #include "doiPETRun.hh" 35 #include "doiPETDetectorConstruction.hh" 36 #include "doiPETPrimaryGeneratorAction.hh" 37 38 #include "G4UnitsTable.hh" 39 #include "G4SystemOfUnits.hh" 40 #include "G4Event.hh" 41 #include "doiPETAnalysis.hh" 42 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 43 44 doiPETRun::doiPETRun() 45 : G4Run() 46 { 47 // 48 totalTime = 0 * s; 49 prev_totalTime = 0 * s; 50 totalEdep = 0*keV; 51 } 52 53 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 54 55 doiPETRun::~doiPETRun() 56 {} 57 58 // 59 void doiPETRun::GetIntractionInfomation(InteractionInformation*step) 60 { 61 //The photon can have many interactions in one block. This stores all the interaction (including multiple copies of the same block.) 62 //Save the interacting block in the associative container including multiple copies. Since there will be multiple interactions in the crystals wthin a single block 63 mapBlockInteraction.insert( std::multimap< G4int, InteractionInformation* >::value_type( step->GetBlockNo(), step) ); 64 65 //The following stores the block number uniquely. Save the interacting block in the associative container without saving multiple copies 66 setBlockInteraction.insert( step->GetBlockNo()); 67 } 68 69 ///////// InteractingCrystal /////////////////////////////////// 70 71 void doiPETRun::FindInteractingCrystal() 72 { 73 std::multimap< G4int, InteractionInformation* >::iterator mitr; 74 std::set<G4int>::iterator sitr; 75 G4int cystalIDTemp; 76 77 for(sitr=setBlockInteraction.begin(); sitr!=setBlockInteraction.end(); sitr++){ 78 blockID = (*sitr); 79 mitr = mapBlockInteraction.find(blockID); 80 numberofInteractions = mapBlockInteraction.count(blockID); //number of interactions in which the interaction deposite energy in a block. 81 totalEdep=edepMax=0.; 82 for(G4int i=0; i<numberofInteractions; i++, mitr++){ 83 edep = (*mitr).second->GetEdep(); 84 cystalIDTemp = (*mitr).second->GetCrystalNo(); 85 edep_AfterCrystalBlurring = doiPETAnalysis::GetInstance()->QuantumEffifciency(edep, blockID, cystalIDTemp); 86 totalEdep += edep_AfterCrystalBlurring; 87 88 crystalID_vec.push_back(cystalIDTemp); 89 posInter_vec.push_back(((*mitr).second->GetInteractionPositionInCrystal())); 90 edepInCry_vec.push_back(edep_AfterCrystalBlurring); 91 92 //This is to identify the crystal with the highest energy deposition. 93 if(edepMax<edep_AfterCrystalBlurring){ 94 edepMax=edep_AfterCrystalBlurring; 95 crystalID = cystalIDTemp; 96 interactionTime = (*mitr).second->GetGlobalTime(); 97 98 //position of interaction based on highest edep. See below for energy weighted position of interaction 99 //interactionPos = (*mitr).second->GetInteractionPositionInCrystal(); 100 } 101 } 102 if(totalEdep==0)continue; 103 104 //The interaction position is based on energy weighted center of mass calculation 105 interactionPos = CenterOfMassInteractionPos(crystalID_vec, edepInCry_vec, totalEdep, posInter_vec); 106 doiPETAnalysis::GetInstance()->ReadOut(blockID,crystalID,interactionTime,timeOfAnnihil,interactionPos,totalEdep); 107 108 crystalID_vec.clear(); 109 posInter_vec.clear(); 110 edepInCry_vec.clear(); 111 } 112 113 } 114 115 //Apply center of mass calculation in determining the interaction position 116 G4ThreeVector doiPETRun::CenterOfMassInteractionPos(const std::vector<G4int> &cryID, const std::vector<G4double> &energyDep, G4double edepTotal, const std::vector<G4ThreeVector> &pos){ 117 posInterInCrystal = G4ThreeVector(); 118 for(std::size_t i=0; i<cryID.size();i++){ 119 posInterInCrystal += pos[i] * energyDep[i]; 120 } 121 return (posInterInCrystal)/edepTotal; 122 } 123 124 125 void doiPETRun::SetAnnihilationTime(G4double t){ 126 timeOfAnnihil = t; 127 } 128 // 129 void doiPETRun::Clear(){ 130 //clear the content of the associative containers, the multimap and the set 131 std::multimap< G4int, InteractionInformation* >::iterator mitr2; 132 for( mitr2= mapBlockInteraction.begin(); mitr2!=mapBlockInteraction.end(); mitr2++){ 133 delete (*mitr2).second; 134 } 135 136 mapBlockInteraction.clear(); 137 setBlockInteraction.clear(); 138 // 139 } 140 // 141 void doiPETRun::GetEventIDRun(G4int evt){ 142 eventID = evt; 143 } 144 // 145 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 146 147 void doiPETRun::Merge(const G4Run* run) 148 { 149 G4Run::Merge(run); 150 } 151