| 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 // 26 //
>> 27 // $Id: GammaRayTelEventAction.cc,v 1.20 2010-06-06 06:18:54 perl Exp $
>> 28 // GEANT4 tag $Name: geant4-09-04-patch-01 $
27 // ------------------------------------------- 29 // ------------------------------------------------------------
28 // GEANT 4 class implementation file 30 // GEANT 4 class implementation file
29 // CERN Geneva Switzerland 31 // CERN Geneva Switzerland
30 // 32 //
31 // 33 //
32 // ------------ GammaRayTelEventAction - 34 // ------------ GammaRayTelEventAction ------
33 // by R.Giannitrapani, F.Longo & G. 35 // by R.Giannitrapani, F.Longo & G.Santin (13 nov 2000)
34 // 36 //
35 // - inclusion of Digits by F.Longo & R.Gianni 37 // - inclusion of Digits by F.Longo & R.Giannitrapani (24 oct 2001)
36 // 38 //
37 // - Modification of analysis management by G. 39 // - Modification of analysis management by G.Santin (18 Nov 2001)
38 // 40 //
39 // ******************************************* 41 // ************************************************************
40 42
41 #include "GammaRayTelEventAction.hh" 43 #include "GammaRayTelEventAction.hh"
42 #include "GammaRayTelTrackerHit.hh" 44 #include "GammaRayTelTrackerHit.hh"
43 #include "GammaRayTelAnticoincidenceHit.hh" 45 #include "GammaRayTelAnticoincidenceHit.hh"
44 #include "GammaRayTelCalorimeterHit.hh" 46 #include "GammaRayTelCalorimeterHit.hh"
>> 47
>> 48 #ifdef G4ANALYSIS_USE
45 #include "GammaRayTelAnalysis.hh" 49 #include "GammaRayTelAnalysis.hh"
46 #include "GammaRayTelDigi.hh" << 50 #endif
47 #include "GammaRayTelDigitizer.hh" <<
48 51
49 #include "G4DigiManager.hh" <<
50 #include "G4Event.hh" 52 #include "G4Event.hh"
51 #include "G4EventManager.hh" 53 #include "G4EventManager.hh"
52 #include "G4HCofThisEvent.hh" 54 #include "G4HCofThisEvent.hh"
53 #include "G4ios.hh" << 55 #include "G4VHitsCollection.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( "GammaRayTelDigitizer" );
>> 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 94
99 if (theRunAction == nullptr) { << 95 G4int evtNb = evt->GetEventID();
100 G4Exception("GammaRayTelEventAction::B << 96 G4cout << "Event: " << evtNb << G4endl;
101 } << 97 G4SDManager * SDman = G4SDManager::GetSDMpointer();
102 << 98
103 #ifdef G4STORE_DATA << 99 if (trackerCollID==-1) {
104 auto *outputFile = theRunAction->GetOutput << 100 trackerCollID = SDman->GetCollectionID("TrackerCollection");
105 #endif << 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 if (THC)
>> 140 {
>> 141 int n_hit = THC->entries();
>> 142 G4cout << "Number of tracker hits in this event = " << n_hit << G4endl;
>> 143 G4double ESil=0;
>> 144 G4int NStrip, NPlane, IsX;
>> 145
>> 146 // This is a cycle on all the tracker hits of this event
>> 147
>> 148 for (int i=0;i<n_hit;i++)
>> 149 {
>> 150 // Here we put the hit data in a an ASCII file for
>> 151 // later analysis
>> 152 ESil = (*THC)[i]->GetEdepSil();
>> 153 NStrip = (*THC)[i]->GetNStrip();
>> 154 NPlane = (*THC)[i]->GetNSilPlane();
>> 155 IsX = (*THC)[i]->GetPlaneType();
>> 156
138 #ifdef G4STORE_DATA 157 #ifdef G4STORE_DATA
139 (*outputFile) << std::setw(7) << 158 outFile << std::setw(7) << event_id << " " <<
140 << " " << depositedEnergy/ << 159 ESil/keV << " " << NStrip <<
141 << " " << stripNumber << 160 " " << NPlane << " " << IsX << " " <<
142 << " " << planeNumber << 161 (*THC)[i]->GetPos().x()/mm <<" "<<
143 << " " << isXPlane << 162 (*THC)[i]->GetPos().y()/mm <<" "<<
144 << " " << (*trackerCollect << 163 (*THC)[i]->GetPos().z()/mm <<" "<<
145 << " " << (*trackerCollect << 164 G4endl;
146 << " " << (*trackerCollect <<
147 << " " << G4endl; <<
148 #else 165 #else
149 G4cout << std::setw(7) << even << 166 G4cout << std::setw(7) << event_id << " " <<
150 << " " << depositedEnergy << 167 ESil/keV << " " << NStrip <<
151 << " " << stripNumber << 168 " " << NPlane << " " << IsX << " " <<
152 << " " << planeNumber << 169 (*THC)[i]->GetPos().x()/mm <<" "<<
153 << " " << isXPlane << 170 (*THC)[i]->GetPos().y()/mm <<" "<<
154 << " " << (*trackerCollect << 171 (*THC)[i]->GetPos().z()/mm <<" "<<
155 << " " << (*trackerCollect << 172 G4endl;
156 << " " << (*trackerCollect << 173 #endif
157 << " " << G4endl; << 174
158 #endif << 175 #ifdef G4ANALYSIS_USE
159 << 176
160 // Here we fill the histograms << 177 // Here we fill the histograms of the Analysis manager
161 if (isXPlane != 0) { << 178 GammaRayTelAnalysis* analysis = GammaRayTelAnalysis::getInstance();
162 if (analysis->GetHisto2DMo << 179
163 analysis->InsertPositi << 180 if(IsX)
164 } else { << 181 {
165 analysis->InsertPositi << 182 if (analysis->GetHisto2DMode()=="position")
166 } << 183 analysis->InsertPositionXZ((*THC)[i]->GetPos().x()/mm,(*THC)[i]->GetPos().z()/mm);
167 << 184 else
168 if (planeNumber == 0) { << 185 analysis->InsertPositionXZ(NStrip, NPlane);
169 analysis->InsertEnergy << 186 if (NPlane == 0) analysis->InsertEnergy(ESil/keV);
170 } << 187 analysis->InsertHits(NPlane);
171 analysis->InsertHits(plane << 188 }
172 } else { << 189 else
173 if (analysis->GetHisto2DMo << 190 {
174 analysis->InsertPositi << 191 if (analysis->GetHisto2DMode()=="position")
175 } else { << 192 analysis->InsertPositionYZ((*THC)[i]->GetPos().y()/mm,(*THC)[i]->GetPos().z()/mm);
176 analysis->InsertPositi << 193 else
177 } << 194 analysis->InsertPositionYZ(NStrip, NPlane);
178 if (planeNumber == 0) { << 195 if (NPlane == 0) analysis->InsertEnergy(ESil/keV);
179 analysis->InsertEnergy << 196 analysis->InsertHits(NPlane);
180 } << 197 }
181 analysis->InsertHits(plane << 198
182 } << 199 #ifdef G4ANALYSIS_USE
183 analysis->setNtuple(depositedE << 200 analysis->setNtuple( ESil/keV, NPlane, (*THC)[i]->GetPos().x()/mm,
184 (*trackerCollection)[i]->G << 201 (*THC)[i]->GetPos().y()/mm,
185 (*trackerCollection)[i]->G << 202 (*THC)[i]->GetPos().z()/mm);
186 (*trackerCollection)[i]->G << 203 #endif
187 ); << 204
188 } << 205 #endif
189 analysis->EndOfEvent(numberOfHits) << 206
190 } << 207 }
191 << 208 // Here we call the analysis manager function for visualization
192 auto *myDM = (GammaRayTelDigitizer*) f << 209 #ifdef G4ANALYSIS_USE
193 myDM->Digitize(); << 210 GammaRayTelAnalysis* analysis = GammaRayTelAnalysis::getInstance();
194 << 211 analysis->EndOfEvent(n_hit);
195 #ifdef G4STORE_DATA << 212 #endif
196 // The whole block is needed only when << 213
197 // protect block to avoid compilations << 214 }
198 << 215
199 auto digitsCollectionIdentifier = fDM->Get << 216 GammaRayTelDigitizer * myDM =
200 // G4cout << "Digits collection: " << << 217 (GammaRayTelDigitizer*)fDM->FindDigitizerModule( "GammaRayTelDigitizer" );
>> 218 myDM->Digitize();
201 219
202 auto *digitsCollection = (GammaRayTelDigit << 220 G4int myDigiCollID = fDM->GetDigiCollectionID("DigitsCollection");
203 221
204 if (digitsCollection != nullptr) { << 222 // G4cout << "digi collecion" << myDigiCollID << G4endl;
205 auto numberOfDigits = digitsCollec << 223
206 // G4cout << "Total number of digi << 224 GammaRayTelDigitsCollection * DC = (GammaRayTelDigitsCollection*)fDM->GetDigiCollection( myDigiCollID );
207 << 225
208 G4int stripNumber; << 226 if(DC) {
209 G4int planeNumber; << 227 // G4cout << "Total Digits " << DC->entries() << G4endl;
210 G4int isXPlane; << 228 G4int n_digi = DC->entries();
211 << 229 G4int NStrip, NPlane, IsX;
212 for (auto i = 0; i < numberOfDigit << 230 for (G4int i=0;i<n_digi;i++) {
213 // Here we put the digi data i << 231 // Here we put the digi data in a an ASCII file for
214 stripNumber = (*digitsCollecti << 232 // later analysis
215 planeNumber = (*digitsCollecti << 233 NStrip = (*DC)[i]->GetStripNumber();
216 isXPlane = (*digitsCollection) << 234 NPlane = (*DC)[i]->GetPlaneNumber();
217 << 235 IsX = (*DC)[i]->GetPlaneType();
218 (*outputFile) << std::setw(7) << 236
219 << eventIdentifier << 237 // outFile << std::setw(7) << event_id << " " << NStrip <<
220 << " " << stripNumber << 238 // " " << NPlane << " " << IsX << " " << G4endl;
221 << " " << planeNumber << 239 }
222 << " " << isXPlane << 240 }
223 << " " << G4endl; <<
224 } <<
225 } <<
226 #endif <<
227 } 241 }
228 } 242 }
>> 243
>> 244
>> 245
>> 246
>> 247
>> 248
>> 249
>> 250
>> 251
>> 252
>> 253
>> 254
>> 255
229 256