Geant4  10.02.p01
B5EventAction.cc
Go to the documentation of this file.
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 // $Id: B5EventAction.cc 94486 2015-11-19 08:33:37Z gcosmo $
27 //
30 
31 #include "B5EventAction.hh"
32 #include "B5HodoscopeHit.hh"
33 #include "B5DriftChamberHit.hh"
34 #include "B5EmCalorimeterHit.hh"
35 #include "B5HadCalorimeterHit.hh"
36 #include "B5Analysis.hh"
37 
38 #include "G4Event.hh"
39 #include "G4RunManager.hh"
40 #include "G4EventManager.hh"
41 #include "G4HCofThisEvent.hh"
42 #include "G4VHitsCollection.hh"
43 #include "G4SDManager.hh"
44 #include "G4SystemOfUnits.hh"
45 #include "G4ios.hh"
46 
47 
48 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
49 
52  fHHC1ID(-1),
53  fHHC2ID(-1),
54  fDHC1ID(-1),
55  fDHC2ID(-1),
56  fECHCID(-1),
57  fHCHCID(-1),
58  fEmCalEdep(),
59  fHadCalEdep()
60 {
61  // set printing per each event
63 
64  // initialize the vectors
65  fEmCalEdep.resize(80, 0.);
66  fHadCalEdep.resize(20, 0.);
67 }
68 
69 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
70 
72 {}
73 
74 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
75 
77 {
78  if (fHHC1ID==-1) {
80  fHHC1ID = sdManager->GetCollectionID("hodoscope1/hodoscopeColl");
81  fHHC2ID = sdManager->GetCollectionID("hodoscope2/hodoscopeColl");
82  fDHC1ID = sdManager->GetCollectionID("chamber1/driftChamberColl");
83  fDHC2ID = sdManager->GetCollectionID("chamber2/driftChamberColl");
84  fECHCID = sdManager->GetCollectionID("EMcalorimeter/EMcalorimeterColl");
85  fHCHCID = sdManager->GetCollectionID("HadCalorimeter/HadCalorimeterColl");
86  }
87 }
88 
89 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
90 
92 {
93  G4HCofThisEvent* hce = event->GetHCofThisEvent();
94  if (!hce)
95  {
97  msg << "No hits collection of this event found." << G4endl;
98  G4Exception("B5EventAction::EndOfEventAction()",
99  "B5Code001", JustWarning, msg);
100  return;
101  }
102 
103 
104  // Get hits collections
106  = static_cast<B5HodoscopeHitsCollection*>(hce->GetHC(fHHC1ID));
107 
109  = static_cast<B5HodoscopeHitsCollection*>(hce->GetHC(fHHC2ID));
110 
112  = static_cast<B5DriftChamberHitsCollection*>(hce->GetHC(fDHC1ID));
113 
115  = static_cast<B5DriftChamberHitsCollection*>(hce->GetHC(fDHC2ID));
116 
118  = static_cast<B5EmCalorimeterHitsCollection*>(hce->GetHC(fECHCID));
119 
121  = static_cast<B5HadCalorimeterHitsCollection*>(hce->GetHC(fHCHCID));
122 
123  if ( (!hHC1) || (!hHC2) || (!dHC1) || (!dHC2) || (!ecHC) || (!hcHC) )
124  {
126  msg << "Some of hits collections of this event not found." << G4endl;
127  G4Exception("B5EventAction::EndOfEventAction()",
128  "B5Code001", JustWarning, msg);
129  return;
130  }
131 
132  //
133  // Fill histograms & ntuple
134  //
135 
136  // Get analysis manager
137  G4AnalysisManager* analysisManager = G4AnalysisManager::Instance();
138 
139  // Fill histograms
140 
141  G4int n_hit = dHC1->entries();
142  analysisManager->FillH1(0, n_hit);
143 
144  for (G4int i=0;i<n_hit;i++)
145  {
146  B5DriftChamberHit* hit = (*dHC1)[i];
147  G4ThreeVector localPos = hit->GetLocalPos();
148  analysisManager->FillH2(0, localPos.x(), localPos.y());
149  }
150 
151  n_hit = dHC2->entries();
152  analysisManager->FillH1(1, n_hit);
153 
154  for (G4int i=0;i<n_hit;i++)
155  {
156  B5DriftChamberHit* hit = (*dHC2)[i];
157  G4ThreeVector localPos = hit->GetLocalPos();
158  analysisManager->FillH2(1, localPos.x(), localPos.y());
159  }
160 
161 
162  // Fill ntuple
163 
164  // Dc1Hits
165  analysisManager->FillNtupleIColumn(0, dHC1->entries());
166  // Dc2Hits
167  analysisManager->FillNtupleIColumn(1, dHC1->entries());
168 
169  // ECEnergy
170  G4int totalEmHit = 0;
171  G4double totalEmE = 0.;
172  for (G4int i=0;i<80;i++)
173  {
174  B5EmCalorimeterHit* hit = (*ecHC)[i];
175  G4double eDep = hit->GetEdep();
176  if (eDep>0.)
177  {
178  totalEmHit++;
179  totalEmE += eDep;
180  }
181  fEmCalEdep[i] = eDep;
182  }
183  analysisManager->FillNtupleDColumn(2, totalEmE);
184 
185  // HCEnergy
186  G4int totalHadHit = 0;
187  G4double totalHadE = 0.;
188  for (G4int i=0;i<20;i++)
189  {
190  B5HadCalorimeterHit* hit = (*hcHC)[i];
191  G4double eDep = hit->GetEdep();
192  if (eDep>0.)
193  {
194  totalHadHit++;
195  totalHadE += eDep;
196  }
197  fHadCalEdep[i] = eDep;
198  }
199  analysisManager->FillNtupleDColumn(3, totalHadE);
200 
201  // Time 1
202  for (G4int i=0;i<hHC1->entries();i++)
203  {
204  analysisManager->FillNtupleDColumn(4,(*hHC1)[i]->GetTime());
205  }
206 
207  // Time 2
208  for (G4int i=0;i<hHC2->entries();i++)
209  {
210  analysisManager->FillNtupleDColumn(5,(*hHC2)[i]->GetTime());
211  }
212 
213  analysisManager->AddNtupleRow();
214 
215  //
216  // Print diagnostics
217  //
218 
220  if ( printModulo==0 || event->GetEventID() % printModulo != 0) return;
221 
222  G4PrimaryParticle* primary = event->GetPrimaryVertex(0)->GetPrimary(0);
223  G4cout << G4endl
224  << ">>> Event " << event->GetEventID() << " >>> Simulation truth : "
225  << primary->GetG4code()->GetParticleName()
226  << " " << primary->GetMomentum() << G4endl;
227 
228  // Hodoscope 1
229  n_hit = hHC1->entries();
230  G4cout << "Hodoscope 1 has " << n_hit << " hits." << G4endl;
231  for (G4int i=0;i<n_hit;i++)
232  {
233  B5HodoscopeHit* hit = (*hHC1)[i];
234  hit->Print();
235  }
236 
237  // Hodoscope 2
238  n_hit = hHC2->entries();
239  G4cout << "Hodoscope 2 has " << n_hit << " hits." << G4endl;
240  for (G4int i=0;i<n_hit;i++)
241  {
242  B5HodoscopeHit* hit = (*hHC2)[i];
243  hit->Print();
244  }
245 
246  // Drift chamber 1
247  n_hit = dHC1->entries();
248  G4cout << "Drift Chamber 1 has " << n_hit << " hits." << G4endl;
249  for (G4int i2=0;i2<5;i2++)
250  {
251  for (G4int i=0;i<n_hit;i++)
252  {
253  B5DriftChamberHit* hit = (*dHC1)[i];
254  if (hit->GetLayerID()==i2) hit->Print();
255  }
256  }
257 
258  // Drift chamber 2
259  n_hit = dHC2->entries();
260  G4cout << "Drift Chamber 2 has " << n_hit << " hits." << G4endl;
261  for (G4int i2=0;i2<5;i2++)
262  {
263  for (G4int i=0;i<n_hit;i++)
264  {
265  B5DriftChamberHit* hit = (*dHC2)[i];
266  if (hit->GetLayerID()==i2) hit->Print();
267  }
268  }
269 
270  // EM calorimeter
271  G4cout << "EM Calorimeter has " << totalEmHit << " hits. Total Edep is "
272  << totalEmE/MeV << " (MeV)" << G4endl;
273 
274  // Had calorimeter
275  G4cout << "Hadron Calorimeter has " << totalHadHit << " hits. Total Edep is "
276  << totalHadE/MeV << " (MeV)" << G4endl;
277 }
278 
279 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
static const double MeV
Definition: G4SIunits.hh:211
G4int GetPrintProgress()
G4ThreeVector GetMomentum() const
Definition of the B5EventAction class.
G4VHitsCollection * GetHC(G4int i)
G4int GetCollectionID(G4String colName)
Definition: G4SDManager.cc:131
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
CLHEP::Hep3Vector G4ThreeVector
void SetPrintProgress(G4int i)
G4int entries() const
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
G4ParticleDefinition * GetG4code() const
G4int GetEventID() const
Definition: G4Event.hh:151
Definition of the B5EmCalorimeterHit class.
G4ThreeVector GetLocalPos() const
G4bool FillNtupleIColumn(G4int id, G4int value)
G4GLOB_DLL std::ostream G4cout
G4double GetEdep() const
G4bool FillNtupleDColumn(G4int id, G4double value)
Drift chamber hit.
EM Calorimeter hit.
virtual void Print()
Hodoscope hit.
Definition of the B5HodoscopeHit class.
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
Selection of the analysis technology.
static G4RunManager * GetRunManager()
Definition: G4RunManager.cc:79
Hadron Calorimeter hit.
G4bool FillH2(G4int id, G4double xvalue, G4double yvalue, G4double weight=1.0)
G4bool FillH1(G4int id, G4double value, G4double weight=1.0)
static G4SDManager * GetSDMpointer()
Definition: G4SDManager.cc:40
virtual void BeginOfEventAction(const G4Event *)
virtual ~B5EventAction()
G4double GetEdep() const
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
G4int GetLayerID() const
virtual void EndOfEventAction(const G4Event *)
std::vector< G4double > fEmCalEdep
Definition of the B5DriftChamberHit class.
Definition of the B5HadCalorimeterHit class.
std::vector< G4double > fHadCalEdep