3 #include <calobase/RawTower.h>
4 #include <calobase/RawTowerContainer.h>
5 #include <calobase/RawTowerGeom.h>
6 #include <calobase/RawTowerGeomContainer.h>
45 std::cout <<
Name() <<
": No Calo Type set" << std::endl;
55 hcalin_energy_eta =
new TH2F(
"hcalin_energy_eta",
"hcalin energy eta", 1000, 0, 10, 240, -1.1, 1.1);
56 hcalin_e_eta_phi =
new TH3F(
"hcalin_e_eta_phi",
"hcalin e eta phi", 50, 0, 10, 24, -1.1, 1.1, 64, -3.14159, 3.14159);
57 for (
int i = 0; i < 24; i++)
59 for (
int j = 0; j < 64; j++)
63 hcal_in_eta_phi[i][j] =
new TH1F(hist_name.c_str(),
"Hcal_in_energy", 1000, 0, 10);
67 for (
int i = 0; i < 24; i++)
71 hcalin_eta[i] =
new TH1F(hist_name.c_str(),
"hcalin eta's", 1000, 0, 10);
76 hcalout_energy_eta =
new TH2F(
"hcalout_energy_eta",
"hcalout energy eta", 10, 0, 10, 24000, -1.1, 1.1);
77 hcalout_e_eta_phi =
new TH3F(
"hcalout_e_eta_phi",
"hcalout e eta phi", 50, 0, 10, 24, -1.1, 1.1, 64, -3.14159, 3.14159);
78 for (
int i = 0; i < 24; i++)
80 for (
int j = 0; j < 64; j++)
84 hcal_out_eta_phi[i][j] =
new TH1F(hist_name.c_str(),
"Hcal_out energy", 1000, 0, 10);
88 for (
int i = 0; i < 24; i++)
92 hcalout_eta[i] =
new TH1F(hist_name.c_str(),
"hcalout eta's", 1000, 0, 10);
97 for (
int i = 0; i < 96; i++)
99 for (
int j = 0; j < 258; j++)
103 cemc_hist_eta_phi[i][j] =
new TH1F(hist_name.c_str(),
"Hist_ieta_phi_leaf(e)", 1000, 0, 10);
108 for (
int i = 0; i < 96; i++)
110 gStyle->SetOptFit(1);
113 eta_hist[i] =
new TH1F(b.c_str(),
"eta and all phi's", 1000, 0, 10);
117 energy_eta_hist =
new TH2F(
"energy_eta_hist",
"energy eta and all phi", 10, 0, 10, 9600, -1, 1);
120 e_eta_phi =
new TH3F(
"e_eta_phi",
"e v eta v phi", 50, 0, 10, 100, -1, 1, 256, -3.14159, 3.14159);
129 if (
_ievent % 100 == 0) std::cout <<
"LiteCaloEval::process_event(PHCompositeNode *topNode) Processing Event " <<
_ievent << std::endl;
131 std::string towernode =
"TOWER_CALIB_" +
_caloname;
132 RawTowerContainer* towers = findNode::getClass<RawTowerContainer>(topNode, towernode);
135 std::cout <<
PHWHERE <<
" ERROR: Can't find " << towernode << std::endl;
139 std::string towergeomnode =
"TOWERGEOM_" +
_caloname;
140 RawTowerGeomContainer* towergeom = findNode::getClass<RawTowerGeomContainer>(topNode, towergeomnode);
143 std::cout <<
PHWHERE <<
" ERROR: Can't find " << towergeomnode << std::endl;
149 for (rtiter = begin_end.first; rtiter != begin_end.second; ++rtiter)
158 std::cout <<
PHWHERE <<
" ERROR: Can't find tower geometry for this tower hit: ";
163 const int towerid = tower->
get_id();
173 std::cout <<
"a towerid was less than 0 " << std::endl;
227 double par_value[24];
228 double eta_value[24];
233 double par_value2[24];
237 for (
int i = 0; i < 24; i++)
246 std::unique_ptr<TF1>
f1(
new TF1(
"f1",
"expo", 0.02, 0.1));
247 std::unique_ptr<TF1>
f2(
new TF1(
"f2",
"expo", 0.1, 1));
253 a = f1->GetParameter(1);
254 b = f1->GetParError(1);
256 c = f2->GetParameter(1);
257 d = f2->GetParError(1);
260 par_value[i] =
abs(a);
263 rel_err[i] = par_err[i] / par_value[i];
265 par_value2[i] =
abs(c);
267 rel_err2[i] = par_err2[i] / par_value2[i];
271 TGraph g1(24, eta_value, par_value);
272 g1.SetTitle(
"HCal In (0.02-0.1 GeV); eta; p1");
273 g1.SetMarkerStyle(20);
275 g1.SetName(
"Fit1_hcalin");
278 TGraph g2(24, eta_value, rel_err);
279 g2.SetTitle(
"HCal In Error (0.02-0.1 GeV); eta; p1 rel error");
280 g2.SetMarkerStyle(20);
282 g2.SetName(
"Fit1_err_hcalin");
285 TGraph g3(24, eta_value, par_value2);
286 g3.SetTitle(
"HCal In (0.1-1 GeV); eta; p1");
287 g3.SetMarkerStyle(20);
289 g3.SetName(
"Fit2_hcalin");
292 TGraph
g4(24, eta_value, rel_err2);
293 g4.SetTitle(
"HCal In Error (0.1-1 GeV); eta; p1 rel error");
294 g4.SetMarkerStyle(20);
296 g4.SetName(
"Fit2_err_hcalin");
301 double par_value[24];
302 double eta_value[24];
307 double par_value2[24];
311 double par_value3[24];
315 for (
int i = 0; i < 24; i++)
326 std::unique_ptr<TF1>
f1(
new TF1(
"f1",
"expo", 0.05, 0.2));
327 std::unique_ptr<TF1>
f2(
new TF1(
"f2",
"expo", 0.2, 1));
328 std::unique_ptr<TF1>
f3(
new TF1(
"f3",
"expo", 1, 2));
336 a = f1->GetParameter(1);
337 b = f1->GetParError(1);
339 c = f2->GetParameter(1);
340 d = f2->GetParError(1);
342 e = f3->GetParameter(1);
343 f = f3->GetParError(1);
346 par_value[i] =
abs(a);
349 rel_err[i] = par_err[i] / par_value[i];
351 par_value2[i] =
abs(c);
353 rel_err2[i] = par_err2[i] / par_value2[i];
355 par_value3[i] =
abs(e);
357 rel_err3[i] = par_err3[i] / par_value3[i];
361 TGraph g1(24, eta_value, par_value);
362 g1.SetTitle(
"HCal Out (0.05-0.2 GeV); eta; p1");
363 g1.SetMarkerStyle(20);
365 g1.SetName(
"Fit1_hcalout");
368 TGraph g2(24, eta_value, rel_err);
369 g2.SetTitle(
"HCal Out Error (0.05-0.2 GeV); eta; p1 rel err");
370 g2.SetMarkerStyle(20);
372 g2.SetName(
"Fit1_err_hcalout");
375 TGraph g3(24, eta_value, par_value2);
376 g3.SetTitle(
"HCal Out (0.2-1 GeV); eta; p1");
377 g3.SetMarkerStyle(20);
379 g3.SetName(
"Fit2_hcalout");
382 TGraph
g4(24, eta_value, rel_err2);
383 g4.SetTitle(
"HCal Out Error (0.2-1 GeV); eta; p1 rel err");
384 g4.SetMarkerStyle(20);
386 g4.SetName(
"Fit2_err_hcalout");
389 TGraph g5(24, eta_value, par_value3);
390 g5.SetTitle(
"HCal Out (1-2 GeV); eta; p1");
391 g5.SetMarkerStyle(20);
393 g5.SetName(
"Fit3_hcalout");
396 TGraph g6(24, eta_value, rel_err3);
397 g6.SetTitle(
"HCal Out Error (1-2 GeV); eta; p1 rel err");
398 g6.SetMarkerStyle(20);
400 g6.SetName(
"Fit3_err_hcalout");
406 double par_value[96];
407 double eta_value[96];
412 double par_value2[96];
416 double par_value3[96];
421 for (
int i = 0; i < 96; i++)
432 std::unique_ptr<TF1>
f1(
new TF1(
"f1",
"expo", 0.04, 0.1));
433 std::unique_ptr<TF1>
f2(
new TF1(
"f2",
"expo", 0.1, 0.4));
434 std::unique_ptr<TF1>
f3(
new TF1(
"f3",
"expo", 0.4, 2));
442 a = f1->GetParameter(1);
443 b = f1->GetParError(1);
445 c = f2->GetParameter(1);
446 d = f2->GetParError(1);
448 e = f3->GetParameter(1);
449 f = f3->GetParError(1);
452 par_value[i] =
abs(a);
455 rel_err[i] = par_err[i] / par_value[i];
457 par_value2[i] =
abs(c);
459 rel_err2[i] = par_err2[i] / par_value2[i];
461 par_value3[i] =
abs(e);
463 rel_err3[i] = par_err3[i] / par_value3[i];
467 TGraph g1(96, eta_value, par_value);
468 g1.SetTitle(
"EMCal (0.04-0.1 GeV); eta; p1");
469 g1.SetMarkerStyle(20);
471 g1.SetName(
"Fit1_emc");
474 TGraph g2(96, eta_value, rel_err);
475 g2.SetTitle(
"EMCal Error (0.04-0.1 GeV); eta; p1 rel error");
476 g2.SetMarkerStyle(20);
478 g2.SetName(
"Fit1_err_emc");
481 TGraph g3(96, eta_value, par_value2);
482 g3.SetTitle(
"EMCal (0.1-0.4 GeV); eta; p1");
483 g3.SetMarkerStyle(20);
485 g3.SetName(
"Fit2_emc");
488 TGraph
g4(96, eta_value, rel_err2);
489 g4.SetTitle(
"EMCal Error (0.1-0.4 GeV); eta; p1 rel error");
490 g4.SetMarkerStyle(20);
492 g4.SetName(
"Fit2_err_emc");
495 TGraph g5(96, eta_value, par_value3);
496 g5.SetTitle(
"EMCal (0.4-2 GeV); eta; p1");
497 g5.SetMarkerStyle(20);
499 g5.SetName(
"Fit3_emc");
502 TGraph g6(96, eta_value, rel_err3);
503 g6.SetTitle(
"EMCal Error (0.4-2 GeV); eta; p1 rel err");
504 g6.SetMarkerStyle(20);
506 g6.SetName(
"Fit3_err_emc");