Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer << 3 // * DISCLAIMER * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th << 5 // * The following disclaimer summarizes all the specific disclaimers * 6 // * the Geant4 Collaboration. It is provided << 6 // * of contributors to this software. The specific disclaimers,which * 7 // * conditions of the Geant4 Software License << 7 // * govern, are listed with their locations in: * 8 // * LICENSE and available at http://cern.ch/ << 8 // * http://cern.ch/geant4/license * 9 // * include a list of copyright holders. << 10 // * 9 // * * 11 // * Neither the authors of this software syst 10 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 11 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 12 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 13 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file << 14 // * use. * 16 // * for the full disclaimer and the limitatio << 17 // * 15 // * * 18 // * This code implementation is the result << 16 // * This code implementation is the intellectual property of the * 19 // * technical work of the GEANT4 collaboratio << 17 // * GEANT4 collaboration. * 20 // * By using, copying, modifying or distri << 18 // * By copying, distributing or modifying the Program (or any work * 21 // * any work based on the software) you ag << 19 // * based on the Program) you indicate your acceptance of this * 22 // * use in resulting scientific publicati << 20 // * statement, and all its terms. * 23 // * acceptance of all terms of the Geant4 Sof << 24 // ******************************************* 21 // ******************************************************************** 25 // 22 // 26 // 23 // >> 24 // $Id: GammaRayTelEventAction.cc,v 1.17 2003/06/16 16:46:25 gunter Exp $ >> 25 // GEANT4 tag $Name: geant4-08-00-patch-01 $ 27 // ------------------------------------------- 26 // ------------------------------------------------------------ 28 // GEANT 4 class implementation file 27 // GEANT 4 class implementation file 29 // CERN Geneva Switzerland 28 // CERN Geneva Switzerland 30 // 29 // 31 // 30 // 32 // ------------ GammaRayTelEventAction - 31 // ------------ GammaRayTelEventAction ------ 33 // by R.Giannitrapani, F.Longo & G. 32 // by R.Giannitrapani, F.Longo & G.Santin (13 nov 2000) 34 // 33 // 35 // - inclusion of Digits by F.Longo & R.Gianni 34 // - inclusion of Digits by F.Longo & R.Giannitrapani (24 oct 2001) 36 // 35 // 37 // - Modification of analysis management by G. 36 // - Modification of analysis management by G.Santin (18 Nov 2001) 38 // 37 // 39 // ******************************************* 38 // ************************************************************ 40 39 41 #include "GammaRayTelEventAction.hh" 40 #include "GammaRayTelEventAction.hh" 42 #include "GammaRayTelTrackerHit.hh" 41 #include "GammaRayTelTrackerHit.hh" 43 #include "GammaRayTelAnticoincidenceHit.hh" 42 #include "GammaRayTelAnticoincidenceHit.hh" 44 #include "GammaRayTelCalorimeterHit.hh" 43 #include "GammaRayTelCalorimeterHit.hh" >> 44 >> 45 #ifdef G4ANALYSIS_USE 45 #include "GammaRayTelAnalysis.hh" 46 #include "GammaRayTelAnalysis.hh" 46 #include "GammaRayTelDigi.hh" << 47 #endif 47 #include "GammaRayTelDigitizer.hh" << 48 48 49 #include "G4DigiManager.hh" << 50 #include "G4Event.hh" 49 #include "G4Event.hh" 51 #include "G4EventManager.hh" 50 #include "G4EventManager.hh" 52 #include "G4HCofThisEvent.hh" 51 #include "G4HCofThisEvent.hh" 53 #include "G4ios.hh" << 52 #include "G4VHitsCollection.hh" >> 53 #include "G4TrajectoryContainer.hh" >> 54 #include "G4Trajectory.hh" >> 55 #include "G4VVisManager.hh" 54 #include "G4SDManager.hh" 56 #include "G4SDManager.hh" 55 #include "G4SystemOfUnits.hh" << 56 #include "G4UImanager.hh" 57 #include "G4UImanager.hh" >> 58 #include "G4ios.hh" 57 #include "G4UnitsTable.hh" 59 #include "G4UnitsTable.hh" 58 #include "G4VHitsCollection.hh" << 59 #include "Randomize.hh" 60 #include "Randomize.hh" 60 61 61 //....oooOO0OOooo........oooOO0OOooo........oo << 62 #include "GammaRayTelDigi.hh" >> 63 #include "GammaRayTelDigitizer.hh" >> 64 #include "G4DigiManager.hh" 62 65 63 // This file is a global variable in which we 66 // This file is a global variable in which we store energy deposition per hit 64 // and other relevant information 67 // and other relevant information 65 68 66 GammaRayTelEventAction::GammaRayTelEventAction << 69 #ifdef G4STORE_DATA 67 auto *digitizer = new GammaRayTelDigitizer << 70 extern std::ofstream outFile; 68 G4DigiManager::GetDMpointer()->AddNewModul << 71 #endif 69 } << 70 72 71 //....oooOO0OOooo........oooOO0OOooo........oo 73 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 72 74 73 GammaRayTelEventAction::~GammaRayTelEventActio << 75 GammaRayTelEventAction::GammaRayTelEventAction() >> 76 :trackerCollID(-1),calorimeterCollID(-1), >> 77 anticoincidenceCollID(-1), drawFlag("all") >> 78 { >> 79 G4DigiManager * fDM = G4DigiManager::GetDMpointer(); >> 80 GammaRayTelDigitizer * myDM = new GammaRayTelDigitizer( "TrackerDigitizer" ); >> 81 fDM->AddNewModule(myDM); 74 } 82 } 75 83 76 //....oooOO0OOooo........oooOO0OOooo........oo 84 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 77 85 78 void GammaRayTelEventAction::BeginOfEventActio << 86 GammaRayTelEventAction::~GammaRayTelEventAction() 79 G4int eventIdentifer = event->GetEventID() << 87 { 80 G4cout << "Event: " << eventIdentifer << G << 81 auto *sensitiveDetectorManager = G4SDManag << 82 << 83 if (trackerCollectionID == -1) { << 84 trackerCollectionID = sensitiveDetecto << 85 } << 86 if (anticoincidenceCollectionID == -1) { << 87 anticoincidenceCollectionID = sensitiv << 88 } << 89 if (calorimeterCollectionID == -1) { << 90 calorimeterCollectionID = sensitiveDet << 91 } << 92 } 88 } 93 89 94 //....oooOO0OOooo........oooOO0OOooo........oo 90 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 95 91 96 void GammaRayTelEventAction::EndOfEventAction( << 92 void GammaRayTelEventAction::BeginOfEventAction(const G4Event* evt) 97 G4int eventIdentifier = event->GetEventID( << 93 { 98 << 99 if (theRunAction == nullptr) { << 100 G4Exception("GammaRayTelEventAction::B << 101 } << 102 94 103 #ifdef G4STORE_DATA << 95 G4int evtNb = evt->GetEventID(); 104 auto *outputFile = theRunAction->GetOutput << 96 G4cout << "Event: " << evtNb << G4endl; 105 #endif << 97 G4SDManager * SDman = G4SDManager::GetSDMpointer(); >> 98 >> 99 if (trackerCollID==-1) { >> 100 trackerCollID = SDman->GetCollectionID("TrackerCollection"); >> 101 } >> 102 if(anticoincidenceCollID==-1) { >> 103 anticoincidenceCollID = >> 104 SDman->GetCollectionID("AnticoincidenceCollection"); >> 105 } >> 106 if(calorimeterCollID==-1) { >> 107 calorimeterCollID = >> 108 SDman->GetCollectionID("CalorimeterCollection"); >> 109 } >> 110 } 106 111 107 auto *collection = event->GetHCofThisEvent << 112 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 108 GammaRayTelTrackerHitsCollection *trackerC << 109 // GammaRayTelCalorimeterHitsCollection* c << 110 // GammaRayTelAnticoincidenceHitsCollectio << 111 << 112 auto *fDM = G4DigiManager::GetDMpointer(); << 113 GammaRayTelAnalysis *analysis = GammaRayTe << 114 << 115 if (collection != nullptr) { << 116 trackerCollection = (GammaRayTelTracke << 117 // calorimeterCollection = (GammaRayTe << 118 // anticoincidenceCollection = (GammaR << 119 << 120 if (trackerCollection != nullptr) { << 121 G4int numberOfHits = trackerCollec << 122 G4cout << "Number of tracker hits << 123 << 124 G4double depositedEnergy{0}; << 125 G4int stripNumber; << 126 G4int planeNumber; << 127 G4int isXPlane; << 128 << 129 // This is a cycle on all the trac << 130 << 131 for (auto i = 0; i < numberOfHits; << 132 // Here we put the hit data in << 133 depositedEnergy = (*trackerCol << 134 stripNumber = (*trackerCollect << 135 planeNumber = (*trackerCollect << 136 isXPlane = (*trackerCollection << 137 113 >> 114 void GammaRayTelEventAction::EndOfEventAction(const G4Event* evt) >> 115 { >> 116 G4int event_id = evt->GetEventID(); >> 117 >> 118 G4TrajectoryContainer * trajectoryContainer = evt->GetTrajectoryContainer(); >> 119 G4int n_trajectories = 0; >> 120 if (trajectoryContainer) n_trajectories = trajectoryContainer->entries(); >> 121 >> 122 >> 123 G4HCofThisEvent* HCE = evt->GetHCofThisEvent(); >> 124 GammaRayTelTrackerHitsCollection* THC = 0; >> 125 GammaRayTelCalorimeterHitsCollection* CHC = 0; >> 126 GammaRayTelAnticoincidenceHitsCollection* AHC = 0; >> 127 >> 128 >> 129 G4DigiManager * fDM = G4DigiManager::GetDMpointer(); >> 130 >> 131 if (HCE) >> 132 { >> 133 THC = (GammaRayTelTrackerHitsCollection*)(HCE->GetHC(trackerCollID)); >> 134 CHC = (GammaRayTelCalorimeterHitsCollection*) >> 135 (HCE->GetHC(calorimeterCollID)); >> 136 AHC = (GammaRayTelAnticoincidenceHitsCollection*) >> 137 (HCE->GetHC(anticoincidenceCollID)); >> 138 } >> 139 >> 140 if (THC) >> 141 { >> 142 int n_hit = THC->entries(); >> 143 G4cout << "Number of tracker hits in this event = " << n_hit << G4endl; >> 144 G4double ESil=0; >> 145 G4int NStrip, NPlane, IsX; >> 146 >> 147 // This is a cycle on all the tracker hits of this event >> 148 >> 149 for (int i=0;i<n_hit;i++) >> 150 { >> 151 // Here we put the hit data in a an ASCII file for >> 152 // later analysis >> 153 ESil = (*THC)[i]->GetEdepSil(); >> 154 NStrip = (*THC)[i]->GetNStrip(); >> 155 NPlane = (*THC)[i]->GetNSilPlane(); >> 156 IsX = (*THC)[i]->GetPlaneType(); >> 157 138 #ifdef G4STORE_DATA 158 #ifdef G4STORE_DATA 139 (*outputFile) << std::setw(7) << 159 outFile << std::setw(7) << event_id << " " << 140 << " " << depositedEnergy/ << 160 ESil/keV << " " << NStrip << 141 << " " << stripNumber << 161 " " << NPlane << " " << IsX << " " << 142 << " " << planeNumber << 162 (*THC)[i]->GetPos().x()/mm <<" "<< 143 << " " << isXPlane << 163 (*THC)[i]->GetPos().y()/mm <<" "<< 144 << " " << (*trackerCollect << 164 (*THC)[i]->GetPos().z()/mm <<" "<< 145 << " " << (*trackerCollect << 165 G4endl; 146 << " " << (*trackerCollect << 147 << " " << G4endl; << 148 #else 166 #else 149 G4cout << std::setw(7) << even << 167 G4cout << std::setw(7) << event_id << " " << 150 << " " << depositedEnergy << 168 ESil/keV << " " << NStrip << 151 << " " << stripNumber << 169 " " << NPlane << " " << IsX << " " << 152 << " " << planeNumber << 170 (*THC)[i]->GetPos().x()/mm <<" "<< 153 << " " << isXPlane << 171 (*THC)[i]->GetPos().y()/mm <<" "<< 154 << " " << (*trackerCollect << 172 (*THC)[i]->GetPos().z()/mm <<" "<< 155 << " " << (*trackerCollect << 173 G4endl; 156 << " " << (*trackerCollect << 174 #endif 157 << " " << G4endl; << 158 #endif << 159 << 160 // Here we fill the histograms << 161 if (isXPlane != 0) { << 162 if (analysis->GetHisto2DMo << 163 analysis->InsertPositi << 164 } else { << 165 analysis->InsertPositi << 166 } << 167 << 168 if (planeNumber == 0) { << 169 analysis->InsertEnergy << 170 } << 171 analysis->InsertHits(plane << 172 } else { << 173 if (analysis->GetHisto2DMo << 174 analysis->InsertPositi << 175 } else { << 176 analysis->InsertPositi << 177 } << 178 if (planeNumber == 0) { << 179 analysis->InsertEnergy << 180 } << 181 analysis->InsertHits(plane << 182 } << 183 analysis->setNtuple(depositedE << 184 (*trackerCollection)[i]->G << 185 (*trackerCollection)[i]->G << 186 (*trackerCollection)[i]->G << 187 ); << 188 } << 189 analysis->EndOfEvent(numberOfHits) << 190 } << 191 175 192 auto *myDM = (GammaRayTelDigitizer*) f << 176 #ifdef G4ANALYSIS_USE 193 myDM->Digitize(); << 194 177 195 #ifdef G4STORE_DATA << 178 // Here we fill the histograms of the Analysis manager 196 // The whole block is needed only when << 179 GammaRayTelAnalysis* analysis = GammaRayTelAnalysis::getInstance(); 197 // protect block to avoid compilations << 198 180 199 auto digitsCollectionIdentifier = fDM->Get << 181 if(IsX) 200 // G4cout << "Digits collection: " << << 182 { 201 << 183 if (analysis->GetHisto2DMode()=="position") 202 auto *digitsCollection = (GammaRayTelDigit << 184 analysis->InsertPositionXZ((*THC)[i]->GetPos().x()/mm,(*THC)[i]->GetPos().z()/mm); >> 185 else >> 186 analysis->InsertPositionXZ(NStrip, NPlane); >> 187 if (NPlane == 0) analysis->InsertEnergy(ESil/keV); >> 188 analysis->InsertHits(NPlane); >> 189 } >> 190 else >> 191 { >> 192 if (analysis->GetHisto2DMode()=="position") >> 193 analysis->InsertPositionYZ((*THC)[i]->GetPos().y()/mm,(*THC)[i]->GetPos().z()/mm); >> 194 else >> 195 analysis->InsertPositionYZ(NStrip, NPlane); >> 196 if (NPlane == 0) analysis->InsertEnergy(ESil/keV); >> 197 analysis->InsertHits(NPlane); >> 198 } >> 199 >> 200 #ifdef G4ANALYSIS_USE >> 201 analysis->setNtuple( ESil/keV, NPlane, (*THC)[i]->GetPos().x()/mm, >> 202 (*THC)[i]->GetPos().y()/mm, >> 203 (*THC)[i]->GetPos().z()/mm); >> 204 #endif >> 205 >> 206 #endif >> 207 >> 208 } >> 209 // Here we call the analysis manager function for visualization >> 210 #ifdef G4ANALYSIS_USE >> 211 GammaRayTelAnalysis* analysis = GammaRayTelAnalysis::getInstance(); >> 212 analysis->EndOfEvent(n_hit); >> 213 #endif >> 214 } >> 215 >> 216 GammaRayTelDigitizer * myDM = >> 217 (GammaRayTelDigitizer*)fDM->FindDigitizerModule( "TrackerDigitizer" ); >> 218 myDM->Digitize(); >> 219 >> 220 G4int myDigiCollID = fDM->GetDigiCollectionID("DigitsCollection"); >> 221 >> 222 // G4cout << "digi collecion" << myDigiCollID << G4endl; >> 223 >> 224 GammaRayTelDigitsCollection * DC = (GammaRayTelDigitsCollection*)fDM->GetDigiCollection( myDigiCollID ); >> 225 >> 226 if(DC) { >> 227 // G4cout << "Total Digits " << DC->entries() << G4endl; >> 228 G4int n_digi = DC->entries(); >> 229 G4int NStrip, NPlane, IsX; >> 230 for (G4int i=0;i<n_digi;i++) { >> 231 // Here we put the digi data in a an ASCII file for >> 232 // later analysis >> 233 NStrip = (*DC)[i]->GetStripNumber(); >> 234 NPlane = (*DC)[i]->GetPlaneNumber(); >> 235 IsX = (*DC)[i]->GetPlaneType(); 203 236 204 if (digitsCollection != nullptr) { << 237 // outFile << std::setw(7) << event_id << " " << NStrip << 205 auto numberOfDigits = digitsCollec << 238 // " " << NPlane << " " << IsX << " " << G4endl; 206 // G4cout << "Total number of digi << 239 } 207 << 240 } 208 G4int stripNumber; << 241 209 G4int planeNumber; << 242 if(G4VVisManager::GetConcreteInstance()) { 210 G4int isXPlane; << 243 for(G4int i=0; i<n_trajectories; i++) { 211 << 244 G4Trajectory* trj = (G4Trajectory *)((*(evt->GetTrajectoryContainer()))[i]); 212 for (auto i = 0; i < numberOfDigit << 245 if (drawFlag == "all") trj->DrawTrajectory(50); 213 // Here we put the digi data i << 246 else if ((drawFlag == "charged")&&(trj->GetCharge() != 0.)) 214 stripNumber = (*digitsCollecti << 247 trj->DrawTrajectory(50); 215 planeNumber = (*digitsCollecti << 216 isXPlane = (*digitsCollection) << 217 << 218 (*outputFile) << std::setw(7) << 219 << eventIdentifier << 220 << " " << stripNumber << 221 << " " << planeNumber << 222 << " " << isXPlane << 223 << " " << G4endl; << 224 } << 225 } << 226 #endif << 227 } 248 } >> 249 } 228 } 250 } >> 251 >> 252 >> 253 >> 254 >> 255 >> 256 >> 257 >> 258 >> 259 >> 260 >> 261 >> 262 >> 263 229 264