ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
LXeEventAction.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file LXeEventAction.cc
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 //
29 //
30 //
31 #include "LXeEventAction.hh"
32 #include "LXeScintHit.hh"
33 #include "LXePMTHit.hh"
34 #include "LXeTrajectory.hh"
35 #include "LXeRun.hh"
36 #include "LXeHistoManager.hh"
38 
39 #include "G4EventManager.hh"
40 #include "G4SDManager.hh"
41 #include "G4RunManager.hh"
42 #include "G4Event.hh"
43 #include "G4EventManager.hh"
44 #include "G4TrajectoryContainer.hh"
45 #include "G4Trajectory.hh"
46 #include "G4VVisManager.hh"
47 #include "G4ios.hh"
48 #include "G4UImanager.hh"
49 #include "G4SystemOfUnits.hh"
50 
51 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
52 
54  : fDetector(det),fScintCollID(-1),fPMTCollID(-1),fVerbose(0),
55  fPMTThreshold(1),fForcedrawphotons(false),fForcenophotons(false)
56 {
58 
59  fHitCount = 0;
62  fAbsorptionCount = 0;
64  fTotE = 0.0;
65 
66  fConvPosSet = false;
67  fEdepMax = 0.0;
68 
70 }
71 
72 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
73 
75 
76 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
77 
79 
80  fHitCount = 0;
83  fAbsorptionCount = 0;
85  fTotE = 0.0;
86 
87  fConvPosSet = false;
88  fEdepMax = 0.0;
89 
91 
93  if(fScintCollID<0)
94  fScintCollID=SDman->GetCollectionID("scintCollection");
95  if(fPMTCollID<0)
96  fPMTCollID=SDman->GetCollectionID("pmtHitCollection");
97 }
98 
99 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
100 
102 
103  G4TrajectoryContainer* trajectoryContainer=anEvent->GetTrajectoryContainer();
104 
105  G4int n_trajectories = 0;
106  if (trajectoryContainer) n_trajectories = trajectoryContainer->entries();
107 
108  // extract the trajectories and draw them
110  for (G4int i=0; i<n_trajectories; i++){
111  LXeTrajectory* trj = (LXeTrajectory*)
112  ((*(anEvent->GetTrajectoryContainer()))[i]);
113  if(trj->GetParticleName()=="opticalphoton"){
116  }
117  trj->DrawTrajectory();
118  }
119  }
120 
121  LXeScintHitsCollection* scintHC = nullptr;
122  LXePMTHitsCollection* pmtHC = nullptr;
123  G4HCofThisEvent* hitsCE = anEvent->GetHCofThisEvent();
124 
125  //Get the hit collections
126  if(hitsCE){
127  if(fScintCollID>=0) {
128  scintHC = (LXeScintHitsCollection*)(hitsCE->GetHC(fScintCollID));
129  }
130  if(fPMTCollID>=0) {
131  pmtHC = (LXePMTHitsCollection*)(hitsCE->GetHC(fPMTCollID));
132  }
133  }
134 
135  //Hits in scintillator
136  if(scintHC){
137  int n_hit = scintHC->entries();
138  G4ThreeVector eWeightPos(0.);
139  G4double edep;
140  G4double edepMax=0;
141 
142  for(int i=0;i<n_hit;i++){ //gather info on hits in scintillator
143  edep=(*scintHC)[i]->GetEdep();
144  fTotE += edep;
145  eWeightPos += (*scintHC)[i]->GetPos()*edep;//calculate energy weighted pos
146  if(edep>edepMax){
147  edepMax=edep;//store max energy deposit
148  G4ThreeVector posMax=(*scintHC)[i]->GetPos();
149  fPosMax = posMax;
150  fEdepMax = edep;
151  }
152  }
153 
155 
156  if(fTotE == 0.){
157  if(fVerbose>0)G4cout<<"No hits in the scintillator this event."<<G4endl;
158  }
159  else{
160  //Finish calculation of energy weighted position
161  eWeightPos /= fTotE;
162  fEWeightPos = eWeightPos;
163  if(fVerbose>0){
164  G4cout << "\tEnergy weighted position of hits in LXe : "
165  << eWeightPos/mm << G4endl;
166  }
167  }
168  if(fVerbose>0){
169  G4cout << "\tTotal energy deposition in scintillator : "
170  << fTotE / keV << " (keV)" << G4endl;
171  }
172  }
173 
174  if(pmtHC){
175  G4ThreeVector reconPos(0.,0.,0.);
176  G4int pmts=pmtHC->entries();
177  //Gather info from all PMTs
178  for(G4int i=0;i<pmts;i++){
179  fHitCount += (*pmtHC)[i]->GetPhotonCount();
180  reconPos+=(*pmtHC)[i]->GetPMTPos()*(*pmtHC)[i]->GetPhotonCount();
181  if((*pmtHC)[i]->GetPhotonCount()>=fPMTThreshold){
183  }
184  else{//wasnt above the threshold, turn it back off
185  (*pmtHC)[i]->SetDrawit(false);
186  }
187  }
188 
191 
192  if(fHitCount > 0) {//dont bother unless there were hits
193  reconPos/=fHitCount;
194  if(fVerbose>0){
195  G4cout << "\tReconstructed position of hits in LXe : "
196  << reconPos/mm << G4endl;
197  }
198  fReconPos = reconPos;
199  }
200  pmtHC->DrawAllHits();
201  }
202 
207 
208  if(fVerbose>0){
209  //End of event output. later to be controlled by a verbose level
210  G4cout << "\tNumber of photons that hit PMTs in this event : "
211  << fHitCount << G4endl;
212  G4cout << "\tNumber of PMTs above threshold("<<fPMTThreshold<<") : "
214  G4cout << "\tNumber of photons produced by scintillation in this event : "
216  G4cout << "\tNumber of photons produced by cerenkov in this event : "
218  G4cout << "\tNumber of photons absorbed (OpAbsorption) in this event : "
219  << fAbsorptionCount << G4endl;
220  G4cout << "\tNumber of photons absorbed at boundaries (OpBoundary) in "
221  << "this event : " << fBoundaryAbsorptionCount << G4endl;
222  G4cout << "Unaccounted for photons in this event : "
225  << G4endl;
226  }
227 
228  // update the run statistics
229  LXeRun* run = static_cast<LXeRun*>(
231 
232  run->IncHitCount(fHitCount);
235  run->IncEDep(fTotE);
239 
240  //If we have set the flag to save 'special' events, save here
241  if(fPhotonCount_Scint + fPhotonCount_Ceren < fDetector->GetSaveThreshold())
242  {
244  }
245 }
246 
247 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......