Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** 25 // 25 // 26 /// \file electromagnetic/TestEm11/src/doiPETR 26 /// \file electromagnetic/TestEm11/src/doiPETRun.cc 27 /// \brief Implementation of the doiPETRun cla 27 /// \brief Implementation of the doiPETRun class 28 // 28 // 29 // $Id: doiPETRun.cc 71376 2013-06-14 07:44:50 << 30 // 29 // 31 //....oooOO0OOooo........oooOO0OOooo........oo 30 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 32 //....oooOO0OOooo........oooOO0OOooo........oo 31 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 33 32 34 #include "doiPETRun.hh" 33 #include "doiPETRun.hh" 35 #include "doiPETDetectorConstruction.hh" 34 #include "doiPETDetectorConstruction.hh" 36 #include "doiPETPrimaryGeneratorAction.hh" 35 #include "doiPETPrimaryGeneratorAction.hh" 37 36 38 #include "G4UnitsTable.hh" 37 #include "G4UnitsTable.hh" 39 #include "G4SystemOfUnits.hh" 38 #include "G4SystemOfUnits.hh" 40 #include "G4Event.hh" 39 #include "G4Event.hh" 41 #include "doiPETAnalysis.hh" 40 #include "doiPETAnalysis.hh" 42 //....oooOO0OOooo........oooOO0OOooo........oo 41 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 43 42 44 doiPETRun::doiPETRun() 43 doiPETRun::doiPETRun() 45 : G4Run() 44 : G4Run() 46 { 45 { 47 // 46 // 48 totalTime = 0 * s; 47 totalTime = 0 * s; 49 prev_totalTime = 0 * s; 48 prev_totalTime = 0 * s; 50 totalEdep = 0*keV; 49 totalEdep = 0*keV; 51 } 50 } 52 51 53 //....oooOO0OOooo........oooOO0OOooo........oo 52 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 54 53 55 doiPETRun::~doiPETRun() 54 doiPETRun::~doiPETRun() 56 {} 55 {} 57 56 58 // 57 // 59 void doiPETRun::GetIntractionInfomation(Intera 58 void doiPETRun::GetIntractionInfomation(InteractionInformation*step) 60 { 59 { 61 //The photon can have many interactions in o 60 //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 associat 61 //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< G 62 mapBlockInteraction.insert( std::multimap< G4int, InteractionInformation* >::value_type( step->GetBlockNo(), step) ); 64 63 65 //The following stores the block number uniq 64 //The following stores the block number uniquely. Save the interacting block in the associative container without saving multiple copies 66 setBlockInteraction.insert( step->GetBlockNo 65 setBlockInteraction.insert( step->GetBlockNo()); 67 } 66 } 68 67 69 ///////// InteractingCrystal ///////////////// 68 ///////// InteractingCrystal /////////////////////////////////// 70 69 71 void doiPETRun::FindInteractingCrystal() 70 void doiPETRun::FindInteractingCrystal() 72 { 71 { >> 72 >> 73 //G4cout<<"totalTime: "<<totalTime/s<<" "<<InitialActivity/becquerel<<G4endl; 73 std::multimap< G4int, InteractionInformation 74 std::multimap< G4int, InteractionInformation* >::iterator mitr; 74 std::set<G4int>::iterator sitr; 75 std::set<G4int>::iterator sitr; 75 G4int cystalIDTemp; 76 G4int cystalIDTemp; 76 77 77 for(sitr=setBlockInteraction.begin(); sitr!= 78 for(sitr=setBlockInteraction.begin(); sitr!=setBlockInteraction.end(); sitr++){ 78 blockID = (*sitr); 79 blockID = (*sitr); 79 mitr = mapBlockInteraction.find(blockID); 80 mitr = mapBlockInteraction.find(blockID); 80 numberofInteractions = mapBlockInteraction 81 numberofInteractions = mapBlockInteraction.count(blockID); //number of interactions in which the interaction deposite energy in a block. 81 totalEdep=edepMax=0.; 82 totalEdep=edepMax=0.; 82 for(G4int i=0; i<numberofInteractions; i++ 83 for(G4int i=0; i<numberofInteractions; i++, mitr++){ 83 edep = (*mitr).second->GetEdep(); 84 edep = (*mitr).second->GetEdep(); 84 cystalIDTemp = (*mitr).second->GetCrysta 85 cystalIDTemp = (*mitr).second->GetCrystalNo(); >> 86 //edep_AfterCrystalBlurring = doiPETAnalysis::GetInstance()->QuantumEffifciency(edep); 85 edep_AfterCrystalBlurring = doiPETAnalys 87 edep_AfterCrystalBlurring = doiPETAnalysis::GetInstance()->QuantumEffifciency(edep, blockID, cystalIDTemp); 86 totalEdep += edep_AfterCrystalBlurring; 88 totalEdep += edep_AfterCrystalBlurring; 87 << 88 crystalID_vec.push_back(cystalIDTemp); << 89 posInter_vec.push_back(((*mitr).second-> << 90 edepInCry_vec.push_back(edep_AfterCrysta << 91 << 92 //This is to identify the crystal with 89 //This is to identify the crystal with the highest energy deposition. 93 if(edepMax<edep_AfterCrystalBlurring){ 90 if(edepMax<edep_AfterCrystalBlurring){ 94 edepMax=edep_AfterCrystalBlurring; 91 edepMax=edep_AfterCrystalBlurring; 95 crystalID = cystalIDTemp; 92 crystalID = cystalIDTemp; 96 interactionTime = (*mitr).second->GetG 93 interactionTime = (*mitr).second->GetGlobalTime(); 97 << 94 interactionPos = (*mitr).second->GetInteractionPositionInCrystal(); 98 //position of interaction based on hig << 99 //interactionPos = (*mitr).second->Get << 100 } 95 } 101 } 96 } 102 if(totalEdep==0)continue; << 103 << 104 //The interaction position is based on ene << 105 interactionPos = CenterOfMassInteractionPo << 106 doiPETAnalysis::GetInstance()->ReadOut(blo 97 doiPETAnalysis::GetInstance()->ReadOut(blockID,crystalID,interactionTime,timeOfAnnihil,interactionPos,totalEdep); 107 << 108 crystalID_vec.clear(); << 109 posInter_vec.clear(); << 110 edepInCry_vec.clear(); << 111 } 98 } 112 << 113 } 99 } 114 << 115 //Apply center of mass calculation in determin << 116 G4ThreeVector doiPETRun::CenterOfMassInteracti << 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 100 125 void doiPETRun::SetAnnihilationTime(G4double t 101 void doiPETRun::SetAnnihilationTime(G4double t){ 126 timeOfAnnihil = t; 102 timeOfAnnihil = t; 127 } 103 } 128 // 104 // 129 void doiPETRun::Clear(){ 105 void doiPETRun::Clear(){ 130 //clear the content of the associative conta 106 //clear the content of the associative containers, the multimap and the set 131 std::multimap< G4int, InteractionInformation 107 std::multimap< G4int, InteractionInformation* >::iterator mitr2; 132 for( mitr2= mapBlockInteraction.begin(); mit 108 for( mitr2= mapBlockInteraction.begin(); mitr2!=mapBlockInteraction.end(); mitr2++){ 133 delete (*mitr2).second; 109 delete (*mitr2).second; 134 } 110 } 135 111 136 mapBlockInteraction.clear(); 112 mapBlockInteraction.clear(); 137 setBlockInteraction.clear(); 113 setBlockInteraction.clear(); 138 // 114 // 139 } 115 } 140 // 116 // 141 void doiPETRun::GetEventIDRun(G4int evt){ 117 void doiPETRun::GetEventIDRun(G4int evt){ 142 eventID = evt; 118 eventID = evt; 143 } 119 } 144 // 120 // 145 //....oooOO0OOooo........oooOO0OOooo........oo 121 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 146 122 147 void doiPETRun::Merge(const G4Run* run) 123 void doiPETRun::Merge(const G4Run* run) 148 { 124 { 149 G4Run::Merge(run); 125 G4Run::Merge(run); 150 } 126 } 151 127