EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
QAG4SimulationEicCalorimeterSum.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file QAG4SimulationEicCalorimeterSum.cc
1 #include "QAG4SimulationEicCalorimeterSum.h"
2 
3 #include <qa_modules/QAHistManagerDef.h>
4 
5 #include <g4eval/CaloEvalStack.h>
6 #include <g4eval/CaloRawClusterEval.h>
7 #include <g4eval/SvtxEvalStack.h>
8 
9 #include <g4main/PHG4Particle.h>
10 #include <g4main/PHG4TruthInfoContainer.h>
11 
12 #include <calobase/RawCluster.h>
13 #include <calobase/RawTower.h>
14 #include <calobase/RawTowerContainer.h>
15 #include <calobase/RawTowerGeomContainer.h>
16 
17 #include <trackbase_historic/SvtxTrack.h>
18 
19 #include <g4eval/SvtxTrackEval.h> // for SvtxTrackEval
20 
21 #include <fun4all/Fun4AllHistoManager.h>
22 #include <fun4all/Fun4AllReturnCodes.h>
23 #include <fun4all/SubsysReco.h>
24 
25 #include <phool/getClass.h>
26 
27 #include <TAxis.h>
28 #include <TH1.h>
29 #include <TH2.h>
30 #include <TNamed.h>
31 #include <TString.h>
32 
33 #include <cassert>
34 #include <cmath>
35 #include <iostream>
36 #include <iterator> // for reverse_iterator
37 #include <utility> // for pair
38 #include <vector>
39 
40 using namespace std;
41 
42 QAG4SimulationEicCalorimeterSum::QAG4SimulationEicCalorimeterSum(
43  QAG4SimulationEicCalorimeterSum::enu_flags flags)
44  : SubsysReco("QAG4SimulationEicCalorimeterSum")
45  , _flags(flags)
46  , m_TrackNodeName("TrackMap")
47  , _calo_name_cemc("CEMC")
48  , _calo_name_hcalin("HCALIN")
49  , _calo_name_hcalout("HCALOUT")
50  , _truth_container(nullptr)
51  , _magField(+1.4)
52 {
53 }
54 
55 int QAG4SimulationEicCalorimeterSum::InitRun(PHCompositeNode *topNode)
56 {
57  _truth_container = findNode::getClass<PHG4TruthInfoContainer>(topNode,
58  "G4TruthInfo");
59  if (!_truth_container)
60  {
61  cout << "QAG4SimulationEicCalorimeterSum::InitRun - Fatal Error - "
62  << "unable to find DST node "
63  << "G4TruthInfo" << endl;
64  assert(_truth_container);
65  }
66 
67  if (flag(kProcessCluster))
68  {
69  if (!_caloevalstack_cemc)
70  {
71  _caloevalstack_cemc.reset(new CaloEvalStack(topNode, _calo_name_cemc));
72  _caloevalstack_cemc->set_strict(true);
73  _caloevalstack_cemc->set_verbosity(Verbosity() + 1);
74  }
75  if (!_caloevalstack_hcalin)
76  {
77  _caloevalstack_hcalin.reset(new CaloEvalStack(topNode, _calo_name_hcalin));
78  _caloevalstack_hcalin->set_strict(true);
79  _caloevalstack_hcalin->set_verbosity(Verbosity() + 1);
80  }
81  if (!_caloevalstack_hcalout)
82  {
83  _caloevalstack_hcalout.reset(new CaloEvalStack(topNode, _calo_name_hcalout));
84  _caloevalstack_hcalout->set_strict(true);
85  _caloevalstack_hcalout->set_verbosity(Verbosity() + 1);
86  }
87  }
88 
89  if (flag(kProcessTrackProj))
90  {
91  if (!_svtxevalstack)
92  {
93  _svtxevalstack.reset(new SvtxEvalStack(topNode));
94  _svtxevalstack->set_track_nodename(m_TrackNodeName);
95  _svtxevalstack->set_strict(true);
96  _svtxevalstack->set_verbosity(Verbosity() + 1);
97  }
98  }
99  return Fun4AllReturnCodes::EVENT_OK;
100 }
101 
102 int QAG4SimulationEicCalorimeterSum::Init(PHCompositeNode *topNode)
103 {
104  Fun4AllHistoManager *hm = QAHistManagerDef::getHistoManager();
105  assert(hm);
106  TH1D *h = new TH1D(TString(get_histo_prefix()) + "Normalization", //
107  TString(get_histo_prefix()) + " Normalization;Items;Count", 10, .5, 10.5);
108  int i = 1;
109  h->GetXaxis()->SetBinLabel(i++, "Event");
110  h->GetXaxis()->SetBinLabel(i++, (_calo_name_cemc + " Tower").c_str());
111  h->GetXaxis()->SetBinLabel(i++, (_calo_name_hcalin + " Tower").c_str());
112  h->GetXaxis()->SetBinLabel(i++, (_calo_name_hcalout + " Tower").c_str());
113  h->GetXaxis()->SetBinLabel(i++, (_calo_name_cemc + " Cluster").c_str());
114  h->GetXaxis()->SetBinLabel(i++, (_calo_name_hcalin + " Cluster").c_str());
115  h->GetXaxis()->SetBinLabel(i++, (_calo_name_hcalout + " Cluster").c_str());
116  h->GetXaxis()->SetBinLabel(i++, "Track");
117  h->GetXaxis()->LabelsOption("v");
118  hm->registerHisto(h);
119 
120  // if (flag(kProcessTower))
121  // {
122  // if (Verbosity() >= 1)
123  // cout << "QAG4SimulationEicCalorimeterSum::Init - Process tower occupancies"
124  // << endl;
125  // Init_Tower(topNode);
126  // }
127  if (flag(kProcessCluster))
128  {
129  if (Verbosity() >= 1)
130  cout << "QAG4SimulationEicCalorimeterSum::Init - Process tower occupancies"
131  << endl;
132  Init_Cluster(topNode);
133  }
134 
135  if (flag(kProcessTrackProj))
136  {
137  if (Verbosity() >= 1)
138  cout << "QAG4SimulationEicCalorimeterSum::Init - Process sampling fraction"
139  << endl;
140  Init_TrackProj(topNode);
141  }
142  return Fun4AllReturnCodes::EVENT_OK;
143 }
144 
145 int QAG4SimulationEicCalorimeterSum::process_event(PHCompositeNode *topNode)
146 {
147  if (Verbosity() > 2)
148  cout << "QAG4SimulationEicCalorimeterSum::process_event() entered" << endl;
149 
150  if (_caloevalstack_cemc)
151  _caloevalstack_cemc->next_event(topNode);
152  if (_caloevalstack_hcalin)
153  _caloevalstack_hcalin->next_event(topNode);
154  if (_caloevalstack_hcalout)
155  _caloevalstack_hcalout->next_event(topNode);
156  if (_svtxevalstack)
157  _svtxevalstack->next_event(topNode);
158 
159  // if (flag(kProcessTower))
160  // {
161  // int ret = process_event_Tower(topNode);
162  //
163  // if (ret != Fun4AllReturnCodes::EVENT_OK)
164  // return ret;
165  // }
166 
167  if (flag(kProcessCluster))
168  {
169  int ret = process_event_Cluster(topNode);
170 
171  if (ret != Fun4AllReturnCodes::EVENT_OK)
172  return ret;
173  }
174 
175  if (flag(kProcessTrackProj))
176  {
177  int ret = process_event_TrackProj(topNode);
178 
179  if (ret != Fun4AllReturnCodes::EVENT_OK)
180  return ret;
181  }
182 
183  // at the end, count success events
184  Fun4AllHistoManager *hm = QAHistManagerDef::getHistoManager();
185  assert(hm);
186  TH1D *h_norm = dynamic_cast<TH1D *>(hm->getHisto(
187  get_histo_prefix() + "Normalization"));
188  assert(h_norm);
189  h_norm->Fill("Event", 1);
190 
191  return Fun4AllReturnCodes::EVENT_OK;
192 }
193 
194 string
195 QAG4SimulationEicCalorimeterSum::get_histo_prefix()
196 {
197  return "h_QAG4Sim_CalorimeterSum_";
198 }
199 
200 PHG4Particle *
201 QAG4SimulationEicCalorimeterSum::get_truth_particle()
202 {
203  // get last primary
204  assert(_truth_container);
205  assert(not _truth_container->GetMap().empty());
206  PHG4Particle *last_primary = _truth_container->GetMap().rbegin()->second;
207  assert(last_primary);
208 
209  return last_primary;
210 }
211 
212 int QAG4SimulationEicCalorimeterSum::Init_TrackProj(PHCompositeNode *topNode)
213 {
214  Fun4AllHistoManager *hm = QAHistManagerDef::getHistoManager();
215  assert(hm);
216 
217  hm->registerHisto(
218  new TH2F(
219  TString(get_histo_prefix()) + TString(_calo_name_cemc.c_str()) + "_TrackProj", //
220  TString(_calo_name_cemc.c_str()) + " Tower Energy Distr. around Track Proj.;Polar distance / Tower width;Azimuthal distance / Tower width",
221  (Max_N_Tower - 1) * 10, -Max_N_Tower / 2, Max_N_Tower / 2,
222  (Max_N_Tower - 1) * 10, -Max_N_Tower / 2, Max_N_Tower / 2));
223 
224  hm->registerHisto(
225  new TH2F(
226  TString(get_histo_prefix()) + TString(_calo_name_hcalin.c_str()) + "_TrackProj", //
227  TString(_calo_name_hcalin.c_str()) + " Tower Energy Distr. around Track Proj.;Polar distance / Tower width;Azimuthal distance / Tower width",
228  (Max_N_Tower - 1) * 10, -Max_N_Tower / 2, Max_N_Tower / 2,
229  (Max_N_Tower - 1) * 10, -Max_N_Tower / 2, Max_N_Tower / 2));
230 
231  hm->registerHisto(
232  new TH2F(
233  TString(get_histo_prefix()) + TString(_calo_name_hcalout.c_str()) + "_TrackProj", //
234  TString(_calo_name_hcalout.c_str()) + " Tower Energy Distr. around Track Proj.;Polar distance / Tower width;Azimuthal distance / Tower width",
235  (Max_N_Tower - 1) * 10, -Max_N_Tower / 2, Max_N_Tower / 2,
236  (Max_N_Tower - 1) * 10, -Max_N_Tower / 2, Max_N_Tower / 2));
237 
238  hm->registerHisto(
239  new TH1F(TString(get_histo_prefix()) + "TrackProj_3x3Tower_EP", //
240  "Tower 3x3 sum /E_{Truth};#Sigma_{3x3}[E_{Tower}] / total truth energy",
241  150, 0, 1.5));
242 
243  hm->registerHisto(
244  new TH1F(TString(get_histo_prefix()) + "TrackProj_5x5Tower_EP", //
245  "Tower 5x5 sum /E_{Truth};#Sigma_{5x5}[E_{Tower}] / total truth energy",
246  150, 0, 1.5));
247 
248  return Fun4AllReturnCodes::EVENT_OK;
249 }
250 
251 int QAG4SimulationEicCalorimeterSum::process_event_TrackProj(PHCompositeNode *topNode)
252 {
253  if (Verbosity() > 2)
254  cout << "QAG4SimulationEicCalorimeterSum::process_event_TrackProj() entered"
255  << endl;
256 
257  PHG4Particle *primary = get_truth_particle();
258  if (!primary)
259  return Fun4AllReturnCodes::DISCARDEVENT;
260 
261  SvtxTrackEval *trackeval = _svtxevalstack->get_track_eval();
262  assert(trackeval);
263  SvtxTrack *track = trackeval->best_track_from(primary);
264  if (!track)
265  return Fun4AllReturnCodes::EVENT_OK; // not through the whole event for missing track.
266 
267  Fun4AllHistoManager *hm = QAHistManagerDef::getHistoManager();
268  assert(hm);
269  TH1D *h_norm = dynamic_cast<TH1D *>(hm->getHisto(
270  get_histo_prefix() + "Normalization"));
271  assert(h_norm);
272  h_norm->Fill("Track", 1);
273 
274  {
275  TH1F *hsum = dynamic_cast<TH1F *>(hm->getHisto(
276  (get_histo_prefix()) + "TrackProj_3x3Tower_EP"));
277  assert(hsum);
278 
279  hsum->Fill(
280  (track->get_cal_energy_3x3(SvtxTrack::CEMC) + track->get_cal_energy_3x3(SvtxTrack::HCALIN) + track->get_cal_energy_3x3(SvtxTrack::HCALOUT)) / (primary->get_e() + 1e-9));
281  }
282  {
283  TH1F *hsum = dynamic_cast<TH1F *>(hm->getHisto(
284  (get_histo_prefix()) + "TrackProj_5x5Tower_EP"));
285  assert(hsum);
286 
287  hsum->Fill(
288  (track->get_cal_energy_5x5(SvtxTrack::CEMC) + track->get_cal_energy_5x5(SvtxTrack::HCALIN) + track->get_cal_energy_5x5(SvtxTrack::HCALOUT)) / (primary->get_e() + 1e-9));
289  }
290 
291  eval_trk_proj(_calo_name_cemc, track, topNode);
292  eval_trk_proj(_calo_name_hcalin, track, topNode);
293  eval_trk_proj(_calo_name_hcalout, track, topNode);
294 
295  return Fun4AllReturnCodes::EVENT_OK;
296 }
297 
298 bool QAG4SimulationEicCalorimeterSum::eval_trk_proj(const string &detector, SvtxTrack *track,
299  PHCompositeNode *topNode)
300 // Track projections
301 {
302  assert(track);
303  assert(topNode);
304 
305  Fun4AllHistoManager *hm = QAHistManagerDef::getHistoManager();
306  assert(hm);
307 
308  TH2F *h2_proj = dynamic_cast<TH2F *>(hm->getHisto(
309  (get_histo_prefix()) + detector + "_TrackProj"));
310  assert(h2_proj);
311 
312  // pull the tower geometry
313  string towergeonodename = "TOWERGEOM_" + detector;
314  RawTowerGeomContainer *towergeo = findNode::getClass<RawTowerGeomContainer>(
315  topNode, towergeonodename.c_str());
316  assert(towergeo);
317 
318  if (Verbosity() > 2)
319  {
320  towergeo->identify();
321  }
322  // pull the towers
323  string towernodename = "TOWER_CALIB_" + detector;
324  RawTowerContainer *towerList = findNode::getClass<RawTowerContainer>(topNode,
325  towernodename.c_str());
326  assert(towerList);
327 
328  if (Verbosity() > 3)
329  {
330  cout << __PRETTY_FUNCTION__ << " - info - handling track track p = ("
331  << track->get_px() << ", " << track->get_py() << ", "
332  << track->get_pz() << "), v = (" << track->get_x() << ", "
333  << track->get_z() << ", " << track->get_z() << ")"
334  << " size_states "
335  << track->size_states() << endl;
336  }
337 
338  // curved tracks inside mag field
339  // straight projections thereafter
340 
341  std::vector<double> point;
342  point.assign(3, NAN);
343 
344  // const double radius = towergeo->get_radius() + towergeo->get_thickness() * 0.5;
345 
346  // PHG4HoughTransform::projectToRadius(track, _magField, radius, point);
347 
348  if (std::isnan(point[0]) or std::isnan(point[1]) or std::isnan(point[2]))
349  {
350  // cout << __PRETTY_FUNCTION__ << "::" << Name()
351  // << " - Error - track extrapolation failure:";
352  // track->identify();
353  return false;
354  }
355 
356  assert(not std::isnan(point[0]));
357  assert(not std::isnan(point[1]));
358  assert(not std::isnan(point[2]));
359 
360  double x = point[0];
361  double y = point[1];
362  double z = point[2];
363 
364  double phi = atan2(y, x);
365  double eta = asinh(z / sqrt(x * x + y * y));
366 
367  // calculate 3x3 tower energy
368  int binphi = towergeo->get_phibin(phi);
369  int bineta = towergeo->get_etabin(eta);
370 
371  double etabin_width = towergeo->get_etabounds(bineta).second - towergeo->get_etabounds(bineta).first;
372  if (bineta > 1 and bineta < towergeo->get_etabins() - 1)
373  etabin_width = (towergeo->get_etacenter(bineta + 1) - towergeo->get_etacenter(bineta - 1)) / 2.;
374 
375  double phibin_width = towergeo->get_phibounds(binphi).second - towergeo->get_phibounds(binphi).first;
376 
377  assert(etabin_width > 0);
378  assert(phibin_width > 0);
379 
380  const double etabin_shift = (towergeo->get_etacenter(bineta) - eta) / etabin_width;
381  const double phibin_shift = (fmod(
382  towergeo->get_phicenter(binphi) - phi + 5 * M_PI, 2 * M_PI) -
383  M_PI) /
384  phibin_width;
385 
386  if (Verbosity() > 3)
387  cout << __PRETTY_FUNCTION__ << " - info - handling track proj. (" << x
388  << ", " << y << ", " << z << ")"
389  << ", eta " << eta << ", phi" << phi
390  << ", bineta " << bineta << " - ["
391  << towergeo->get_etabounds(bineta).first << ", "
392  << towergeo->get_etabounds(bineta).second << "], etabin_shift = "
393  << etabin_shift << ", binphi " << binphi << " - ["
394  << towergeo->get_phibounds(binphi).first << ", "
395  << towergeo->get_phibounds(binphi).second << "], phibin_shift = "
396  << phibin_shift << endl;
397 
398  const int bin_search_range = (Max_N_Tower - 1) / 2;
399  for (int iphi = binphi - bin_search_range; iphi <= binphi + bin_search_range;
400  ++iphi)
401 
402  for (int ieta = bineta - bin_search_range;
403  ieta <= bineta + bin_search_range; ++ieta)
404  {
405  // wrap around
406  int wrapphi = iphi;
407  if (wrapphi < 0)
408  {
409  wrapphi = towergeo->get_phibins() + wrapphi;
410  }
411  if (wrapphi >= towergeo->get_phibins())
412  {
413  wrapphi = wrapphi - towergeo->get_phibins();
414  }
415 
416  // edges
417  if (ieta < 0)
418  continue;
419  if (ieta >= towergeo->get_etabins())
420  continue;
421 
422  if (Verbosity() > 3)
423  cout << __PRETTY_FUNCTION__ << " - info - processing tower"
424  << ", bineta " << ieta << " - ["
425  << towergeo->get_etabounds(ieta).first << ", "
426  << towergeo->get_etabounds(ieta).second << "]"
427  << ", binphi "
428  << wrapphi << " - [" << towergeo->get_phibounds(wrapphi).first
429  << ", " << towergeo->get_phibounds(wrapphi).second << "]" << endl;
430 
431  RawTower *tower = towerList->getTower(ieta, wrapphi);
432  double energy = 0;
433  if (tower)
434  {
435  if (Verbosity() > 1)
436  cout << __PRETTY_FUNCTION__ << " - info - tower " << ieta << " "
437  << wrapphi << " energy = " << tower->get_energy() << endl;
438 
439  energy = tower->get_energy();
440  }
441 
442  h2_proj->Fill(ieta - bineta + etabin_shift,
443  iphi - binphi + phibin_shift, energy);
444 
445  } // for (int ieta = bineta-1; ieta < bineta+2; ++ieta) {
446 
447  return true;
448 } // // Track projections
449 
450 //int
451 //QAG4SimulationEicCalorimeterSum::Init_Tower(PHCompositeNode *topNode)
452 //{
453 //
454 // Fun4AllHistoManager *hm = QAHistManagerDef::getHistoManager();
455 // assert(hm);
456 //
458 //
459 // return Fun4AllReturnCodes::EVENT_OK;
460 //}
461 //
462 //int
463 //QAG4SimulationEicCalorimeterSum::process_event_Tower(PHCompositeNode *topNode)
464 //{
465 // if (Verbosity() > 2)
466 // cout << "QAG4SimulationEicCalorimeterSum::process_event_Tower() entered"
467 // << endl;
468 //
469 // return Fun4AllReturnCodes::EVENT_OK;
470 //}
471 
472 int QAG4SimulationEicCalorimeterSum::Init_Cluster(PHCompositeNode *topNode)
473 {
474  Fun4AllHistoManager *hm = QAHistManagerDef::getHistoManager();
475  assert(hm);
476 
477  hm->registerHisto(
478  new TH2F(
479  TString(get_histo_prefix()) + "Cluster_" + _calo_name_cemc.c_str() + "_" + _calo_name_hcalin.c_str(), //
480  TString(_calo_name_hcalin.c_str()) + " VS " + TString(_calo_name_cemc.c_str()) + ": best cluster energy;" + TString(_calo_name_cemc.c_str()) + " cluster energy (GeV);" + TString(_calo_name_hcalin.c_str()) + " cluster energy (GeV)",
481  70, 0, 70, 70, 0, 70));
482 
483  hm->registerHisto(
484  new TH2F(
485  TString(get_histo_prefix()) + "Cluster_" + _calo_name_cemc.c_str() + "_" + _calo_name_hcalin.c_str() + "_" + _calo_name_hcalout.c_str(), //
486  TString(_calo_name_cemc.c_str()) + " + " + TString(_calo_name_hcalin.c_str()) + " VS " + TString(_calo_name_hcalout.c_str()) + ": best cluster energy;" + TString(_calo_name_cemc.c_str()) + " + " + TString(_calo_name_hcalin.c_str()) + " cluster energy (GeV);" + TString(_calo_name_hcalout.c_str()) + " cluster energy (GeV)",
487  70, 0, 70, 70, 0, 70));
488 
489  hm->registerHisto(
490  new TH1F(TString(get_histo_prefix()) + "Cluster_EP", //
491  "Total Cluster E_{Reco}/E_{Truth};Reco cluster energy sum / total truth energy",
492  150, 0, 1.5));
493 
494  hm->registerHisto(
495  new TH1F(
496  TString(get_histo_prefix()) + "Cluster_Ratio_" + _calo_name_cemc.c_str() + "_" + _calo_name_hcalin.c_str(), //
497  "Energy ratio " + TString(_calo_name_cemc.c_str()) + " VS " + TString(_calo_name_hcalin.c_str()) + ";Best cluster " + TString(_calo_name_cemc.c_str()) + " / (" + TString(_calo_name_cemc.c_str()) + " + " + TString(_calo_name_hcalin.c_str()) + ")", 110, 0, 1.1));
498 
499  hm->registerHisto(
500  new TH1F(
501  TString(get_histo_prefix()) + "Cluster_Ratio_" + _calo_name_cemc.c_str() + "_" + _calo_name_hcalin.c_str() + "_" + TString(_calo_name_hcalout.c_str()), //
502  "Energy ratio " + TString(_calo_name_cemc.c_str()) + " + " + TString(_calo_name_hcalin.c_str()) + " VS " + TString(_calo_name_hcalout.c_str()) + ";Best cluster (" + TString(_calo_name_cemc.c_str()) + " + " + TString(_calo_name_hcalin.c_str()) + ") / (" + TString(_calo_name_cemc.c_str()) + " + " + TString(_calo_name_hcalin.c_str()) + " + " + TString(_calo_name_hcalout.c_str()) + ")", 110, 0, 1.1));
503 
504  return Fun4AllReturnCodes::EVENT_OK;
505 }
506 
507 int QAG4SimulationEicCalorimeterSum::process_event_Cluster(PHCompositeNode *topNode)
508 {
509  if (Verbosity() > 2)
510  cout << "QAG4SimulationEicCalorimeterSum::process_event_Cluster() entered"
511  << endl;
512 
513  Fun4AllHistoManager *hm = QAHistManagerDef::getHistoManager();
514  assert(hm);
515 
516  PHG4Particle *primary = get_truth_particle();
517  if (!primary)
518  return Fun4AllReturnCodes::DISCARDEVENT;
519 
520  RawCluster *cluster_cemc =
521  _caloevalstack_cemc->get_rawcluster_eval()->best_cluster_from(primary);
522  RawCluster *cluster_hcalin =
523  _caloevalstack_hcalin->get_rawcluster_eval()->best_cluster_from(primary);
524  RawCluster *cluster_hcalout =
525  _caloevalstack_hcalout->get_rawcluster_eval()->best_cluster_from(primary);
526 
527  const double cluster_cemc_e = cluster_cemc ? cluster_cemc->get_energy() : 0;
528  const double cluster_hcalin_e =
529  cluster_hcalin ? cluster_hcalin->get_energy() : 0;
530  const double cluster_hcalout_e =
531  cluster_hcalout ? cluster_hcalout->get_energy() : 0;
532 
533  if (cluster_cemc_e + cluster_hcalin_e > 0)
534  {
535  TH2F *h2 = dynamic_cast<TH2F *>(hm->getHisto(
536  (get_histo_prefix()) + "Cluster_" + _calo_name_cemc + "_" + _calo_name_hcalin));
537  assert(h2);
538 
539  h2->Fill(cluster_cemc_e, cluster_hcalin_e);
540 
541  TH1F *hr = dynamic_cast<TH1F *>(hm->getHisto(
542  (get_histo_prefix()) + "Cluster_Ratio_" + _calo_name_cemc + "_" + _calo_name_hcalin));
543  assert(hr);
544 
545  hr->Fill(cluster_cemc_e / (cluster_cemc_e + cluster_hcalin_e));
546  }
547 
548  if (cluster_cemc_e + cluster_hcalin_e + cluster_hcalout_e > 0)
549  {
550  if (Verbosity())
551  {
552  cout << "QAG4SimulationEicCalorimeterSum::process_event_Cluster - "
553  << " cluster_cemc_e = " << cluster_cemc_e
554  << " cluster_hcalin_e = " << cluster_hcalin_e
555  << " cluster_hcalout_e = " << cluster_hcalout_e << " hr = "
556  << (cluster_cemc_e + cluster_hcalin_e) / (cluster_cemc_e + cluster_hcalin_e + cluster_hcalout_e)
557  << endl;
558  }
559 
560  TH2F *h2 = dynamic_cast<TH2F *>(hm->getHisto(
561  (get_histo_prefix()) + "Cluster_" + _calo_name_cemc + "_" + _calo_name_hcalin + "_" + _calo_name_hcalout));
562  assert(h2);
563 
564  h2->Fill((cluster_cemc_e + cluster_hcalin_e), cluster_hcalout_e);
565 
566  TH1F *hr = dynamic_cast<TH1F *>(hm->getHisto(
567  (get_histo_prefix()) + "Cluster_Ratio_" + _calo_name_cemc + "_" + _calo_name_hcalin + "_" + _calo_name_hcalout));
568  assert(hr);
569 
570  hr->Fill(
571  (cluster_cemc_e + cluster_hcalin_e) / (cluster_cemc_e + cluster_hcalin_e + cluster_hcalout_e));
572 
573  TH1F *hsum = dynamic_cast<TH1F *>(hm->getHisto(
574  (get_histo_prefix()) + "Cluster_EP"));
575  assert(hsum);
576 
577  hsum->Fill(
578  (cluster_cemc_e + cluster_hcalin_e + cluster_hcalout_e) / (primary->get_e() + 1e-9));
579  }
580 
581  return Fun4AllReturnCodes::EVENT_OK;
582 }