Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DrawPrototype4EMCalTower.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file DrawPrototype4EMCalTower.C
1 // $Id: $
2 
11 #include <TFile.h>
12 #include <TGraphErrors.h>
13 #include <TLatex.h>
14 #include <TLine.h>
15 #include <TString.h>
16 #include <TTree.h>
17 #include <cassert>
18 #include <cmath>
19 #include "SaveCanvas.C"
20 #include "SetOKStyle.C"
21 using namespace std;
22 
23 //#include "Prototype4_DSTReader.h"
24 
25 TFile *_file0 = NULL;
26 TTree *T = NULL;
27 TString cuts = "";
28 double beam_momentum_selection = -16;
29 
30 void DrawPrototype4EMCalTower( //
31  // const TString infile = "/gpfs/mnt/gpfs04/sphenix/user/jinhuang/Prototype_2018/test_production_5/beam_00000683-0000_DSTReader.root", //
32  // const TString infile = "/gpfs/mnt/gpfs04/sphenix/user/jinhuang/Prototype_2018/Scan2Block18/beam_0000068ALL-0000_DSTReader.root", //
33  // const TString infile = "/gpfs/mnt/gpfs04/sphenix/user/jinhuang/Prototype_2018/Scan1Block36/2GeV/beam_00000ALL-0000_DSTReader.root", //
34  const TString infile = "data/beam_00001214.root_DSTReader.root", //
35  bool plot_all = false, const double momentum = -16)
36 {
37  beam_momentum_selection = momentum;
38 
39  SetOKStyle();
40  gStyle->SetOptStat(0);
41  gStyle->SetOptFit(1111);
42  TVirtualFitter::SetDefaultFitter("Minuit2");
43  gSystem->Load("libg4eval.so");
44  gSystem->Load("libqa_modules.so");
45  gSystem->Load("libPrototype4.so");
46 
47  // gROOT->LoadMacro("Prototype4_DSTReader.C+");
48 
49  if (!_file0)
50  {
51  TString chian_str = infile;
52  chian_str.ReplaceAll("ALL", "*");
53 
54  TChain *t = new TChain("T");
55  const int n = t->Add(chian_str);
56 
57  cout << "Loaded " << n << " root files with " << chian_str << endl;
58  assert(n > 0);
59 
60  T = t;
61 
62  _file0 = new TFile;
63  _file0->SetName(infile);
64  }
65 
66  assert(_file0);
67 
68  T->SetAlias("ActiveTower_LG",
69  "TOWER_LG_CEMC[].get_binphi()<8 && TOWER_LG_CEMC[].get_bineta()<8");
70  T->SetAlias("EnergySum_LG",
71  "1*Sum$(TOWER_LG_CEMC[].get_energy() * ActiveTower_LG)");
72 
73  T->SetAlias("ActiveTower_HG",
74  "TOWER_HG_CEMC[].get_binphi()<8 && TOWER_HG_CEMC[].get_bineta()<8");
75  T->SetAlias("EnergySum_HG",
76  "1*Sum$(TOWER_HG_CEMC[].get_energy() * ActiveTower_HG)");
77 
78  // T->SetAlias("C2_Inner_e", "1*TOWER_RAW_C2[0].energy");
79  T->SetAlias("C2_Inner_e", "1*(TOWER_CALIB_C2[0].energy)");
80  T->SetAlias("C2_Outer_e", "1*(TOWER_CALIB_C2[1].energy)");
81  T->SetAlias("C2_Sum_e", "C2_Inner_e + C2_Outer_e");
82  T->SetAlias("C1", "0 + TOWER_CALIB_C1[0].energy");
83 
84  // "TOWER_CALIB_CEMC.energy * ( Sum$( TOWER_CALIB_CEMC.get_column()==2 && TOWER_CALIB_CEMC.get_row()==1
85 
86  T->SetAlias("Average_column",
87  "Sum$(TOWER_CALIB_CEMC.get_column() * TOWER_CALIB_CEMC.get_energy())/Sum$(TOWER_CALIB_CEMC.get_energy())");
88  T->SetAlias("Average_row",
89  "Sum$(TOWER_CALIB_CEMC.get_row() * TOWER_CALIB_CEMC.get_energy())/Sum$(TOWER_CALIB_CEMC.get_energy())");
90 
91  T->SetAlias("Average_HODO_VERTICAL",
92  "Sum$(TOWER_CALIB_HODO_VERTICAL.towerid * (abs(TOWER_CALIB_HODO_VERTICAL.energy)>0.5) * abs(TOWER_CALIB_HODO_VERTICAL.energy))/Sum$((abs(TOWER_CALIB_HODO_VERTICAL.energy)>0.5) * abs(TOWER_CALIB_HODO_VERTICAL.energy))");
93  T->SetAlias("Valid_HODO_VERTICAL", "Sum$((TOWER_CALIB_HODO_VERTICAL.energy)>0.5) > 0");
94 
95  T->SetAlias("Average_HODO_HORIZONTAL",
96  "Sum$(TOWER_CALIB_HODO_HORIZONTAL.towerid * (abs(TOWER_CALIB_HODO_HORIZONTAL.energy)>0.5) * abs(TOWER_CALIB_HODO_HORIZONTAL.energy))/Sum$((abs(TOWER_CALIB_HODO_HORIZONTAL.energy)>0.5) * abs(TOWER_CALIB_HODO_HORIZONTAL.energy))");
97  T->SetAlias("Valid_HODO_HORIZONTAL", "Sum$((TOWER_CALIB_HODO_HORIZONTAL.energy)>0.5) > 0");
98 
99  T->SetAlias("No_Triger_VETO",
100  "Sum$((TOWER_CALIB_TRIGGER_VETO.energy)>0.2)==0");
101 
102  T->SetAlias("Energy_Sum_col1_row2_3x3",
103  "Sum$( (abs(TOWER_CALIB_CEMC.get_column()-1)<=1 && abs(TOWER_CALIB_CEMC.get_row()-2)<=1 ) * TOWER_CALIB_CEMC.get_energy())");
104  T->SetAlias("Energy_Sum_col1_row2_5x5",
105  "Sum$( (abs(TOWER_CALIB_CEMC.get_column()-1)<=2 && abs(TOWER_CALIB_CEMC.get_row()-2)<=2 ) * TOWER_CALIB_CEMC.get_energy())");
106  T->SetAlias("Energy_Sum_col2_row2_5x5",
107  "Sum$( (abs(TOWER_CALIB_CEMC.get_column()-2)<=2 && abs(TOWER_CALIB_CEMC.get_row()-2)<=2 ) * TOWER_CALIB_CEMC.get_energy())");
108  T->SetAlias("Energy_Sum_CEMC", "1*Sum$(TOWER_CALIB_CEMC.get_energy())");
109  T->SetAlias("Energy_Sum_RAW_CEMC", "1*Sum$(TOWER_RAW_CEMC.get_energy())");
110 
111  // 12 GeV calibration
112  // EDM=9.83335e-18 STRATEGY= 1 ERROR MATRIX ACCURATE
113  //EXT PARAMETER STEP FIRST
114  //NO. NAME VALUE ERROR SIZE DERIVATIVE
115  //1 p0 1.19768e+01 7.30605e-02 3.76799e-05 4.24290e-09
116  //2 p1 8.71776e+00 6.82987e-02 3.52240e-05 6.80808e-08
117 
118  // T->SetAlias("Energy_Sum_Hadron_CEMC",
119  // "1*Sum$(TOWER_CALIB_CEMC.get_energy())"); // full bias
120  T->SetAlias("Energy_Sum_Hadron_CEMC",
121  "1.14*12./8.71776e+00*Sum$(TOWER_CALIB_CEMC.get_energy())"); // full bias
122  // T->SetAlias("Energy_Sum_Hadron_CEMC",
123  // "1.14*12./8.71776e+00*(16./6.93250e+00)*(28/33.3405)*Sum$(TOWER_CALIB_CEMC.get_energy())"); // half-gain bias
124  T->SetAlias("CEMC_MIP", "Energy_Sum_Hadron_CEMC<0.7");
125 
126  // 12 GeV calibration
127  // FCN=9.63681 FROM HESSE STATUS=OK 14 CALLS 56 TOTAL
128  // EDM=1.49963e-17 STRATEGY= 1 ERROR MATRIX ACCURATE
129  // EXT PARAMETER STEP FIRST
130  // NO. NAME VALUE ERROR SIZE DERIVATIVE
131  // 1 p0 9.50430e+00 8.42691e-02 3.83666e-06 1.41105e-08
132  // 2 p1 6.99727e+00 1.06583e-01 4.85258e-06 5.85631e-08
133 
134  // T->SetAlias("Energy_Sum_Hadron_HCALIN",
135  // "1*Sum$(TOWER_CALIB_LG_HCALIN.get_energy())");
136  // T->SetAlias("HCALIN_MIP", "Energy_Sum_Hadron_HCALIN<0.5");
137  // T->SetAlias("Energy_Sum_Hadron_HCALOUT",
138  // "1*Sum$(TOWER_CALIB_LG_HCALOUT.get_energy())");
139  T->SetAlias("Energy_Sum_Hadron_HCALIN",
140  "12./6.99727e+00*Sum$(TOWER_CALIB_LG_HCALIN.get_energy())");
141  T->SetAlias("HCALIN_MIP", "Energy_Sum_Hadron_HCALIN<0.5");
142  T->SetAlias("Energy_Sum_Hadron_HCALOUT",
143  "12./9.50430e+00*Sum$(TOWER_CALIB_LG_HCALOUT.get_energy())");
144 
145  T->SetAlias("MIP_Count_Col2",
146  "Sum$( abs( TOWER_RAW_CEMC.get_energy() )>20 && abs( TOWER_RAW_CEMC.get_energy() )<200 && TOWER_CALIB_CEMC.get_column() == 2 )");
147  T->SetAlias("Pedestal_Count_AllCEMC",
148  "Sum$( abs( TOWER_RAW_CEMC.get_energy() )<20)");
149 
150  //
151  TCut event_sel = "1*1";
152 
153  if (plot_all)
154  {
155  // event_sel = "1*1";
156  // cuts = "_all_event";
157  // event_sel = "Valid_HODO_HORIZONTAL && Valid_HODO_VERTICAL";
158  // cuts = "_Valid_HODO";
159  event_sel =
160  "Valid_HODO_HORIZONTAL && Valid_HODO_VERTICAL && No_Triger_VETO";
161  cuts = "_Valid_HODO_Trigger_VETO";
162  }
163  else
164  {
165  if (0)
166  { // energy selection
167  event_sel = Form("(beam_MTNRG_GeV == %f)", beam_momentum_selection);
168  cuts = Form("_%.0fGeV", beam_momentum_selection);
169 
170  cout << "Build event selection of " << (const char *) event_sel
171  << endl;
172 
173  T->Draw(">>EventListRunCut", event_sel);
174  TEventList *elist = gDirectory->GetObjectChecked("EventListRunCut",
175  "TEventList");
176  cout << elist->GetN() << " / " << T->GetEntriesFast()
177  << " events selected" << endl;
178  T->SetEventList(elist);
179  }
180 
181  // event_sel = "1*1";
182  // cuts = "_all_event";
183 
184  // event_sel = "Valid_HODO_HORIZONTAL && Valid_HODO_VERTICAL";
185  // cuts = "_Valid_HODO";
186 
187  event_sel =
188  "Valid_HODO_HORIZONTAL && Valid_HODO_VERTICAL && No_Triger_VETO";
189  cuts = "_Valid_HODO_Trigger_VETO";
190 
191 // event_sel =
192 // "Valid_HODO_HORIZONTAL && Valid_HODO_VERTICAL && No_Triger_VETO && abs(Average_HODO_HORIZONTAL - 3)<.2 && abs(Average_HODO_VERTICAL - 4)<.2 ";
193 // cuts = "_Valid_HODO_Trigger_VETO_h3_v4";
194  // event_sel =
195  // "Valid_HODO_HORIZONTAL && Valid_HODO_VERTICAL && No_Triger_VETO && abs(Average_HODO_HORIZONTAL - 3)<.2 && abs(Average_HODO_VERTICAL - 4)<.2 ";
196  // cuts = "_Valid_HODO_Trigger_VETO_h3_v4";
197  // event_sel =
198  // "Valid_HODO_HORIZONTAL && Valid_HODO_VERTICAL && No_Triger_VETO && abs(Average_HODO_HORIZONTAL - 3)<.2 && abs(Average_HODO_VERTICAL - 4)<.2 && (C2_Outer_e>2000 || C2_Inner_e>1000)";
199  // cuts = "_Valid_HODO_Trigger_VETO_good_e_h3_v4";
200  // event_sel =
201  // "Entry$==1643 && Valid_HODO_HORIZONTAL && Valid_HODO_VERTICAL && No_Triger_VETO && abs(Average_HODO_HORIZONTAL - 5)<.2 && abs(Average_HODO_VERTICAL - 2)<.2 && (C2_Outer_e>2000 || C2_Inner_e>1000)";
202  // cuts = "_Valid_HODO_Trigger_VETO_good_e_h5_v2_ev1643";
203  }
204 
205  cout << "Build event selection of " << (const char *) event_sel << endl;
206 
207  T->Draw(">>EventList", event_sel);
208  TEventList *elist = gDirectory->GetObjectChecked("EventList", "TEventList");
209  cout << elist->GetN() << " / " << T->GetEntriesFast() << " events selected"
210  << endl;
211 
212  T->SetEventList(elist);
213 
214  // int rnd = rand();
215  // gDirectory->mkdir(Form("dir_%d", rnd));
216  // gDirectory->cd(Form("dir_%d", rnd));
217  // if (plot_all)
218  // EMCDistribution_SUM("Energy_Sum_col1_row2_5x5");
219 
220  int rnd = rand();
221  gDirectory->mkdir(Form("dir_%d", rnd));
222  gDirectory->cd(Form("dir_%d", rnd));
223 // if (plot_all)
224  EMC_HodoScope_Calibration("C2_Sum_e");
225 
226  int rnd = rand();
227  gDirectory->mkdir(Form("dir_%d", rnd));
228  gDirectory->cd(Form("dir_%d", rnd));
229  if (plot_all)
230  EMCDistribution_ShowShape("C2_Sum_e");
231 
232  int rnd = rand();
233  gDirectory->mkdir(Form("dir_%d", rnd));
234  gDirectory->cd(Form("dir_%d", rnd));
235  if (plot_all)
236  EMCDistribution_SUM("Energy_Sum_CEMC", "C2_Sum_e");
237 
238  int rnd = rand();
239  gDirectory->mkdir(Form("dir_%d", rnd));
240  gDirectory->cd(Form("dir_%d", rnd));
241  // if (plot_all)
242  EMCDistribution_SUM_RawADC("Energy_Sum_RAW_CEMC", "C2_Sum_e");
243 
244  int rnd = rand();
245  gDirectory->mkdir(Form("dir_%d", rnd));
246  gDirectory->cd(Form("dir_%d", rnd));
247  if (plot_all)
248  EMCDistribution_HCalCalibration();
249 
250  int rnd = rand();
251  gDirectory->mkdir(Form("dir_%d", rnd));
252  gDirectory->cd(Form("dir_%d", rnd));
253  if (plot_all)
254  EMCDistribution_Fast();
255 
256  int rnd = rand();
257  gDirectory->mkdir(Form("dir_%d", rnd));
258  gDirectory->cd(Form("dir_%d", rnd));
259  if (plot_all)
260  EMCDistribution_Fast("RAW");
261 
262  int rnd = rand();
263  gDirectory->mkdir(Form("dir_%d", rnd));
264  gDirectory->cd(Form("dir_%d", rnd));
265  if (plot_all)
266  EMCDistribution_Fast("RAW", true);
267  //
268  // int rnd = rand();
269  // gDirectory->mkdir(Form("dir_%d", rnd));
270  // gDirectory->cd(Form("dir_%d", rnd));
271  // if (plot_all)
272  // EMCDistribution_PeakSample_Fast();
273 
274  int rnd = rand();
275  gDirectory->mkdir(Form("dir_%d", rnd));
276  gDirectory->cd(Form("dir_%d", rnd));
277  if (plot_all)
278  EMCDistribution_Fast("CALIB", true);
279 
280  int rnd = rand();
281  gDirectory->mkdir(Form("dir_%d", rnd));
282  gDirectory->cd(Form("dir_%d", rnd));
283 // if (plot_all)
284  EMCDistribution_ADC();
285 
286  // if (!plot_all)
287  // T->Process("Prototype4_DSTReader.C+",
288  // TString(_file0->GetName())
289  // + TString("_DrawPrototype4EMCalTower_Prototype4_DSTReader") + cuts
290  // + TString(".dat"));
291 }
292 
293 void EMC_HodoScope_Calibration(TString Cherenkov_Choice)
294 {
295  TString cut_C(Cherenkov_Choice + " > 1000");
296  TString cuts = Cherenkov_Choice;
297 
298  gStyle->SetOptStat(0);
299  gStyle->SetOptFit(1111);
300 
301  TH1F *h1_HCalOut = new TH1F("h1_HCalOut",
302  "(CEMC MIP && HCal_{IN} MIP);EMCal + HCal_{IN} + HCal_{OUT} (GeV)", 100,
303  0, abs(beam_momentum_selection) * 2.);
304 
305  TH2 *h2_hodo_v_h =
306  new TH2F("h2_hodo_v_h",
307  "Hodoscope hit distribution;Average Horizontal Hodoscope Idx;Average Vertical Hodoscope Idx",
308  8, -.5, 7.5, 8, -.5, 7.5);
309 
310  TH2 *h2_hodoE_v_h =
311  new TH2F("h2_hodoE_v_h",
312  "Event-Averaged Tower-Summed EMCal Energy [A.U.];Average Horizontal Hodoscope Idx;Average Vertical Hodoscope Idx",
313  8, -.5, 7.5, 8, -.5, 7.5);
314  TH2 *h2_hodoE_v_h_norm =
315  new TH2F("h2_hodoE_v_h_norm",
316  "Event-Averaged Tower-Summed EMCal Energy [A.U.];Average Horizontal Hodoscope Idx;Average Vertical Hodoscope Idx",
317  8, -.5, 7.5, 8, -.5, 7.5);
318 
319  TH2 *h2_hodo_c_h =
320  new TH2F("h2_hodo_c_h",
321  "Column VS Horizontal Hodoscope;Average Horizontal Hodoscope Idx;Average EMCal Column",
322  24, -.5, 7.5, 24, -.5, 7.5);
323  TH2 *h2_hodo_r_v =
324  new TH2F("h2_hodo_r_v",
325  "Row VS Vertical Hodoscope;Average Vertical Hodoscope Idx;Average EMCal Row",
326  24, -.5, 7.5, 24, -.5, 7.5);
327 
328  TH2 *h2_hodoE_h =
329  new TH2F("h2_hodoE_h",
330  "Total EMCal Energy VS Horizontal Hodoscope;Average Horizontal Hodoscope Idx;Total EMCal Energy [A.U.]",
331  8, -.5, 7.5, 80, -.5, 70);
332  TH2 *h2_hodoE_v =
333  new TH2F("h2_hodoE_v",
334  "Total EMCal Energy VS Vertical Hodoscope;Average Vertical Hodoscope Idx;Total EMCal Energy [A.U.]",
335  8, -.5, 7.5, 80, -.5, 70);
336 
337  TText *t;
338  TCanvas *c1 = new TCanvas("EMC_HodoScope_Calibration" + cuts,
339  "EMC_HodoScope_Calibration" + cuts, 1800, 1100);
340  c1->Divide(3, 2);
341  int idx = 1;
342  TPad *p;
343 
344  p = (TPad *) c1->cd(idx++);
345  c1->Update();
346 
347  T->Draw("Average_HODO_VERTICAL:Average_HODO_HORIZONTAL>>h2_hodo_v_h", cut_C, "colz");
348 
349  p = (TPad *) c1->cd(idx++);
350  c1->Update();
351  p->SetLogz();
352 
353  T->Draw("Average_column:Average_HODO_HORIZONTAL>>h2_hodo_c_h", cut_C, "colz");
354 
355  p = (TPad *) c1->cd(idx++);
356  c1->Update();
357 
358  T->Draw("Energy_Sum_RAW_CEMC:Average_HODO_HORIZONTAL>>h2_hodoE_h", cut_C, "colz");
359 
360  p = (TPad *) c1->cd(idx++);
361  c1->Update();
362 
363  T->Draw("Average_HODO_VERTICAL:Average_HODO_HORIZONTAL>>h2_hodoE_v_h", "(" + cut_C + ")*(Energy_Sum_RAW_CEMC) * (Energy_Sum_RAW_CEMC>0.1)", "goff");
364  T->Draw("Average_HODO_VERTICAL:Average_HODO_HORIZONTAL>>h2_hodoE_v_h_norm", "(" + cut_C + " ) && (Energy_Sum_RAW_CEMC>0.1)", "goff");
365  h2_hodoE_v_h->Divide(h2_hodoE_v_h_norm);
366  h2_hodoE_v_h->Draw("colz");
367 
368  p = (TPad *) c1->cd(idx++);
369  c1->Update();
370  p->SetLogz();
371 
372  T->Draw("Average_row:Average_HODO_VERTICAL>>h2_hodo_r_v", cut_C, "colz");
373 
374  p = (TPad *) c1->cd(idx++);
375  c1->Update();
376 
377  T->Draw("Energy_Sum_RAW_CEMC:Average_HODO_VERTICAL>>h2_hodoE_v", cut_C, "colz");
378 
379  SaveCanvas(c1,
380  TString(_file0->GetName()) + TString("_DrawPrototype4EMCalTower_") + TString(c1->GetName()), false);
381 }
382 
383 void EMCDistribution_HCalCalibration()
384 {
385  gStyle->SetOptStat(0);
386  gStyle->SetOptFit(1111);
387 
388  TH1F *h1_HCalOut = new TH1F("h1_HCalOut",
389  "(CEMC MIP && HCal_{IN} MIP);EMCal + HCal_{IN} + HCal_{OUT} (GeV)", 100,
390  0, abs(beam_momentum_selection) * 2.);
391 
392  TH2 *h2_HCalOut_HCalIn =
393  new TH2F("h2_HCalOut_HCalIn",
394  "Energy Balance HCal_{OUT} and HCal_{IN} (CEMC MIP && NOT HCal_{IN} MIP);(HCal_{OUT} - HCal_{IN})/HCal Sum;EMCal + HCal_{IN} + HCal_{OUT} (GeV)",
395  10, -1, 1, 100, 0, abs(beam_momentum_selection) * 2.);
396 
397  TH2 *h2_HCal_EMCal =
398  new TH2F("h2_HCal_EMCal",
399  "Energy Balance HCal_{SUM} and EMCal (NOT CEMC MIP);(HCal_{SUM} - EMCal)/Sum;EMCal + HCal_{IN} + HCal_{OUT} (GeV)",
400  10, -1, 1, 100, 0, abs(beam_momentum_selection) * 2.);
401 
402  TText *t;
403  TCanvas *c1 = new TCanvas("EMCDistribution_HCalCalibration" + cuts,
404  "EMCDistribution_HCalCalibration" + cuts, 1800, 600);
405  c1->Divide(3, 1);
406  int idx = 1;
407  TPad *p;
408 
409  p = (TPad *) c1->cd(idx++);
410  c1->Update();
411 
412  T->Draw(
413  "Energy_Sum_Hadron_CEMC + Energy_Sum_Hadron_HCALIN + Energy_Sum_Hadron_HCALOUT:(Energy_Sum_Hadron_HCALOUT - Energy_Sum_Hadron_HCALIN)/(Energy_Sum_Hadron_HCALIN+Energy_Sum_Hadron_HCALOUT)>>h2_HCalOut_HCalIn",
414  "(!HCALIN_MIP) && CEMC_MIP", "goff");
415  h2_HCalOut_HCalIn->FitSlicesY();
416 
417  TH1 *h2_HCalOut_HCalIn_1 = (TH1 *) gDirectory->GetObjectChecked(
418  "h2_HCalOut_HCalIn_1", "TH1D");
419  assert(h2_HCalOut_HCalIn_1);
420  h2_HCalOut_HCalIn_1->SetMarkerStyle(kFullCircle);
421 
422  TF1 *f_HCalOut_HCalIn = new TF1("f_HCalOut_HCalIn",
423  "[0] * (x+1)/2 + [1]*(1-x)/2", -.8, 0.8);
424  f_HCalOut_HCalIn->SetLineColor(kMagenta);
425  f_HCalOut_HCalIn->SetLineWidth(3);
426  h2_HCalOut_HCalIn_1->Fit(f_HCalOut_HCalIn, "MR0");
427 
428  h2_HCalOut_HCalIn->DrawClone("colz");
429  h2_HCalOut_HCalIn_1->DrawClone("SAME");
430  f_HCalOut_HCalIn->DrawClone("same");
431 
432  p = (TPad *) c1->cd(idx++);
433  c1->Update();
434 
435  T->Draw(
436  "Energy_Sum_Hadron_CEMC + Energy_Sum_Hadron_HCALIN + Energy_Sum_Hadron_HCALOUT:(Energy_Sum_Hadron_HCALIN + Energy_Sum_Hadron_HCALOUT - Energy_Sum_Hadron_CEMC)/(Energy_Sum_Hadron_CEMC + Energy_Sum_Hadron_HCALIN + Energy_Sum_Hadron_HCALOUT)>>h2_HCal_EMCal",
437  "! CEMC_MIP", "goff");
438  h2_HCal_EMCal->FitSlicesY();
439 
440  TH1 *h2_HCal_EMCal_1 = (TH1 *) gDirectory->GetObjectChecked(
441  "h2_HCal_EMCal_1", "TH1D");
442  assert(h2_HCal_EMCal_1);
443  h2_HCal_EMCal_1->SetMarkerStyle(kFullCircle);
444 
445  TF1 *f_HCal_EMCal = new TF1("f_HCal_EMCal", "[0] * (x+1)/2 + [1]*(1-x)/2",
446  -.8, 0.8);
447  f_HCal_EMCal->SetLineColor(kMagenta);
448  f_HCal_EMCal->SetLineWidth(3);
449  h2_HCal_EMCal_1->Fit(f_HCal_EMCal, "MR0");
450 
451  h2_HCal_EMCal->DrawClone("colz");
452  h2_HCal_EMCal_1->DrawClone("SAME");
453  f_HCal_EMCal->DrawClone("same");
454 
455  p = (TPad *) c1->cd(idx++);
456  c1->Update();
457 
458  T->Draw(
459  "Energy_Sum_Hadron_CEMC + Energy_Sum_Hadron_HCALIN + Energy_Sum_Hadron_HCALOUT>>h1_HCalOut",
460  "HCALIN_MIP && CEMC_MIP", "goff");
461 
462  TH1 *h2_HCal_EMCal_ProjY = (TH1 *) h2_HCal_EMCal->ProjectionY(
463  "h2_HCal_EMCal_ProjY", 2, 9)
464  ->DrawClone();
465  TH1 *h2_HCalOut_HCalIn_ProjY = (TH1 *) h2_HCalOut_HCalIn->ProjectionY(
466  "h2_HCalOut_HCalIn_ProjY")
467  ->DrawClone("same");
468  h1_HCalOut->Draw("same");
469 
470  h2_HCal_EMCal_ProjY->SetTitle(
471  "Red: 3 Calo shower, Blue: MIP EMCal, Green: MIP EMCal & HCal_{IN}");
472 
473  h2_HCal_EMCal_ProjY->SetLineColor(kRed + 2);
474  h2_HCal_EMCal_ProjY->SetLineWidth(3);
475  h2_HCalOut_HCalIn_ProjY->SetLineColor(kBlue + 2);
476  h2_HCalOut_HCalIn_ProjY->SetLineWidth(2);
477  h1_HCalOut->SetLineColor(kGreen + 3);
478  h1_HCalOut->SetLineWidth(3);
479 
480  SaveCanvas(c1,
481  TString(_file0->GetName()) + TString("_DrawPrototype4EMCalTower_") + TString(c1->GetName()), false);
482 }
483 
484 void EMCDistribution_ShowShape(TString CherenkovSignal = "C2_Inner",
485  const double che_cut = 2000)
486 {
487  TString cut_pass = CherenkovSignal + Form(">%.1f", che_cut);
488  TString cut_rej = CherenkovSignal + Form("<%.1f", che_cut);
489 
490  const double event_pass = T->GetEntries(cut_pass);
491  const double event_rej = T->GetEntries(cut_rej);
492 
493  TH2 *EnergyDist_pass = new TH2F("EnergyDist_pass",
494  cut_pass + ";Column;Row;<Energy> / Event / Tower", 8, -.5, 7.5, 8, -.5,
495  7.5);
496  TH2 *EnergyDist_rej = new TH2F("EnergyDist_rej",
497  cut_rej + ";Column;Row;<Energy> / Event / Tower", 8, -.5, 7.5, 8, -.5,
498  7.5);
499  TH2 *Hodoscope_dist = new TH2F("Hodoscope_dist",
500  cut_rej + ";7 - Horizontal Hodoscope (5 mm);7 - Vertical Hodoscope (5 mm); Event / finger^2", 8, -.5, 7.5, 8, -.5,
501  7.5);
502 
503  TH1 *Che_full = new TH1F("Che_full",
504  ";" + CherenkovSignal + " Signal (ADC);Count / bin", 200, 0, 17000);
505  TH1 *Che_pass = new TH1F("Che_pass",
506  ";" + CherenkovSignal + " Signal (ADC);Count / bin", 200, 0, 17000);
507  TH1 *Che_rej = new TH1F("Che_rej",
508  ";" + CherenkovSignal + " Signal (ADC);Count / bin", 200, 0, 17000);
509 
510  Che_full->SetLineColor(kBlue + 3);
511  Che_full->SetLineWidth(2);
512  Che_pass->SetLineColor(kGreen + 3);
513  Che_pass->SetLineWidth(2);
514  Che_pass->SetFillColor(kGreen + 3);
515  Che_pass->SetFillStyle(1);
516  Che_rej->SetLineColor(kBlue + 3);
517  Che_rej->SetLineWidth(2);
518  Che_rej->SetFillColor(kBlue + 3);
519  Che_rej->SetFillStyle(1);
520 
521  T->Draw(
522  "TOWER_CALIB_CEMC[].get_binphi():TOWER_CALIB_CEMC[].get_bineta()>>EnergyDist_pass",
523  Form("(%s) * (TOWER_CALIB_CEMC[].get_energy())", cut_pass.Data()),
524  "goff");
525  T->Draw(
526  "TOWER_CALIB_CEMC[].get_binphi():TOWER_CALIB_CEMC[].get_bineta()>>EnergyDist_rej",
527  Form("(%s) * (TOWER_CALIB_CEMC[].get_energy())", cut_rej.Data()), "goff");
528  T->Draw(
529  "7 - Average_HODO_VERTICAL:7 - Average_HODO_HORIZONTAL>>Hodoscope_dist", "1", "goff");
530 
531  T->Draw(CherenkovSignal + ">>Che_full", NULL, "goff");
532  T->Draw(CherenkovSignal + ">>Che_pass", cut_pass, "goff");
533  T->Draw(CherenkovSignal + ">>Che_rej", cut_rej, "goff");
534 
535  EnergyDist_pass->Scale(1. / event_pass);
536  EnergyDist_rej->Scale(1. / event_rej);
537 
538  TText *t;
539  TCanvas *c1 = new TCanvas("EMCDistribution_ShowShape" + cuts,
540  "EMCDistribution_ShowShape" + cuts, 1100, 950);
541  c1->Divide(2, 2);
542  int idx = 1;
543  TPad *p;
544 
545  p = (TPad *) c1->cd(idx++);
546  c1->Update();
547  p->SetLogy();
548  p->SetGridx(0);
549  p->SetGridy(0);
550 
551  Che_full->DrawClone();
552  Che_pass->DrawClone("same");
553  Che_rej->DrawClone("same");
554 
555  p = (TPad *) c1->cd(idx++);
556  c1->Update();
557  // p->SetLogy();
558 
559  EnergyDist_pass->DrawClone("colz");
560 
561  p = (TPad *) c1->cd(idx++);
562  c1->Update();
563  // p->SetLogy();
564  // p->SetGridx(0);
565  // p->SetGridy(0);
566 
567  EnergyDist_rej->DrawClone("colz");
568 
569  p = (TPad *) c1->cd(idx++);
570  c1->Update();
571  // p->SetLogy();
572  // p->SetGridx(0);
573  // p->SetGridy(0);
574 
575  Hodoscope_dist->DrawClone("colz");
576 
577  SaveCanvas(c1,
578  TString(_file0->GetName()) + TString("_DrawPrototype4EMCalTower_") + TString(c1->GetName()), false);
579 }
580 
581 void EMCDistribution_SUM(TString sTOWER = "Energy_Sum_col1_row2_5x5",
582  TString CherenkovSignal = "C2_Inner")
583 {
584  TH1 *EnergySum_LG_full = new TH1F("EnergySum_LG_full",
585  ";Full range Tower Energy Sum (GeV);Count / bin", 900, -5, 120);
586  TH1 *EnergySum_LG = new TH1F("EnergySum_LG",
587  ";Full range Tower Energy Sum (GeV);Count / bin", 900, -5, 120);
588  // TH1 * EnergySum_HG = new TH1F("EnergySum_HG",
589  // ";Low range Tower Energy Sum (ADC);Count / bin", 50, 0, 500);
590 
591  TH1 *C2_Inner_full = new TH1F("C2_Inner_full",
592  CherenkovSignal + ";Cherenkov Signal (ADC);Count / bin", 180, -1000, 17000);
593  TH1 *C2_Inner = new TH1F("C2_Inner",
594  CherenkovSignal + ";Cherenkov Inner Signal (ADC);Count / bin", 180, -1000,
595  17000);
596 
597  EnergySum_LG_full->SetLineColor(kBlue + 3);
598  EnergySum_LG_full->SetLineWidth(2);
599 
600  EnergySum_LG->SetLineColor(kGreen + 3);
601  EnergySum_LG->SetLineWidth(3);
602  EnergySum_LG->SetMarkerColor(kGreen + 3);
603 
604  C2_Inner_full->SetLineColor(kBlue + 3);
605  C2_Inner_full->SetLineWidth(2);
606 
607  C2_Inner->SetLineColor(kGreen + 3);
608  C2_Inner->SetLineWidth(3);
609  C2_Inner->SetMarkerColor(kGreen + 3);
610 
611  TCut c2 = CherenkovSignal + ">2000";
612 
613  T->Draw(sTOWER + ">>EnergySum_LG_full", "", "goff");
614  T->Draw(sTOWER + ">>EnergySum_LG", c2, "goff");
615  T->Draw(CherenkovSignal + ">>C2_Inner_full", "", "goff");
616  T->Draw(CherenkovSignal + ">>C2_Inner", c2, "goff");
617 
618  TText *t;
619  TCanvas *c1 = new TCanvas(
620  "EMCDistribution_SUM_" + sTOWER + "_" + CherenkovSignal + cuts,
621  "EMCDistribution_SUM_" + sTOWER + "_" + CherenkovSignal + cuts, 1800,
622  600);
623  c1->Divide(3, 1);
624  int idx = 1;
625  TPad *p;
626 
627  p = (TPad *) c1->cd(idx++);
628  c1->Update();
629  p->SetLogy();
630  p->SetGridx(0);
631  p->SetGridy(0);
632 
633  C2_Inner_full->DrawClone();
634  C2_Inner->DrawClone("same");
635 
636  p = (TPad *) c1->cd(idx++);
637  c1->Update();
638  p->SetLogy();
639  p->SetGridx(0);
640  p->SetGridy(0);
641 
642  TH1 *h = (TH1 *) EnergySum_LG_full->DrawClone();
643  h->GetXaxis()->SetRangeUser(-1, h->GetMean() + 5 * h->GetRMS());
644  (TH1 *) EnergySum_LG->DrawClone("same");
645 
646  p = (TPad *) c1->cd(idx++);
647  c1->Update();
648  // p->SetLogy();
649  p->SetGridx(0);
650  p->SetGridy(0);
651 
652  TH1 *h_full = (TH1 *) EnergySum_LG_full->DrawClone();
653  // h_full->GetXaxis()->SetRangeUser(0.5,32);
654  TH1 *h = (TH1 *) EnergySum_LG->DrawClone("same");
655 
656  TF1 *fgaus_g = new TF1("fgaus_LG_g", "gaus", h->GetMean() - 1 * h->GetRMS(),
657  h->GetMean() + 4 * h->GetRMS());
658  fgaus_g->SetParameters(1, h->GetMean() - 2 * h->GetRMS(),
659  h->GetMean() + 2 * h->GetRMS());
660  h->Fit(fgaus_g, "MR0N");
661 
662  TF1 *fgaus = new TF1("fgaus_LG", "gaus",
663  fgaus_g->GetParameter(1) - 1 * fgaus_g->GetParameter(2),
664  fgaus_g->GetParameter(1) + 4 * fgaus_g->GetParameter(2));
665  fgaus->SetParameters(fgaus_g->GetParameter(0), fgaus_g->GetParameter(1),
666  fgaus_g->GetParameter(2));
667  h->Fit(fgaus, "MR");
668 
669  h->Sumw2();
670  h_full->Sumw2();
671  h_full->GetXaxis()->SetRangeUser(h->GetMean() - 4 * h->GetRMS(),
672  h->GetMean() + 4 * h->GetRMS());
673  h_full->GetYaxis()->SetRangeUser(0,
674  fgaus_g->GetParameter(0) * 4);
675 
676  h->SetLineWidth(2);
677  h->SetMarkerStyle(kFullCircle);
678 
679  h_full->SetTitle(
680  Form("#DeltaE/<E> = %.1f%%",
681  100 * fgaus->GetParameter(2) / fgaus->GetParameter(1)));
682 
683  // p = (TPad *) c1->cd(idx++);
684  // c1->Update();
685  // p->SetLogy();
686  // p->SetGridx(0);
687  // p->SetGridy(0);
688  //
689  // TH1 * h = (TH1 *) EnergySum_LG->DrawClone();
690  // h->GetXaxis()->SetRangeUser(0,500);
691  // h->SetLineWidth(2);
692  // h->SetLineColor(kBlue + 3);
694  // h->GetXaxis()->SetRangeUser(0, h->GetMean() + 5 * h->GetRMS());
695  //
696  // p = (TPad *) c1->cd(idx++);
697  // c1->Update();
699  // p->SetGridx(0);
700  // p->SetGridy(0);
701  //
702  // TH1 * h = (TH1 *) EnergySum_LG->DrawClone();
703  // h->GetXaxis()->SetRangeUser(0,500);
704  //
705  // TF1 * fgaus = new TF1("fgaus_HG", "gaus", 0, 100);
706  // fgaus->SetParameters(1, h->GetMean() - 2 * h->GetRMS(),
707  // h->GetMean() + 2 * h->GetRMS());
708  // h->Fit(fgaus, "M");
709  //
710  // h->Sumw2();
711  // h->GetXaxis()->SetRangeUser(h->GetMean() - 4 * h->GetRMS(),
712  // h->GetMean() + 4 * h->GetRMS());
713  //
714  // h->SetLineWidth(2);
715  // h->SetMarkerStyle(kFullCircle);
716  //
717  // h->SetTitle(
718  // Form("#DeltaE/<E> = %.1f%%",
719  // 100 * fgaus->GetParameter(2) / fgaus->GetParameter(1)));
720 
721  SaveCanvas(c1,
722  TString(_file0->GetName()) + TString("_DrawPrototype4EMCalTower_") + TString(c1->GetName()), false);
723 }
724 
725 void EMCDistribution_SUM_RawADC(TString sTOWER = "Energy_Sum_col1_row2_5x5",
726  TString CherenkovSignal = "C2_Inner")
727 {
728  TH1 *EnergySum_LG_full = new TH1F("EnergySum_LG_full",
729  ";Tower Energy Sum (ADC);Count / bin", 210, -1000, 20000);
730  TH1 *EnergySum_LG = new TH1F("EnergySum_LG",
731  ";Tower Energy Sum (ADC);Count / bin", 210, -1000, 20000);
732  // TH1 * EnergySum_HG = new TH1F("EnergySum_HG",
733  // ";Low range Tower Energy Sum (ADC);Count / bin", 50, 0, 500);
734 
735  TH1 *C2_Inner_full = new TH1F("C2_Inner_full",
736  CherenkovSignal + ";Cherenkov Signal (ADC);Count / bin", 180, -1000, 17000);
737  TH1 *C2_Inner = new TH1F("C2_Inner",
738  CherenkovSignal + ";Cherenkov Inner Signal (ADC);Count / bin", 180, -1000,
739  17000);
740 
741  EnergySum_LG_full->SetLineColor(kBlue + 3);
742  EnergySum_LG_full->SetLineWidth(2);
743 
744  EnergySum_LG->SetLineColor(kGreen + 3);
745  EnergySum_LG->SetLineWidth(3);
746  EnergySum_LG->SetMarkerColor(kGreen + 3);
747 
748  C2_Inner_full->SetLineColor(kBlue + 3);
749  C2_Inner_full->SetLineWidth(2);
750 
751  C2_Inner->SetLineColor(kGreen + 3);
752  C2_Inner->SetLineWidth(3);
753  C2_Inner->SetMarkerColor(kGreen + 3);
754 
755  TCut c2 = CherenkovSignal + ">2000";
756 
757  T->Draw(sTOWER + ">>EnergySum_LG_full", "", "goff");
758  T->Draw(sTOWER + ">>EnergySum_LG", c2, "goff");
759  T->Draw(CherenkovSignal + ">>C2_Inner_full", "", "goff");
760  T->Draw(CherenkovSignal + ">>C2_Inner", c2, "goff");
761 
762  TText *t;
763  TCanvas *c1 = new TCanvas(
764  "EMCDistribution_SUM_RawADC_" + sTOWER + "_" + CherenkovSignal + cuts,
765  "EMCDistribution_SUM_RawADC_" + sTOWER + "_" + CherenkovSignal + cuts, 1800,
766  600);
767  c1->Divide(3, 1);
768  int idx = 1;
769  TPad *p;
770 
771  p = (TPad *) c1->cd(idx++);
772  c1->Update();
773  p->SetLogy();
774  p->SetGridx(0);
775  p->SetGridy(0);
776 
777  C2_Inner_full->DrawClone();
778  C2_Inner->DrawClone("same");
779 
780  p = (TPad *) c1->cd(idx++);
781  c1->Update();
782  p->SetLogy();
783  p->SetGridx(0);
784  p->SetGridy(0);
785 
786  TH1 *h = (TH1 *) EnergySum_LG_full->DrawClone();
787  // h->GetXaxis()->SetRangeUser(-1, h->GetMean() + 5 * h->GetRMS());
788  (TH1 *) EnergySum_LG->DrawClone("same");
789 
790  p = (TPad *) c1->cd(idx++);
791  c1->Update();
792  // p->SetLogy();
793  p->SetGridx(0);
794  p->SetGridy(0);
795 
796  TH1 *h_full = (TH1 *) EnergySum_LG_full->DrawClone();
797  // h_full->GetXaxis()->SetRangeUser(0.5,32);
798  TH1 *h = (TH1 *) EnergySum_LG->DrawClone("same");
799 
800  TF1 *fgaus_g = new TF1("fgaus_LG_g", "gaus", h->GetMean() - 1 * h->GetRMS(),
801  h->GetMean() + 4 * h->GetRMS());
802  fgaus_g->SetParameters(1, h->GetMean() - 2 * h->GetRMS(),
803  h->GetMean() + 2 * h->GetRMS());
804  h->Fit(fgaus_g, "MR0N");
805 
806  TF1 *fgaus = new TF1("fgaus_LG", "gaus",
807  fgaus_g->GetParameter(1) - 1 * fgaus_g->GetParameter(2),
808  fgaus_g->GetParameter(1) + 4 * fgaus_g->GetParameter(2));
809  fgaus->SetParameters(fgaus_g->GetParameter(0), fgaus_g->GetParameter(1),
810  fgaus_g->GetParameter(2));
811  h->Fit(fgaus, "MR");
812 
813  h->Sumw2();
814  h_full->Sumw2();
815  h_full->GetXaxis()->SetRangeUser(h->GetMean() - 4 * h->GetRMS(),
816  h->GetMean() + 4 * h->GetRMS());
817  h_full->GetYaxis()->SetRangeUser(0,
818  fgaus_g->GetParameter(0) * 4);
819 
820  h->SetLineWidth(2);
821  h->SetMarkerStyle(kFullCircle);
822 
823  h_full->SetTitle(
824  Form("#DeltaE/<E> = %.1f%%",
825  100 * fgaus->GetParameter(2) / fgaus->GetParameter(1)));
826 
827  // p = (TPad *) c1->cd(idx++);
828  // c1->Update();
829  // p->SetLogy();
830  // p->SetGridx(0);
831  // p->SetGridy(0);
832  //
833  // TH1 * h = (TH1 *) EnergySum_LG->DrawClone();
834  // h->GetXaxis()->SetRangeUser(0,500);
835  // h->SetLineWidth(2);
836  // h->SetLineColor(kBlue + 3);
838  // h->GetXaxis()->SetRangeUser(0, h->GetMean() + 5 * h->GetRMS());
839  //
840  // p = (TPad *) c1->cd(idx++);
841  // c1->Update();
843  // p->SetGridx(0);
844  // p->SetGridy(0);
845  //
846  // TH1 * h = (TH1 *) EnergySum_LG->DrawClone();
847  // h->GetXaxis()->SetRangeUser(0,500);
848  //
849  // TF1 * fgaus = new TF1("fgaus_HG", "gaus", 0, 100);
850  // fgaus->SetParameters(1, h->GetMean() - 2 * h->GetRMS(),
851  // h->GetMean() + 2 * h->GetRMS());
852  // h->Fit(fgaus, "M");
853  //
854  // h->Sumw2();
855  // h->GetXaxis()->SetRangeUser(h->GetMean() - 4 * h->GetRMS(),
856  // h->GetMean() + 4 * h->GetRMS());
857  //
858  // h->SetLineWidth(2);
859  // h->SetMarkerStyle(kFullCircle);
860  //
861  // h->SetTitle(
862  // Form("#DeltaE/<E> = %.1f%%",
863  // 100 * fgaus->GetParameter(2) / fgaus->GetParameter(1)));
864 
865  SaveCanvas(c1,
866  TString(_file0->GetName()) + TString("_DrawPrototype4EMCalTower_") + TString(c1->GetName()), false);
867 }
868 
869 void EMCDistribution_Fast(TString gain = "CALIB", bool full_gain = false)
870 {
871  TString hname = "EMCDistribution_" + gain + TString(full_gain ? "_FullGain" : "") + cuts;
872 
873  TH2 *h2 = NULL;
874  if (gain.Contains("CALIB"))
875  {
876  if (full_gain)
877  {
878  h2 = new TH2F(hname,
879  Form(";Calibrated Tower Energy Sum (GeV);Count / bin"), 100, -1,
880  30, 64, -.5, 63.5);
881  // QAHistManagerDef::useLogBins(h2->GetXaxis());
882  }
883  else
884  {
885  h2 = new TH2F(hname,
886  Form(";Calibrated Tower Energy Sum (GeV);Count / bin"), 260, -.2,
887  5, 64, -.5, 63.5);
888  }
889  T->Draw(
890  "TOWER_" + gain + "_CEMC[].get_bineta() + 8* TOWER_" + gain + "_CEMC[].get_binphi():TOWER_" + gain + "_CEMC[].get_energy()>>" + hname, "", "goff");
891  }
892  else if (gain.Contains("RAW"))
893  {
894  if (full_gain)
895  {
896  h2 = new TH2F(hname,
897  Form(";Tower Peak Amp (ADC);Count / bin"), 1 << 8,
898  -3000, 17000, 64, -.5, 63.5);
899  // QAHistManagerDef::useLogBins(h2->GetXaxis());
900  }
901  else
902  {
903  h2 = new TH2F(hname,
904  Form(";Tower Peak Amp (ADC);Count / bin"), 1 << 8,
905  0, 1 << 14, 64, -.5, 63.5);
906  }
907  T->Draw(
908  "TOWER_" + gain + "_CEMC[].get_bineta() + 8* TOWER_" + gain + "_CEMC[].get_binphi():TOWER_" + gain + "_CEMC[].get_energy()>>" + hname, "", "goff");
909  }
910 
911  TText *t;
912  TCanvas *c1 = new TCanvas(
913  "EMCDistribution_" + gain + TString(full_gain ? "_FullGain" : "") + cuts,
914  "EMCDistribution_" + gain + TString(full_gain ? "_FullGain" : "") + cuts,
915  1800, 950);
916  c1->Divide(8, 8, 0., 0.01);
917  int idx = 1;
918  TPad *p;
919 
920  for (int iphi = 8 - 1; iphi >= 0; iphi--)
921  {
922  for (int ieta = 0; ieta < 8; ieta++)
923  {
924  p = (TPad *) c1->cd(idx++);
925  c1->Update();
926 
927  p->SetLogy();
928  // if (full_gain)
929  // {
930  // p->SetLogx();
931  // }
932  p->SetGridx(0);
933  p->SetGridy(0);
934 
935  TString hname = Form("hEnergy_ieta%d_iphi%d", ieta, iphi) + TString(full_gain ? "_FullGain" : "");
936 
937  TH1 *h = h2->ProjectionX(hname, ieta + 8 * iphi + 1,
938  ieta + 8 * iphi + 1); // axis bin number is encoded as ieta+8*iphi+1
939 
940  h->SetLineWidth(0);
941  h->SetLineColor(kBlue + 3);
942  h->SetFillColor(kBlue + 3);
943 
944  h->GetXaxis()->SetTitleSize(.09);
945  h->GetXaxis()->SetLabelSize(.08);
946  h->GetYaxis()->SetLabelSize(.08);
947 
948  h->Draw();
949 
950  if (full_gain)
951  h->Fit("x*gaus", "M");
952  else
953  h->Fit("landau", "M");
954 
955  double peak = -1;
956 
957  TF1 *fit = ((TF1 *) (h->GetListOfFunctions()->At(0)));
958  if (fit)
959  {
960  fit->SetLineColor(kRed);
961  peak = fit->GetParameter(1);
962  }
963 
964  cout << Form("Finished <Col%d Row%d> = %.1f", ieta, iphi, peak)
965  << endl;
966 
967  TText *t = new TText(.9, .9,
968  Form("<Col%d Row%d> = %.1f", ieta, iphi, peak));
969  t->SetTextAlign(33);
970  t->SetTextSize(.15);
971  t->SetNDC();
972  t->Draw();
973  }
974  }
975 
976  SaveCanvas(c1,
977  TString(_file0->GetName()) + TString("_DrawPrototype4EMCalTower_") + TString(c1->GetName()), false);
978 }
979 
980 void EMCDistribution_PeakSample_Fast(bool full_gain = false)
981 {
982  const TString gain = "RAW";
983 
984  TString hname = "EMCDistribution_" + gain + TString(full_gain ? "_FullGain" : "") + cuts;
985 
986  TH2 *h2 = NULL;
987  {
988  if (full_gain)
989  {
990  h2 = new TH2F(hname,
991  Form(";Calibrated Tower Energy Sum (ADC);Count / bin"), 100,
992  .05 * 100, 25 * 100, 64, -.5, 63.5);
993  QAHistManagerDef::useLogBins(h2->GetXaxis());
994  }
995  else
996  {
997  h2 = new TH2F(hname,
998  Form(";Calibrated Tower Energy Sum (ADC);Count / bin"), 260,
999  -.2 * 100, 5 * 100, 64, -.5, 63.5);
1000  }
1001  T->Draw(
1002  "TOWER_" + gain + "_CEMC[].get_bineta() + 8* TOWER_" + gain + "_CEMC[].get_binphi():(TOWER_RAW_CEMC[].signal_samples[10] - TOWER_RAW_CEMC[].signal_samples[0])*(-1)>>" + hname, "", "goff");
1003  }
1004 
1005  TText *t;
1006  TCanvas *c1 = new TCanvas(
1007  "EMCDistribution_PeakSample_Fast_" + TString(full_gain ? "_FullGain" : "") + cuts,
1008  "EMCDistribution_PeakSample_Fast_" + TString(full_gain ? "_FullGain" : "") + cuts, 1800, 950);
1009  c1->Divide(8, 8, 0., 0.01);
1010  int idx = 1;
1011  TPad *p;
1012 
1013  for (int iphi = 8 - 1; iphi >= 0; iphi--)
1014  {
1015  for (int ieta = 0; ieta < 8; ieta++)
1016  {
1017  p = (TPad *) c1->cd(idx++);
1018  c1->Update();
1019 
1020  p->SetLogy();
1021  if (full_gain)
1022  {
1023  p->SetLogx();
1024  }
1025  p->SetGridx(0);
1026  p->SetGridy(0);
1027 
1028  TString hname = Form("hEnergy_ieta%d_iphi%d", ieta, iphi) + TString(full_gain ? "_FullGain" : "");
1029 
1030  TH1 *h = h2->ProjectionX(hname, ieta + 8 * iphi + 1,
1031  ieta + 8 * iphi + 1); // axis bin number is encoded as ieta+8*iphi+1
1032 
1033  h->SetLineWidth(0);
1034  h->SetLineColor(kBlue + 3);
1035  h->SetFillColor(kBlue + 3);
1036 
1037  h->GetXaxis()->SetTitleSize(.09);
1038  h->GetXaxis()->SetLabelSize(.08);
1039  h->GetYaxis()->SetLabelSize(.08);
1040 
1041  h->Draw();
1042 
1043  if (full_gain)
1044  h->Fit("x*gaus", "M");
1045  else
1046  h->Fit("landau", "M");
1047 
1048  double peak = -1;
1049 
1050  TF1 *fit = ((TF1 *) (h->GetListOfFunctions()->At(0)));
1051  if (fit)
1052  {
1053  fit->SetLineColor(kRed);
1054  peak = fit->GetParameter(1);
1055  }
1056 
1057  cout << Form("Finished <Col%d Row%d> = %.1f", ieta, iphi, peak)
1058  << endl;
1059 
1060  TText *t = new TText(.9, .9,
1061  Form("<Col%d Row%d> = %.1f", ieta, iphi, peak));
1062  t->SetTextAlign(33);
1063  t->SetTextSize(.15);
1064  t->SetNDC();
1065  t->Draw();
1066  }
1067  }
1068 
1069  SaveCanvas(c1,
1070  TString(_file0->GetName()) + TString("_DrawPrototype4EMCalTower_") + TString(c1->GetName()), false);
1071 }
1072 
1073 void EMCDistribution(TString gain = "CALIB", bool log_scale = false)
1074 {
1075  TText *t;
1076  TCanvas *c1 = new TCanvas(
1077  "EMCDistribution_" + gain + TString(log_scale ? "_Log" : "") + cuts,
1078  "EMCDistribution_" + gain + TString(log_scale ? "_Log" : "") + cuts, 1800,
1079  1000);
1080  c1->Divide(8, 8, 0., 0.01);
1081  int idx = 1;
1082  TPad *p;
1083 
1084  for (int iphi = 8 - 1; iphi >= 0; iphi--)
1085  {
1086  for (int ieta = 0; ieta < 8; ieta++)
1087  {
1088  p = (TPad *) c1->cd(idx++);
1089  c1->Update();
1090 
1091  p->SetLogy();
1092  p->SetGridx(0);
1093  p->SetGridy(0);
1094 
1095  TString hname = Form("hEnergy_ieta%d_iphi%d", ieta, iphi) + TString(log_scale ? "_Log" : "");
1096 
1097  TH1 *h = NULL;
1098 
1099  if (log_scale)
1100  h = new TH1F(hname,
1101  Form(";Calibrated Tower Energy Sum (GeV);Count / bin"), 300,
1102  5e-3, 3096);
1103  else
1104  // h = new TH1F(hname,
1105  // Form(";Calibrated Tower Energy Sum (GeV);Count / bin"), 196,
1106  // 1900, 2096);
1107  h = new TH1F(hname,
1108  Form(";Calibrated Tower Energy Sum (GeV);Count / bin"), 596,
1109  -96, 500);
1110 
1111  h->SetLineWidth(0);
1112  h->SetLineColor(kBlue + 3);
1113  h->SetFillColor(kBlue + 3);
1114  h->GetXaxis()->SetTitleSize(.09);
1115  h->GetXaxis()->SetLabelSize(.08);
1116  h->GetYaxis()->SetLabelSize(.08);
1117 
1118  if (log_scale)
1119  QAHistManagerDef::useLogBins(h->GetXaxis());
1120 
1121  T->Draw(
1122  "TOWER_" + gain + "_CEMC[].get_energy_power_law_exp()>>" + hname,
1123  Form(
1124  "TOWER_%s_CEMC[].get_bineta()==%d && TOWER_%s_CEMC[].get_binphi()==%d",
1125  gain.Data(), ieta, gain.Data(), iphi),
1126  "");
1127 
1128  TText *t = new TText(.9, .9, Form("Col%d Row%d", ieta, iphi));
1129  t->SetTextAlign(33);
1130  t->SetTextSize(.15);
1131  t->SetNDC();
1132  t->Draw();
1133 
1134  // return;
1135  }
1136  }
1137 
1138  SaveCanvas(c1,
1139  TString(_file0->GetName()) + TString("_DrawPrototype4EMCalTower_") + TString(c1->GetName()), false);
1140 }
1141 
1142 void EMCDistribution_ADC(bool log_scale = true)
1143 {
1144  TString gain = "RAW";
1145 
1146  TText *t;
1147  TCanvas *c1 = new TCanvas(
1148  "EMCDistribution_ADC_" + gain + TString(log_scale ? "_Log" : "") + cuts,
1149  "EMCDistribution_ADC_" + gain + TString(log_scale ? "_Log" : "") + cuts,
1150  1800, 1000);
1151  c1->Divide(8, 8, 0., 0.01);
1152  int idx = 1;
1153  TPad *p;
1154 
1155  for (int iphi = 8 - 1; iphi >= 0; iphi--)
1156  {
1157  for (int ieta = 0; ieta < 8; ieta++)
1158  {
1159  p = (TPad *) c1->cd(idx++);
1160  c1->Update();
1161 
1162  if (log_scale)
1163  {
1164  p->SetLogz();
1165  }
1166  p->SetGridx(0);
1167  p->SetGridy(0);
1168 
1169  TString hname = Form("hEnergy_ieta%d_iphi%d", ieta, iphi) + TString(log_scale ? "_Log" : "");
1170 
1171  TH1 *h = NULL;
1172 
1173  if (log_scale)
1174  h = new TH2F(hname,
1175  Form(";Sample ID;ADC"), 31, -.5,
1176  30.5,
1177  // 128+64, 0, 3096);
1178  1 << 10, 0, 1 << 14);
1179  // else
1180  // h = new TH2F(hname,
1181  // Form(";Calibrated Tower Energy Sum (GeV);Count / bin"), 100,
1182  // -.050, .5,128,0,2048);
1183 
1184  h->SetLineWidth(0);
1185  h->SetLineColor(kBlue + 3);
1186  h->SetFillColor(kBlue + 3);
1187  h->GetXaxis()->SetTitleSize(.09);
1188  h->GetXaxis()->SetLabelSize(.08);
1189  h->GetYaxis()->SetLabelSize(.08);
1190  // h->GetYaxis()->SetRangeUser(2000,3000);
1191 
1192  // if (log_scale)
1193  // QAHistManagerDef::useLogBins(h->GetYaxis());
1194 
1195  TString sdraw = "TOWER_" + gain + "_CEMC[].signal_samples[]:fmod(Iteration$,31)>>" + hname;
1196  TString scut =
1197  Form(
1198  "TOWER_%s_CEMC[].get_bineta()==%d && TOWER_%s_CEMC[].get_binphi()==%d",
1199  gain.Data(), ieta, gain.Data(), iphi);
1200 
1201  cout << "T->Draw(\"" << sdraw << "\",\"" << scut << "\");" << endl;
1202 
1203  T->Draw(sdraw, scut, "colz");
1204 
1205  TText *t = new TText(.9, .9, Form("Col%d Row%d", ieta, iphi));
1206  t->SetTextAlign(33);
1207  t->SetTextSize(.15);
1208  t->SetNDC();
1209  t->Draw();
1210 
1211  // return;
1212  }
1213  }
1214 
1215  SaveCanvas(c1,
1216  TString(_file0->GetName()) + TString("_DrawPrototype4EMCalTower_") + TString(c1->GetName()), false);
1217 }