ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
JetTruthEval.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file JetTruthEval.cc
1 #include "JetTruthEval.h"
2 
3 #include "CaloTruthEval.h"
4 #include "SvtxTruthEval.h"
5 
6 #include <g4jets/Jet.h>
7 #include <g4jets/JetMap.h>
8 
10 
11 #include <phool/getClass.h>
12 #include <phool/phool.h>
13 
14 #include <cassert>
15 #include <cstdlib>
16 #include <iostream>
17 #include <map>
18 #include <utility>
19 
20 using namespace std;
21 
23  const std::string& truthjetname)
24  : _truthjetname(truthjetname)
25  , _svtxevalstack(topNode)
26  , _cemcevalstack(topNode, "CEMC")
27  , _hcalinevalstack(topNode, "HCALIN")
28  , _hcaloutevalstack(topNode, "HCALOUT")
29  , _femcevalstack(topNode, "FEMC")
30  , _fhcalevalstack(topNode, "FHCAL")
31  , _eemcevalstack(topNode, "EEMC")
32  , _truthinfo(nullptr)
33  , _truthjets(nullptr)
34  , _strict(false)
35  , _verbosity(1)
36  , _errors(0)
37  , _do_cache(true)
38  , _cache_all_truth_particles()
39  , _cache_all_truth_showers()
40  , _cache_all_truth_hits()
41  , _cache_get_truth_jet()
42 {
43  get_node_pointers(topNode);
44 }
45 
47 {
48  if (_verbosity > 0)
49  {
50  if ((_errors > 0) || (_verbosity > 1))
51  {
52  cout << "JetTruthEval::~JetTruthEval() - Error Count: " << _errors << endl;
53  }
54  }
55 }
56 
58 {
59  _svtxevalstack.next_event(topNode);
60  _cemcevalstack.next_event(topNode);
63  _femcevalstack.next_event(topNode);
64  _fhcalevalstack.next_event(topNode);
65 
68  _cache_all_truth_hits.clear();
69  _cache_get_truth_jet.clear();
70 
71  get_node_pointers(topNode);
72 }
73 
74 std::set<PHG4Particle*> JetTruthEval::all_truth_particles(Jet* truthjet)
75 {
76  if (_strict)
77  {
78  assert(truthjet);
79  }
80  else if (!truthjet)
81  {
82  ++_errors;
83  return std::set<PHG4Particle*>();
84  }
85 
86  if (_do_cache)
87  {
88  std::map<Jet*, std::set<PHG4Particle*> >::iterator iter =
89  _cache_all_truth_particles.find(truthjet);
90  if (iter != _cache_all_truth_particles.end())
91  {
92  return iter->second;
93  }
94  }
95 
96  std::set<PHG4Particle*> truth_particles;
97 
98  // loop over all the entries in the truthjet
99  for (Jet::ConstIter iter = truthjet->begin_comp();
100  iter != truthjet->end_comp();
101  ++iter)
102  {
103  Jet::SRC source = iter->first;
104  unsigned int index = iter->second;
105  if (source != Jet::PARTICLE)
106  {
107  cout << PHWHERE << " truth jet contains something other than particles!" << endl;
108  exit(-1);
109  }
110 
111  PHG4Particle* truth_particle = _truthinfo->GetParticle(index);
112 
113  if (_strict)
114  {
115  assert(truth_particle);
116  }
117  else if (!truth_particle)
118  {
119  ++_errors;
120  continue;
121  }
122 
123  truth_particles.insert(truth_particle);
124  }
125 
126  if (_do_cache) _cache_all_truth_particles.insert(make_pair(truthjet, truth_particles));
127 
128  return truth_particles;
129 }
130 
131 std::set<PHG4Shower*> JetTruthEval::all_truth_showers(Jet* truthjet)
132 {
133  if (_strict)
134  {
135  assert(truthjet);
136  }
137  else if (!truthjet)
138  {
139  ++_errors;
140  return std::set<PHG4Shower*>();
141  }
142 
143  if (_do_cache)
144  {
145  std::map<Jet*, std::set<PHG4Shower*> >::iterator iter =
146  _cache_all_truth_showers.find(truthjet);
147  if (iter != _cache_all_truth_showers.end())
148  {
149  return iter->second;
150  }
151  }
152 
153  std::set<PHG4Shower*> truth_showers;
154 
155  std::set<PHG4Particle*> truth_particles = all_truth_particles(truthjet);
156 
157  // loop over all the entries in the truthjet
158  for (std::set<PHG4Particle*>::iterator iter = truth_particles.begin();
159  iter != truth_particles.end();
160  ++iter)
161  {
162  PHG4Particle* particle = *iter;
163 
164  if (_strict)
165  {
166  assert(particle);
167  }
168  else if (!particle)
169  {
170  ++_errors;
171  continue;
172  }
173 
174  // any calo truth eval module would work here...
175  CaloTruthEval* cemc_truth_eval = _cemcevalstack.get_truth_eval();
176  PHG4Shower* shower = cemc_truth_eval->get_primary_shower(particle);
177  if (shower)
178  {
179  truth_showers.insert(shower);
180  }
181  }
182 
183  if (_do_cache) _cache_all_truth_showers.insert(make_pair(truthjet, truth_showers));
184 
185  return truth_showers;
186 }
187 
188 std::set<PHG4Hit*> JetTruthEval::all_truth_hits(Jet* truthjet)
189 {
190  if (_strict)
191  {
192  assert(truthjet);
193  }
194  else if (!truthjet)
195  {
196  ++_errors;
197  return std::set<PHG4Hit*>();
198  }
199 
200  if (_do_cache)
201  {
202  std::map<Jet*, std::set<PHG4Hit*> >::iterator iter =
203  _cache_all_truth_hits.find(truthjet);
204  if (iter != _cache_all_truth_hits.end())
205  {
206  return iter->second;
207  }
208  }
209 
210  std::set<PHG4Hit*> truth_hits;
211 
212  std::set<PHG4Particle*> truth_particles = all_truth_particles(truthjet);
213 
214  // loop over all the entries in the truthjet
215  for (std::set<PHG4Particle*>::iterator iter = truth_particles.begin();
216  iter != truth_particles.end();
217  ++iter)
218  {
219  PHG4Particle* particle = *iter;
220 
221  if (_strict)
222  {
223  assert(particle);
224  }
225  else if (!particle)
226  {
227  ++_errors;
228  continue;
229  }
230 
231  // ask the svtx truth eval to backtrack the particles to g4hits
232  SvtxTruthEval* svtx_truth_eval = _svtxevalstack.get_truth_eval();
233  std::set<PHG4Hit*> svtx_g4hits = svtx_truth_eval->all_truth_hits(particle);
234 
235  for (std::set<PHG4Hit*>::iterator jter = svtx_g4hits.begin();
236  jter != svtx_g4hits.end();
237  ++jter)
238  {
239  PHG4Hit* g4hit = *jter;
240 
241  if (_strict)
242  {
243  assert(g4hit);
244  }
245  else if (!g4hit)
246  {
247  ++_errors;
248  continue;
249  }
250 
251  truth_hits.insert(g4hit);
252  }
253  }
254 
255  if (_do_cache) _cache_all_truth_hits.insert(make_pair(truthjet, truth_hits));
256 
257  return truth_hits;
258 }
259 
261 {
262  if (_strict)
263  {
264  assert(particle);
265  }
266  else if (!particle)
267  {
268  ++_errors;
269  return nullptr;
270  }
271 
272  if (_do_cache)
273  {
274  std::map<PHG4Particle*, Jet*>::iterator iter =
275  _cache_get_truth_jet.find(particle);
276  if (iter != _cache_get_truth_jet.end())
277  {
278  return iter->second;
279  }
280  }
281 
282  Jet* truth_jet = nullptr;
283 
284  // loop over all jets and look for this particle...
285  for (JetMap::Iter iter = _truthjets->begin();
286  iter != _truthjets->end();
287  ++iter)
288  {
289  Jet* candidate = iter->second;
290 
291  // loop over all consituents and look for this particle
292  for (Jet::ConstIter jter = candidate->begin_comp();
293  jter != candidate->end_comp();
294  ++jter)
295  {
296  unsigned int index = jter->second;
297 
298  PHG4Particle* constituent = _truthinfo->GetParticle(index);
299  if (_strict)
300  {
301  assert(constituent);
302  }
303  else if (!constituent)
304  {
305  ++_errors;
306  continue;
307  }
308 
309  if (get_svtx_eval_stack()->get_truth_eval()->are_same_particle(constituent, particle))
310  {
311  truth_jet = candidate;
312  break;
313  }
314  }
315 
316  if (truth_jet) break;
317  }
318 
319  if (_do_cache) _cache_get_truth_jet.insert(make_pair(particle, truth_jet));
320 
321  return truth_jet;
322 }
323 
325 {
326  _truthinfo = findNode::getClass<PHG4TruthInfoContainer>(topNode, "G4TruthInfo");
327  if (!_truthinfo)
328  {
329  cerr << PHWHERE << " ERROR: Can't find G4TruthInfo" << endl;
330  exit(-1);
331  }
332 
333  _truthjets = findNode::getClass<JetMap>(topNode, _truthjetname.c_str());
334  if (!_truthjets)
335  {
336  cerr << PHWHERE << " ERROR: Can't find " << _truthjetname << endl;
337  exit(-1);
338  }
339 
340  return;
341 }