EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
EICG4ZDCHitTree.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file EICG4ZDCHitTree.cc
1 // - 1/June/2021 TTree production for ZDC G4 hits Shima Shimizu
2 // - 14/Dec/2021 Track info is added
3 //
4 #include "EICG4ZDCHitTree.h"
5 
6 #include <g4main/PHG4Hit.h>
7 #include <g4main/PHG4HitContainer.h>
8 #include <g4main/PHG4InEvent.h>
9 #include <g4main/PHG4Particle.h>
10 
11 #include <fun4all/Fun4AllHistoManager.h>
12 #include <fun4all/SubsysReco.h> // for SubsysReco
13 
14 #include <phool/getClass.h>
15 
16 #include <TFile.h>
17 #include <TH1.h>
18 #include <TTree.h>
19 
20 #include <sstream>
21 #include <utility>
22 
23 using namespace std;
24 
25 EICG4ZDCHitTree::EICG4ZDCHitTree(const std::string &name, const std::string &filename)
26  : SubsysReco(name)
27  , nblocks(0)
28  , hm(nullptr)
29  , _filename(filename)
30  , tree(nullptr)
31  , outfile(nullptr)
32  , Nhit(0)
33  , Ntrack(0)
34 {
35 }
36 
37 EICG4ZDCHitTree::~EICG4ZDCHitTree(){
38  // delete ntup;
39  delete hm;
40 }
41 
42 int EICG4ZDCHitTree::Init(PHCompositeNode *)
43  {
44  hm = new Fun4AllHistoManager(Name());
45  outfile = new TFile(_filename.c_str(), "RECREATE");
46  tree = new TTree("zdchit", "Collection of EICG4ZDC G4Hits");
47 
48  tree->Branch("Nhit", &Nhit, "Nhit/I");
49  tree->Branch("layerType", &layerType);
50  tree->Branch("layerID", &layerID);
51  tree->Branch("xID", &xID);
52  tree->Branch("yID", &yID);
53  tree->Branch("x0", &x0);
54  tree->Branch("y0", &y0);
55  tree->Branch("z0", &z0);
56  tree->Branch("x1", &x1);
57  tree->Branch("y1", &y1);
58  tree->Branch("z1", &z1);
59  tree->Branch("time0", &time0);
60  tree->Branch("time1", &time1);
61  tree->Branch("edep", &edep);
62 
63  tree->Branch("Ntrack",&Ntrack, "Ntrack/I");
64  tree->Branch("trk_px", &trk_px);
65  tree->Branch("trk_py", &trk_py);
66  tree->Branch("trk_pz", &trk_pz);
67  tree->Branch("trk_e", &trk_e);
68  tree->Branch("trk_pid", &trk_pid);
69 
70  return 0;
71  }
72 
73  int EICG4ZDCHitTree::process_event(PHCompositeNode *topNode)
74  {
75  ostringstream nodename;
76  set<string>::const_iterator iter;
77  for (iter = _node_postfix.begin(); iter != _node_postfix.end(); ++iter)
78  {
79 
80  nodename.str("");
81  nodename << "G4HIT_" << *iter;
82  PHG4HitContainer *hits = findNode::getClass<PHG4HitContainer>(topNode, nodename.str());
83  if (!hits) return 0;
84 
85  Nhit = 0;
86 
87  PHG4HitContainer::ConstRange hit_range = hits->getHits();
88  for (PHG4HitContainer::ConstIterator hit_iter = hit_range.first; hit_iter != hit_range.second; hit_iter++){
89  if(hit_iter->second->get_hit_type()<0) continue;
90  if(hit_iter->second->get_hit_type()>5) continue; //to skip absorber hits
91  Nhit++;
92 
93  layerType.push_back(hit_iter->second->get_hit_type());
94  layerID.push_back(hit_iter->second->get_layer());
95  xID.push_back(hit_iter->second->get_index_i());
96  yID.push_back(hit_iter->second->get_index_j());
97  x0.push_back(hit_iter->second->get_x(0));
98  y0.push_back(hit_iter->second->get_y(0));
99  z0.push_back(hit_iter->second->get_z(0));
100  x1.push_back(hit_iter->second->get_x(1));
101  y1.push_back(hit_iter->second->get_y(1));
102  z1.push_back(hit_iter->second->get_z(1));
103  time0.push_back(hit_iter->second->get_t(0));
104  time1.push_back(hit_iter->second->get_t(1));
105  edep.push_back(hit_iter->second->get_edep());
106  }
107 
108  }
109 
110  PHG4InEvent *track_truthinfo
111  = findNode::getClass<PHG4InEvent>(topNode, "PHG4INEVENT");
112 
113  if(track_truthinfo) {
114 
115  Ntrack=0;
116 
117  pair<multimap<int, PHG4Particle *>::const_iterator, multimap<int, PHG4Particle *>::const_iterator > particlebegin_end = track_truthinfo->GetParticles();
118  multimap<int,PHG4Particle *>::const_iterator particle_iter;
119 
120  for (particle_iter = particlebegin_end.first; particle_iter != particlebegin_end.second; ++particle_iter){
121 
123  const PHG4Particle *truth = particle_iter->second;
124 
125  Ntrack++;
127  trk_px.push_back(truth->get_px());
128  trk_py.push_back(truth->get_py());
129  trk_pz.push_back(truth->get_pz());
130  trk_e.push_back(truth->get_e());
131  trk_pid.push_back(truth->get_pid());
132  }
133  }
134 
135  tree->Fill();
136 
137  layerType.clear();
138  layerID.clear();
139  xID.clear();
140  yID.clear();
141  x0.clear();
142  y0.clear();
143  z0.clear();
144  x1.clear();
145  y1.clear();
146  z1.clear();
147  time0.clear();
148  time1.clear();
149  edep.clear();
150  trk_px.clear();
151  trk_py.clear();
152  trk_pz.clear();
153  trk_e.clear();
154  trk_pid.clear();
155 
156  return 0;
157  }
158 
159  int EICG4ZDCHitTree::End(PHCompositeNode *topNode)
160  {
161  outfile->cd();
162  tree->Write();
163  outfile->Write();
164  outfile->Close();
165  delete outfile;
166  hm->dumpHistos(_filename, "UPDATE");
167  return 0;
168  }
169 
170  void EICG4ZDCHitTree::AddNode(const std::string &name, const int detid)
171  {
172  _node_postfix.insert(name);
173  _detid[name] = detid;
174  return;
175  }