source: trunk/examples/advanced/gammaray_telescope/src/GammaRayTelEventAction.cc@ 1302

Last change on this file since 1302 was 1230, checked in by garnier, 16 years ago

update to geant4.9.3

File size: 8.5 KB
Line 
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// $Id: GammaRayTelEventAction.cc,v 1.19 2006/06/29 15:56:39 gunter Exp $
28// GEANT4 tag $Name: geant4-09-03-cand-01 $
29// ------------------------------------------------------------
30// GEANT 4 class implementation file
31// CERN Geneva Switzerland
32//
33//
34// ------------ GammaRayTelEventAction ------
35// by R.Giannitrapani, F.Longo & G.Santin (13 nov 2000)
36//
37// - inclusion of Digits by F.Longo & R.Giannitrapani (24 oct 2001)
38//
39// - Modification of analysis management by G.Santin (18 Nov 2001)
40//
41// ************************************************************
42
43#include "GammaRayTelEventAction.hh"
44#include "GammaRayTelTrackerHit.hh"
45#include "GammaRayTelAnticoincidenceHit.hh"
46#include "GammaRayTelCalorimeterHit.hh"
47
48#ifdef G4ANALYSIS_USE
49#include "GammaRayTelAnalysis.hh"
50#endif
51
52#include "G4Event.hh"
53#include "G4EventManager.hh"
54#include "G4HCofThisEvent.hh"
55#include "G4VHitsCollection.hh"
56#include "G4TrajectoryContainer.hh"
57#include "G4Trajectory.hh"
58#include "G4VVisManager.hh"
59#include "G4SDManager.hh"
60#include "G4UImanager.hh"
61#include "G4ios.hh"
62#include "G4UnitsTable.hh"
63#include "Randomize.hh"
64
65#include "GammaRayTelDigi.hh"
66#include "GammaRayTelDigitizer.hh"
67#include "G4DigiManager.hh"
68
69// This file is a global variable in which we store energy deposition per hit
70// and other relevant information
71
72#ifdef G4STORE_DATA
73extern std::ofstream outFile;
74#endif
75
76//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
77
78GammaRayTelEventAction::GammaRayTelEventAction()
79 :trackerCollID(-1),calorimeterCollID(-1),
80 anticoincidenceCollID(-1), drawFlag("all")
81{
82 G4DigiManager * fDM = G4DigiManager::GetDMpointer();
83 GammaRayTelDigitizer * myDM = new GammaRayTelDigitizer( "GammaRayTelDigitizer" );
84 fDM->AddNewModule(myDM);
85}
86
87//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
88
89GammaRayTelEventAction::~GammaRayTelEventAction()
90{
91}
92
93//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
94
95void GammaRayTelEventAction::BeginOfEventAction(const G4Event* evt)
96{
97
98 G4int evtNb = evt->GetEventID();
99 G4cout << "Event: " << evtNb << G4endl;
100 G4SDManager * SDman = G4SDManager::GetSDMpointer();
101
102 if (trackerCollID==-1) {
103 trackerCollID = SDman->GetCollectionID("TrackerCollection");
104 }
105 if(anticoincidenceCollID==-1) {
106 anticoincidenceCollID =
107 SDman->GetCollectionID("AnticoincidenceCollection");
108 }
109 if(calorimeterCollID==-1) {
110 calorimeterCollID =
111 SDman->GetCollectionID("CalorimeterCollection");
112 }
113}
114
115//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
116
117void GammaRayTelEventAction::EndOfEventAction(const G4Event* evt)
118{
119 G4int event_id = evt->GetEventID();
120
121 G4TrajectoryContainer * trajectoryContainer = evt->GetTrajectoryContainer();
122 G4int n_trajectories = 0;
123 if (trajectoryContainer) n_trajectories = trajectoryContainer->entries();
124
125
126 G4HCofThisEvent* HCE = evt->GetHCofThisEvent();
127 GammaRayTelTrackerHitsCollection* THC = 0;
128 GammaRayTelCalorimeterHitsCollection* CHC = 0;
129 GammaRayTelAnticoincidenceHitsCollection* AHC = 0;
130
131
132 G4DigiManager * fDM = G4DigiManager::GetDMpointer();
133
134 if (HCE)
135 {
136 THC = (GammaRayTelTrackerHitsCollection*)(HCE->GetHC(trackerCollID));
137 CHC = (GammaRayTelCalorimeterHitsCollection*)
138 (HCE->GetHC(calorimeterCollID));
139 AHC = (GammaRayTelAnticoincidenceHitsCollection*)
140 (HCE->GetHC(anticoincidenceCollID));
141
142 if (THC)
143 {
144 int n_hit = THC->entries();
145 G4cout << "Number of tracker hits in this event = " << n_hit << G4endl;
146 G4double ESil=0;
147 G4int NStrip, NPlane, IsX;
148
149 // This is a cycle on all the tracker hits of this event
150
151 for (int i=0;i<n_hit;i++)
152 {
153 // Here we put the hit data in a an ASCII file for
154 // later analysis
155 ESil = (*THC)[i]->GetEdepSil();
156 NStrip = (*THC)[i]->GetNStrip();
157 NPlane = (*THC)[i]->GetNSilPlane();
158 IsX = (*THC)[i]->GetPlaneType();
159
160#ifdef G4STORE_DATA
161 outFile << std::setw(7) << event_id << " " <<
162 ESil/keV << " " << NStrip <<
163 " " << NPlane << " " << IsX << " " <<
164 (*THC)[i]->GetPos().x()/mm <<" "<<
165 (*THC)[i]->GetPos().y()/mm <<" "<<
166 (*THC)[i]->GetPos().z()/mm <<" "<<
167 G4endl;
168#else
169 G4cout << std::setw(7) << event_id << " " <<
170 ESil/keV << " " << NStrip <<
171 " " << NPlane << " " << IsX << " " <<
172 (*THC)[i]->GetPos().x()/mm <<" "<<
173 (*THC)[i]->GetPos().y()/mm <<" "<<
174 (*THC)[i]->GetPos().z()/mm <<" "<<
175 G4endl;
176#endif
177
178#ifdef G4ANALYSIS_USE
179
180 // Here we fill the histograms of the Analysis manager
181 GammaRayTelAnalysis* analysis = GammaRayTelAnalysis::getInstance();
182
183 if(IsX)
184 {
185 if (analysis->GetHisto2DMode()=="position")
186 analysis->InsertPositionXZ((*THC)[i]->GetPos().x()/mm,(*THC)[i]->GetPos().z()/mm);
187 else
188 analysis->InsertPositionXZ(NStrip, NPlane);
189 if (NPlane == 0) analysis->InsertEnergy(ESil/keV);
190 analysis->InsertHits(NPlane);
191 }
192 else
193 {
194 if (analysis->GetHisto2DMode()=="position")
195 analysis->InsertPositionYZ((*THC)[i]->GetPos().y()/mm,(*THC)[i]->GetPos().z()/mm);
196 else
197 analysis->InsertPositionYZ(NStrip, NPlane);
198 if (NPlane == 0) analysis->InsertEnergy(ESil/keV);
199 analysis->InsertHits(NPlane);
200 }
201
202#ifdef G4ANALYSIS_USE
203 analysis->setNtuple( ESil/keV, NPlane, (*THC)[i]->GetPos().x()/mm,
204 (*THC)[i]->GetPos().y()/mm,
205 (*THC)[i]->GetPos().z()/mm);
206#endif
207
208#endif
209
210 }
211 // Here we call the analysis manager function for visualization
212#ifdef G4ANALYSIS_USE
213 GammaRayTelAnalysis* analysis = GammaRayTelAnalysis::getInstance();
214 analysis->EndOfEvent(n_hit);
215#endif
216
217 }
218
219 GammaRayTelDigitizer * myDM =
220 (GammaRayTelDigitizer*)fDM->FindDigitizerModule( "GammaRayTelDigitizer" );
221 myDM->Digitize();
222
223 G4int myDigiCollID = fDM->GetDigiCollectionID("DigitsCollection");
224
225 // G4cout << "digi collecion" << myDigiCollID << G4endl;
226
227 GammaRayTelDigitsCollection * DC = (GammaRayTelDigitsCollection*)fDM->GetDigiCollection( myDigiCollID );
228
229 if(DC) {
230 // G4cout << "Total Digits " << DC->entries() << G4endl;
231 G4int n_digi = DC->entries();
232 G4int NStrip, NPlane, IsX;
233 for (G4int i=0;i<n_digi;i++) {
234 // Here we put the digi data in a an ASCII file for
235 // later analysis
236 NStrip = (*DC)[i]->GetStripNumber();
237 NPlane = (*DC)[i]->GetPlaneNumber();
238 IsX = (*DC)[i]->GetPlaneType();
239
240 // outFile << std::setw(7) << event_id << " " << NStrip <<
241 // " " << NPlane << " " << IsX << " " << G4endl;
242 }
243 }
244 }
245
246 if(G4VVisManager::GetConcreteInstance()) {
247 for(G4int i=0; i<n_trajectories; i++) {
248 G4Trajectory* trj = (G4Trajectory *)((*(evt->GetTrajectoryContainer()))[i]);
249 if (drawFlag == "all") trj->DrawTrajectory(50);
250 else if ((drawFlag == "charged")&&(trj->GetCharge() != 0.))
251 trj->DrawTrajectory(50);
252 }
253 }
254}
255
256
257
258
259
260
261
262
263
264
265
266
267
Note: See TracBrowser for help on using the repository browser.