AliPhysics  8dc8609 (8dc8609)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
CreateEMCALRunQA.C
Go to the documentation of this file.
1 #if !defined(__CINT__) || defined(__MAKECINT__)
2 #include <Riostream.h>
3 #include "AliEMCALGeometry.h"
4 #include <TColor.h>
5 #include <TStyle.h>
6 #include <TSystem.h>
7 #include <TDirectory.h>
8 #include <TF1.h>
9 #include <TFile.h>
10 #include <TH2F.h>
11 #include <TCanvas.h>
12 #include <TGraphErrors.h>
13 #include <TLegend.h>
14 #include <TTree.h>
15 #include <TPRegexp.h>
16 #include <TList.h>
17 #include <TObjString.h>
18 #include <TDatime.h>
19 #include <TError.h>
20 #include <AliLog.h>
21 #endif
22 
37 
38 Int_t DrawOccupancy(Long_t run, TString period, TString pass, TString fTrigger, TString system, TFile* f, TFile* fout, AliEMCALGeometry* geom, Int_t SavePlots);
39 Int_t DrawRun(Long_t run, TString period, TString pass, TString fTrigger, TFile *f, TFile* fout, Int_t SavePlots, Int_t nSM, Bool_t kFilter);
40 Int_t TrendingEMCALTree(Long_t RunId, TString fCalorimeter, TString system, TString period, TString pass, const int n, TList* TriggersList, TFile* f, TFile *fout, Int_t SavePlots);
41 
42 TH2F* FormatRunHisto(TH2F* aHisto, const char* title, const char* YTitle="");
43 TH2F* HistoPerMod(TH2F* name, const char* title);
44 TH2* AutoZoom(TH2* H, Option_t* aType="all", Int_t EntryMin=0);
45 int FindNumberOfSM(TFile* f, TString fTrigger, TString period);
46 
47 TString QAPATH;
48 TString QAPATHF = "./";
49 
50 //-----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
51 void set_plot_style(){
52 
53  const Int_t NRGBs = 5;
54  const Int_t NCont = 255;
55 
56  Double_t stops[NRGBs] = { 0.00, 0.34, 0.61, 0.84, 1.00 };
57  Double_t red[NRGBs] = { 0.00, 0.00, 0.87, 1.00, 0.51 };
58  Double_t green[NRGBs] = { 0.00, 0.81, 1.00, 0.20, 0.00 };
59  Double_t blue[NRGBs] = { 0.51, 1.00, 0.12, 0.00, 0.00 };
60  TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont);
61  gStyle->SetNumberContours(NCont);
62 
63 }
64 
65 //-----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
66 Double_t pi0massP2(Double_t *x, Double_t *par){
67 
68  Double_t gaus;
69 
70  if(par[2] != 0.)
71  gaus = par[0] * TMath::Exp( -(x[0]-par[1])*(x[0]-par[1]) / (2*par[2]*par[2]) );
72  else
73  gaus = 99999999.;
74 
75  Double_t back = par[3] + par[4]*x[0] + par[5]*x[0]*x[0];
76 
77  return gaus+back;
78 
79 }
80 
81 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------
82 Double_t pi0massP1(Double_t *x, Double_t *par){
83 
84  Double_t gaus = par[0] * TMath::Exp( -(x[0]-par[1])*(x[0]-par[1]) / (2*par[2]*par[2]) );
85  Double_t back = par[3] + par[4]*x[0];
86 
87  return gaus+back;
88 
89 }
90 
91 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------
92 Double_t fitE(Double_t *x, Double_t *par){
93 
94  Double_t levy;
95 
96  levy = par[0] * TMath::Exp( -par[1]/x[0]) * TMath::Power(x[0], -par[2]) ;
97 
98  return levy;
99 
100 }
101 
102 //--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
103 int CreateEMCALRunQA(const char* filename, TString RunId, TString period, TString pass, Int_t SavePlots = 0, Bool_t filter=0, TString fTrigger = "", TString system = "", TString fCalorimeter = "EMCAL"){
104 
105  QAPATH = TString(gSystem->Getenv("QAPATH"));
106  if(QAPATH.IsNull())
107  QAPATH = QAPATHF;
108 
109  if(! QAPATH.BeginsWith("./"))
110  QAPATH = QAPATH + RunId + "/";
111 
112  AliLog::SetGlobalLogLevel(AliLog::kError);
113  TFile *f = new TFile(filename);
114  AliLog::SetGlobalLogLevel(AliLog::kInfo);
115 
116  if(f->IsZombie()){
117  Error(__FUNCTION__, Form("Error openning the input file %s", filename));
118  return -1;
119  }
120 
121  TList* TriggersList = new TList();
122 
123  if(fTrigger==""){
124  TPMERegexp name_re("CaloQA_\\w+");
125  TObjLink* link = f->GetListOfKeys()->FirstLink();
126 
127  while(link){
128  TString name = link->GetObject()->GetName();
129  if(name_re.Match(name)){
130  TriggersList->Add(link->GetObject());
131  if(TString(filename).Contains("barrel") && ! name.Contains("default")) TriggersList->Remove(link->GetObject());
132  if(TString(filename).Contains("outer") && ! name.Contains("EMC")) TriggersList->Remove(link->GetObject());
133  }
134  link = link->Next();
135  }
136  }
137  else{
138  TriggersList->Add(new TObjString(fTrigger.Data()));
139  }
140 
141  if(!TriggersList->GetEntries()){
142  Error(__FUNCTION__, "No trigger found!");
143  return -2;
144  }
145 
146  int nSM = FindNumberOfSM(f, ((TObjString*)TriggersList->Last())->GetString(), period);
147 
148  if(nSM<0){
149  Error(__FUNCTION__, "Could not find the number of super modules!");
150  return -3;
151  }
152  Info(__FUNCTION__, Form("%i super modules were discovered", nSM));
153 
154  TString GeomName;
155  if(nSM <= 4) { nSM = 4; GeomName = "EMCAL_FIRSTYEARv1"; }
156  else if(nSM <= 10){ nSM = 10; GeomName = "EMCAL_COMPLETEv1"; }
157  else if(nSM <= 12){ nSM = 12; GeomName = "EMCAL_COMPLETE12SMv1"; }
158  else { nSM = 20; GeomName = "EMCAL_COMPLETE12SMv1_DCAL_8SM";}
159 
160  AliEMCALGeometry *geom = new AliEMCALGeometry(GeomName.Data(), "EMCAL");
161  Info(__FUNCTION__, Form("Using %i super modules and the Geometry %s", nSM, GeomName.Data()));
162 
163  TFile *fout = new TFile(TString( QAPATH + period+"_"+pass + fTrigger+"_"+ (Long_t)RunId.Atoi() +"_QAplots.root").Data(), "RECREATE");
164 
165  if((system.IsNull()) && (period.EndsWith("h")))
166  system = "PbPb";
167 
168  Int_t ret=0;
169  TIter next(TriggersList);
170  while(TObject *obj = next()){
171  fTrigger= TString(obj->GetName());
172  ret -= DrawOccupancy(RunId.Atoi(), period, pass, fTrigger, system, f, fout, geom, SavePlots);
173  ret -= DrawRun(RunId.Atoi(), period, pass, fTrigger, f, fout, SavePlots, nSM, filter);
174  }
175  ret-= TrendingEMCALTree(RunId.Atoi(), fCalorimeter, system, period, pass, nSM, TriggersList, f, fout, SavePlots);
176 
177  f->Close();
178  fout->Close();
179  delete f;
180  delete geom;
181 
182  return ret;
183 }
184 
185 //----------------------------------------------------------------------------------------------------------------------------------------------------------
186 Int_t DrawOccupancy(Long_t run, TString period, TString pass, TString fTrigger, TString system, TFile* f, TFile* fout, AliEMCALGeometry* geom, Int_t SavePlots){
187 
188  set_plot_style();
189  gStyle->SetOptStat(0);
190 
191  TH1::AddDirectory(kFALSE);
192  TH2D *hEnergyMapReal = new TH2D("hEnergyMapReal", "", 96, -48, 48, 120, -0, 120);
193  TH2D *hOccupancyMapReal = new TH2D("hOccupancyMapReal", "", 96, -48, 48, 120, -0, 120);
194 
195  hEnergyMapReal->SetXTitle("#eta (bin)");
196  hEnergyMapReal->SetYTitle("#varphi (bin)");
197  hEnergyMapReal->SetZTitle("E (GeV)/event");
198  hEnergyMapReal->GetYaxis()->SetTitleOffset(1.2);
199  hEnergyMapReal->GetZaxis()->SetLabelSize(0.02);
200  hEnergyMapReal->GetZaxis()->SetTitleOffset(1.36);
201 
202  hOccupancyMapReal->SetXTitle("#eta (bin)");
203  hOccupancyMapReal->SetYTitle("#varphi (bin)");
204  hOccupancyMapReal->GetYaxis()->SetTitleOffset(1.2);
205  hOccupancyMapReal->GetZaxis()->SetLabelSize(0.02);
206 
207  Int_t nSupMod, nModule, nIphi, nIeta;
208  Int_t iphi, ieta;
209  Int_t realbineta=0;
210  Int_t realbinphi=0;
211 
212  // NO MASK
213  Int_t mask[1] = {2222222};
214 
215  TH2F *hCellAmplitude;
216  TH1F *hNEvents;
217  Int_t Events;
218  Int_t n=0;
219 
220  TString direct;
221  if(!fTrigger.Contains("QA")){
222  direct = "CaloQA_";
223  }
224  direct += fTrigger;
225 
226  Bool_t dirok = f->cd(direct);
227  if(!dirok){
228  Error(__FUNCTION__, Form("No input drectory %s", direct.Data()));
229  return -1;
230  }
231 
232  TList *outputList = (TList*)gDirectory->Get(direct);
233  if(!outputList){
234  Error(__FUNCTION__, "No input list! ");
235  return -1;
236  }
237  outputList->SetOwner();
238 
239  fout->mkdir(Form("%s/%s/%ld/%s/%s", period.Data(), pass.Data(), run, "RunLevelQA", fTrigger.Data()));
240  fout->cd();
241  fout->Cd(Form("%s/%s/%ld/%s/%s", period.Data(), pass.Data(), run, "RunLevelQA", fTrigger.Data()));
242 
243  hNEvents = (TH1F *)outputList->FindObject("hNEvents");
244  if(!hNEvents){
245  Error(__FUNCTION__, Form("hNEvent histogram not found for trigger %s ", fTrigger.Data()));
246  return -2;
247  }
248 
249  Events = (Int_t)hNEvents->GetEntries();
250  if(Events==0){
251  Error(__FUNCTION__, Form("No event in trigger %s", fTrigger.Data()));
252  return -3;
253  }
254 
255  TH1F* hE = (TH1F *)outputList->FindObject("EMCAL_hE");
256  if(!hE || (!hE->GetEntries())){
257  Error(__FUNCTION__, Form("hE Histogram not found or it is empty for trigger %s, EMCal was not in this run?", fTrigger.Data()));
258  return -4;
259  }
260 
261  Double_t Eth=1;
262  if(system=="PbPb"){
263  Eth = 5.;
264  if(fTrigger.Contains("EMC"))
265  Eth=20.;
266  }
267  if(system=="pp"){
268  Eth = 1.;
269  if(fTrigger.Contains("EMC"))
270  Eth=5.;
271  }
272 
273  hCellAmplitude = (TH2F *)outputList->FindObject("EMCAL_hAmpId");
274 
275  for(Int_t i = 0; i < geom->GetNCells() ; i++){
276  Double_t Esum = 0;
277  Double_t Nsum = 0;
278 
279  for(Int_t j = 1; j <= hCellAmplitude->GetNbinsX(); j++){
280  Double_t E = hCellAmplitude->GetXaxis()->GetBinCenter(j);
281  Double_t N = hCellAmplitude->GetBinContent(j, i+1);
282 
283  if(E < 0.3)
284  continue;
285 
286  if(E <= Eth){
287  Esum += E*N;
288  Nsum += N;
289  }
290  }
291 
292  Int_t absId = i;
293  if(n!=0){
294  if(mask[n]<=mask[n-1])
295  Warning(__FUNCTION__, "The list of bad cells is not sorted !!");
296  }
297 
298  if(i==mask[n]){
299  n++ ;
300  continue;
301  } // skip bad cells
302 
303  geom->GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
304  geom->GetCellPhiEtaIndexInSModule(nSupMod, nModule, nIphi, nIeta, iphi, ieta);
305 
306  realbinphi = 120-(nSupMod/2)*24 -iphi -1;
307  if(nSupMod%2==0)
308  realbineta= 48-ieta -1;
309  if(nSupMod%2==1)
310  realbineta= -ieta -1;
311 
312  hEnergyMapReal->Fill(realbineta, realbinphi, Esum/(Double_t)Events);
313  // hOccupancyMapReal->Fill(realbineta, realbinphi, Nsum/(Double_t)Events);
314  }
315 
316  cout <<" Run: " << run << " trigger: " << fTrigger << " N_events: "<<Events<<endl;
317 
318  TPMERegexp r("_\\w+");
319  TString Energy; Energy = QAPATH + "MapEnergy" + fTrigger(r) + ".pdf";
320  TString Energy2; Energy2 = QAPATH + "MapEnergy" + fTrigger(r) + ".png";
321  TString Entries; Entries = QAPATH + "MapEntries" + fTrigger(r) + ".pdf";
322  TString Entries2; Entries2 = QAPATH + "MapEntries" + fTrigger(r) + ".png";
323 
324  TCanvas *c1 = new TCanvas("Energymap", "Energy Map", 600, 600);
325  c1->SetFillColor(0);
326  c1->SetGrid();
327  c1->SetRightMargin(0.14);
328 
329  TString title = "run ";
330  title += run ;
331  if(fTrigger.Contains("EMC"))
332  title += " EMC ";
333  else
334  title += " MB ";
335  title += " Summed energy map";
336 
337  hEnergyMapReal->SetTitle(title);
338  AutoZoom(hEnergyMapReal, "miny")->DrawCopy("colz");
339 
340  if(nSupMod<=12){
341  if(SavePlots==2)
342  c1->SaveAs(Energy);
343  if(SavePlots)
344  c1->SaveAs(Energy2);
345  c1->Write();
346  }
347  delete c1;
348 
349  // TCanvas *c2 = new TCanvas("Occupancy", "Occupancy Map", 600, 600);
350  // c2->SetFillColor(0);
351  // c2->SetGrid();
352  // c2->SetRightMargin(0.14);
353  // TString title2 = "run ";
354  // title2 += run ;
355  // if(fTrigger.Contains("EMC")){ title2 += " EMC ";} else{ title2 += " MB ";}
356  // title2 += " Occupancy map";
357  // hOccupancyMapReal->SetTitle(title2);
358  // AutoZoom(hOccupancyMapReal, "miny")->DrawCopy("colz");
359  // if(SavePlots==2) c2->SaveAs(Entries);
360  // if(SavePlots) c2->SaveAs(Entries2);
361  // c2->Write();
362  // delete c2;
363 
364  if(outputList)
365  outputList->Delete();
366  delete hEnergyMapReal;
367  delete hOccupancyMapReal;
368 
369  return 0;
370 }
371 
372 //-----------------------------------------------------------------------------------------------------------------------
373 Int_t DrawRun(const Long_t run, TString period, TString pass, TString fTrigger, TFile *f, TFile *fout, Int_t SavePlots, Int_t nSM, Bool_t kFilter){
374 
375  TString direct;
376  if(!fTrigger.Contains("QA")){
377  direct = "CaloQA_";
378  }
379  direct += fTrigger;
380 
381  f->cd(direct);
382  if(!direct){
383  Error(__FUNCTION__, Form("No input directory %s", direct.Data()));
384  return -1;
385  }
386 
387  TList *outputList = (TList*)gDirectory->Get(direct);
388  if(!outputList){
389  Error(__FUNCTION__, Form("No input list! %s", direct.Data()));
390  return -2;
391  }
392  outputList->SetOwner();
393 
394  if(kFilter){
395  fout->mkdir(Form("%s/%s/%ld/%s", period.Data(), pass.Data(), run, fTrigger.Data()));
396  fout->cd();
397  fout->Cd(Form("%s/%s/%ld/%s", period.Data(), pass.Data(), run, fTrigger.Data()));
398  outputList->Write();
399  }
400  fout->cd();
401  fout->Cd(Form("%s/%s/%ld/%s/%s", period.Data(), pass.Data(), run, "RunLevelQA", fTrigger.Data()));
402 
403  set_plot_style();
404  gStyle->SetPalette(1);
405  set_plot_style();
406  // gStyle->SetOptStat(1);
407 
408  TH1::AddDirectory(kFALSE);
409  TString outfilename;
410  TString outfilename2;
411  const char* legend="";
412  TPMERegexp r("_\\w+");
413 
414  if(fTrigger.Contains("EMC"))
415  legend = Form(" Run %d EMC ", (int)run);
416  else
417  legend = Form(" Run %d MB ", (int)run);
418 
419  TH1F* hNEvents =(TH1F *)outputList->FindObject("hNEvents");
420  if(!hNEvents){
421  Error(__FUNCTION__, Form("hNEvent histogram not found for trigger %s ", fTrigger.Data()));
422  return -3;
423  }
424 
425  Int_t Events = (Int_t)hNEvents->GetEntries();
426  if(Events==0){
427  Error(__FUNCTION__, Form("No event in trigger %s", fTrigger.Data()));
428  return -4 ;
429  }
430 
431  TH1F* hE =(TH1F *)outputList->FindObject("EMCAL_hE");
432  if(!hE || (!hE->GetEntries())){
433  Error(__FUNCTION__, Form("hE Histogram not found or it is empty for trigger %s, EMCal was not in this run?", fTrigger.Data()));
434  return -4;
435  }
436 
437  TCanvas* c00 = new TCanvas("Occupancy map", "Occupancy", 600, 600);
438  c00->SetLogz();
439  c00->SetFillColor(0);
440  c00->SetBorderSize(0);
441  c00->SetFrameBorderMode(0);
442 
443  TH2F* hOccupancyMapReal =(TH2F *)outputList->FindObject("EMCAL_hEtaPhi");
444  if(!hOccupancyMapReal){
445  Error(__FUNCTION__, Form("EMCAL_hEtaPhi: Histogram for trigger %s not found!", fTrigger.Data()));
446  return -5;
447  }
448  FormatRunHisto(hOccupancyMapReal, Form("Occupancy%s", legend), "#varphi (bin)");
449  hOccupancyMapReal->SetXTitle("#eta (bin)");
450 
451  AutoZoom(hOccupancyMapReal, "all")->DrawCopy("colz");
452  outfilename = QAPATH + "OccupancyMap" + fTrigger(r) + ".pdf" ;
453  outfilename2 = QAPATH + "OccupancyMap" + fTrigger(r) + ".png" ;
454 
455  if(SavePlots==2)
456  c00->SaveAs(outfilename);
457  if(SavePlots)
458  c00->SaveAs(outfilename2);
459  c00->Write();
460  delete c00;
461 
462  TCanvas* c0 = new TCanvas("Grid Cell", "Occupancy from Grid Cells", 600, 600);
463  c0->SetLogz();
464  c0->SetFillColor(0);
465  c0->SetBorderSize(0);
466  c0->SetFrameBorderMode(0);
467 
468  TH2F* hGridCells =(TH2F *)outputList->FindObject("EMCAL_hGridCells");
469  if(!hGridCells){
470  Error(__FUNCTION__, Form("EMCAL_hGridCells: Histogram for trigger %s not found!", fTrigger.Data()));
471  return -5;
472  }
473  FormatRunHisto(hGridCells, Form("Occupancy from Grid Cells%s", legend), "#varphi row columns");
474  hGridCells->SetXTitle("#eta row columns");
475  AutoZoom(hGridCells, "all")->DrawCopy("colz");
476 
477  outfilename = QAPATH + "GridCells" + fTrigger(r) + ".pdf" ;
478  outfilename2 = QAPATH + "GridCells" + fTrigger(r) + ".png" ;
479 
480  if(SavePlots==2)
481  c0->SaveAs(outfilename);
482  if(SavePlots)
483  c0->SaveAs(outfilename2);
484  c0->Write();
485  delete c0;
486 
487  TCanvas* c1 = new TCanvas("TimeVsE", "Cluster Time vs. Energy", 600, 600);
488  c1->SetLogz();
489  c1->SetFillColor(0);
490  c1->SetBorderSize(0);
491  c1->SetFrameBorderMode(0);
492 
493  TH2F* hClusterTimeEnergy =(TH2F *)outputList->FindObject("EMCAL_hClusterTimeEnergy");
494  if(!hClusterTimeEnergy){
495  Error(__FUNCTION__, Form("EMCAL_hClusterTimeEnergy: Histogram for trigger %s not found!", fTrigger.Data()));
496  return -5;
497  }
498  FormatRunHisto(hClusterTimeEnergy, Form("Time vs. Energy%s", legend), "EMCal ToF (ns)");
499  AutoZoom(hClusterTimeEnergy, "maxx")->DrawCopy("colz");
500 
501  outfilename = QAPATH + "TimeRun" + fTrigger(r) + ".pdf" ;
502  outfilename2 = QAPATH + "TimeRun" + fTrigger(r) + ".png" ;
503 
504  if(SavePlots==2)
505  c1->SaveAs(outfilename);
506  if(SavePlots)
507  c1->SaveAs(outfilename2);
508  c1->Write();
509  delete c1;
510 
511  if(pass!="simu"){
512  TCanvas* c1a = new TCanvas("TimeVsAbsId", "Cell Id vs. time ", 600, 600);
513  c1a->SetLogz();
514  c1a->SetFillColor(0);
515  c1a->SetBorderSize(0);
516  c1a->SetFrameBorderMode(0);
517 
518  TH2F* hTimeAbsId =(TH2F *)outputList->FindObject("EMCAL_hTimeId");
519  if(!hTimeAbsId){
520  Error(__FUNCTION__, Form("EMCAL_hTimeId: Histogram for trigger %s not found!", fTrigger.Data()));
521  return -5;
522  }
523  FormatRunHisto(hTimeAbsId, Form("Cell ID vs. Time%s", legend), "Cell AbsID");
524  hTimeAbsId->RebinY(16);
525  hTimeAbsId->SetMinimum(32);
526  AutoZoom(hTimeAbsId, "maxx")->DrawCopy("colz");
527 
528  outfilename = QAPATH + "TimeAbsIdRun" + fTrigger(r) + ".pdf" ;
529  outfilename2 = QAPATH + "TimeAbsRun" + fTrigger(r) + ".png" ;
530 
531  if(SavePlots==2)
532  c1a->SaveAs(outfilename);
533  if(SavePlots)
534  c1a->SaveAs(outfilename2);
535  c1a->Write();
536  delete c1a;
537  }
538 
539  // TCanvas* c1aLG = new TCanvas("TimeVsAbsIdLG", "Cell Id vs. time LG ", 600, 600);
540  // c1aLG->SetLogz();
541  // c1aLG->SetFillColor(0);
542  // c1aLG->SetBorderSize(0);
543  // c1aLG->SetFrameBorderMode(0);
544 
545  // TH2F* hTimeAbsIdLG =(TH2F *)outputList->FindObject("EMCAL_hTimeIdLG");
546  // if(!hTimeAbsIdLG){
547  // Error(__FUNCTION__, Form("EMCAL_hTimeIdLG: Histogram for trigger %s not found!", fTrigger.Data()));
548  // return -5;
549  // }
550  // FormatRunHisto(hTimeAbsIdLG, Form("Cell ID vs. Time LG%s", legend), "cell AbsID");
551 
552  // AutoZoom(hTimeAbsIdLG, "maxx")->DrawCopy("colz");
553  // outfilename = QAPATH + "TimeAbsIdRunLG" + fTrigger(r) + ".pdf" ;
554  // outfilename2 = QAPATH + "TimeAbsRunLG" + fTrigger(r) + ".png" ;
555 
556  // if(SavePlots==2) c1a->SaveAs(outfilename);
557  // if(SavePlots) c1a->SaveAs(outfilename2);
558  // c1a->Write();
559  // delete c1a;
560 
561  TCanvas* c1b = new TCanvas("EVsAbsId", "Cell Id vs. Energy", 600, 600);
562  c1b->SetLogz();
563  c1b->SetFillColor(0);
564  c1b->SetBorderSize(0);
565  c1b->SetFrameBorderMode(0);
566 
567  TH2F* hEAbsId =(TH2F *)outputList->FindObject("EMCAL_hAmpId");
568  if(!hEAbsId){
569  Error(__FUNCTION__, Form("EMCAL_hAmpId: Histogram for trigger %s not found!", fTrigger.Data()));
570  return -5;
571  }
572  FormatRunHisto(hEAbsId, Form("Cell ID vs. Energy%s", legend), "Cell AbsID");
573  hEAbsId->RebinY(16);
574  hEAbsId->SetMinimum(32);
575  AutoZoom(hEAbsId, "maxx")->DrawCopy("colz");
576 
577  outfilename = QAPATH + "EAbsIdRun" + fTrigger(r) + ".pdf" ;
578  outfilename2 = QAPATH + "EAbsIdRun" + fTrigger(r) + ".png" ;
579 
580  if(SavePlots==2)
581  c1b->SaveAs(outfilename);
582  if(SavePlots)
583  c1b->SaveAs(outfilename2);
584  c1b->Write();
585  delete c1b;
586 
587  // TCanvas* c1bLG = new TCanvas("EVsAbsIdLG", "Cell Id vs. E LG", 600, 600);
588  // c1bLG->SetLogz();
589  // c1bLG->SetFillColor(0);
590  // c1bLG->SetBorderSize(0);
591  // c1bLG->SetFrameBorderMode(0);
592 
593  // TH2F* hEAbsIdLG =(TH2F *)outputList->FindObject("EMCAL_hAmpIdLG");
594  // if(!hEAbsIdLG){
595  // Error(__FUNCTION__, Form("EMCAL_hAmpIdLG: Histogram for trigger %s not found!", fTrigger.Data()));
596  // return -5;
597  // }
598  // FormatRunHisto(hEAbsIdLG, Form("Cell ID vs. Energy%s", legend), "Cell AbsID ");
599 
600  // AutoZoom(hEAbsIdLG, "maxx")->DrawCopy("colz");
601  // outfilename = QAPATH + "EAbsIdRunLG" + fTrigger(r) + ".pdf" ;
602  // outfilename2 = QAPATH + "EAbsIdRunLG" + fTrigger(r) + ".png" ;
603 
604  // if(SavePlots==2) c1bLG->SaveAs(outfilename);
605  // if(SavePlots) c1bLG->SaveAs(outfilename2);
606  // c1bLG->Write();
607  // delete c1bLG;
608 
609  if(pass!="muon_calo_pass1"){
610  TCanvas * c2 = new TCanvas("ClusterVsTrack ", "Correlation calo Mult vs. Track Multiplicity", 600, 600);
611  c2->SetLogz();
612  c2->SetFillColor(0);
613  c2->SetBorderSize(0);
614  c2->SetFrameBorderMode(0);
615 
616  TH2F* hClusterVsTrack =(TH2F *)outputList->FindObject("EMCAL_hCaloTrackMNClusters");
617  FormatRunHisto(hClusterVsTrack, Form("#it{N} clusters vs. #it{N} tracks%s", legend), "#it{N} clusters");
618  AutoZoom(hClusterVsTrack, "maxx, maxy", 1)->DrawCopy("colz");
619 
620  outfilename = QAPATH + "CaloTrackMult" + fTrigger(r) + ".pdf";
621  outfilename2 = QAPATH + "CaloTrackMult" + fTrigger(r) + ".png";
622 
623  // if(pass!="muon_calo_pass1"){
624  if(SavePlots==2)
625  c2->SaveAs(outfilename);
626  if(SavePlots)
627  c2->SaveAs(outfilename2);
628  c2->Write();
629  delete c2;
630 
631  TCanvas* c3 = new TCanvas("ClusterEVsTrack ", "Correlation #it{E} calo vs. Track Multiplicity", 600, 600);
632  c3->SetLogz();
633  c3->SetFillColor(0);
634  c3->SetBorderSize(0);
635  c3->SetFrameBorderMode(0);
636 
637  TH2F* hClusterEVsTrack =(TH2F*)outputList->FindObject("EMCAL_hCaloTrackMEClusters");
638  FormatRunHisto(hClusterEVsTrack, Form("Sum #it{E} clusters vs. #it{N} tracks%s", legend), "#Sigma#it{E} clusters");
639  AutoZoom(hClusterEVsTrack, "maxx, maxy", 1)->DrawCopy("colz");
640 
641  outfilename = QAPATH + "ETrackMult" + fTrigger(r) + ".pdf";
642  outfilename2 = QAPATH + "ETrackMult" + fTrigger(r) + ".png";
643 
644  // if(pass!="muon_calo_pass1"){
645  if(SavePlots==2)
646  c3->SaveAs(outfilename);
647  if(SavePlots)
648  c3->SaveAs(outfilename2);
649  c3->Write();
650  // }
651  delete c3;
652  }
653 
654  TCanvas* c4 = new TCanvas("ClusterEVsV0 ", "Correlation #it{E} calo vs. V0 signal", 600, 600);
655  c4->SetLogz();
656  c4->SetFillColor(0);
657  c4->SetBorderSize(0);
658  c4->SetFrameBorderMode(0);
659 
660  TH2F* hClusterEVsV0S =(TH2F*)outputList->FindObject("EMCAL_hCaloV0SEClusters");
661  FormatRunHisto(hClusterEVsV0S, Form("Sum #it{E} clusters vs. V0 signal%s", legend), "#Sigma#it{E} clusters");
662  AutoZoom(hClusterEVsV0S, "maxx, maxy", 1)->DrawCopy("colz");
663 
664  outfilename = QAPATH +"EVsV0s" + fTrigger(r) + ".pdf";
665  outfilename2 = QAPATH +"EVsV0s" + fTrigger(r) + ".png";
666 
667  if(SavePlots==2)
668  c4->SaveAs(outfilename);
669  if(SavePlots)
670  c4->SaveAs(outfilename2);
671  c4->Write();
672  delete c4;
673 
674  TCanvas* c6 = new TCanvas("SumCellEvsSM", "sum cell energy vs. SM", 600, 600);
675  c6->SetLogz();
676  c6->SetFillColor(0);
677  c6->SetBorderSize(0);
678  c6->SetFrameBorderMode(0);
679 
680  TH2F* hSumCellEVsSM =(TH2F*)outputList->FindObject("EMCAL_hSumCellsAmp_Mod");
681  FormatRunHisto(hSumCellEVsSM, Form("Sum #it{E} cells vs. SM%s", legend), "SM number");
682  AutoZoom(hSumCellEVsSM, "maxx, maxy", 1)->DrawCopy("colz");
683 
684  outfilename = QAPATH +"SumCellEVsSM" + fTrigger(r) + ".pdf";
685  outfilename2 = QAPATH +"SumCellEVsSM" + fTrigger(r) + ".png";
686 
687  if(SavePlots==2)
688  c6->SaveAs(outfilename);
689  if(SavePlots)
690  c6->SaveAs(outfilename2);
691  c6->Write();
692  delete c6;
693 
694  TCanvas* c5 = new TCanvas("CellsperCluster", "Nb of cells per cluster for each SM", 600, 600);
695  c5->SetLogz();
696  c5->SetFillColor(0);
697  c5->SetBorderSize(0);
698  c5->SetFrameBorderMode(0);
699 
700  Bool_t mod3=0;
701  if(nSM%3)
702  mod3=1;
703  c5->Divide(3, (nSM/3)+mod3);
704 
705  for(int ism = 0; ism < nSM; ism++){
706  c5->cd(ism+1);
707  gPad->SetLogz();
708 
709  if(TString(Form("Nb cells per clus%s Mod %d", legend, ism)).Length() > 60){
710  Error(__FUNCTION__, "Title too long!");
711  return -6;
712  }
713 
714  AutoZoom(HistoPerMod((TH2F*)outputList->FindObject(Form("EMCAL_hNCellsPerCluster_Mod%i", ism)), Form("Nb of cells per cluster%s Mod %d", legend, ism)), "all", 1)->DrawCopy("colz");
715  }
716 
717  outfilename = QAPATH + "CellsperClusterSM" + fTrigger(r) + ".pdf";
718  outfilename2 = QAPATH + "CellsperClusterSM" + fTrigger(r) + ".png";
719 
720  if(SavePlots==2)
721  c5->SaveAs(outfilename);
722  if(SavePlots)
723  c5->SaveAs(outfilename2);
724  c5->Write();
725  delete c5;
726  if(outputList)
727  outputList->Delete();
728 
729  return 0;
730 }
731 
732 //----------------------------------------------------------------------------------------------------------------------------------
733 Int_t TrendingEMCALTree(Long_t RunId, TString fCalorimeter, TString system, TString period, TString pass, int n, TList* TriggersList, TFile* f, TFile *fout, Int_t SavePlots){
734 
735  TString fTrigger="";
736  TString aCalorimeter;
737  if(n<=12)
738  aCalorimeter = fCalorimeter;
739  else
740  aCalorimeter = TString("EMCAL_and_DCAL");
741 
742  TDatime now;
743 
744  Double_t Nevent=0 ;
745  Double_t xe=0.5;
746  Double_t NTClusters=0; // total number of clusters
747  Double_t NTClustersRMS=0;
748  Double_t NCClusters=0; // number of charged clusters
749  Double_t NCClustersRMS=0;
750  Double_t NMatchClustersP=0;
751  Double_t NMatchClustersPRMS=0;
752 
753  Double_t CellMean=0;
754  Double_t CellRMS=0;
755  Double_t ClusterMean=0;
756  Double_t ClusterRMS=0;
757  Double_t EtotalMean=0;
758  Double_t EtotalRMS=0;
759 
760  Double_t CellPerClusterMean=0;
761  Double_t CellPerClusterRMS=0;
762 
763  Double_t mPDG = 134.9766;
764  Double_t Npi0EMCAL=0;
765  Double_t Npi0EMCALErr=0;
766  Double_t MeanPosEMCAL=0;
767  Double_t MeanPosEMCALErr=0;
768  Double_t WidthEMCAL=0;
769  Double_t WidthEMCALErr=0;
770  Double_t Chi2NdfPi0EMCAL=0;
771  Double_t NggEMCAL=0;
772  Double_t NggEMCALErr=0;
773  Double_t SignifEMCAL=0; // !S/(S+B)
774  Double_t SignifEMCALErr=0; // !S/(S+B)
775 
776  Double_t Npi0DCAL=0;
777  Double_t Npi0DCALErr=0;
778  Double_t MeanPosDCAL=0;
779  Double_t MeanPosDCALErr=0;
780  Double_t WidthDCAL=0;
781  Double_t WidthDCALErr=0;
782  Double_t Chi2NdfPi0DCAL=0;
783  Double_t NggDCAL=0;
784  Double_t NggDCALErr=0;
785  Double_t SignifDCAL=0; // !S/(S+B)
786  Double_t SignifDCALErr=0; // !S/(S+B)
787 
788  TFile* ftree = new TFile(Form("%s/trending.root", QAPATH.Data()), "RECREATE");
789 
790  TTree *tree = new TTree("trending", "Trending QA Tree");
791  tree->Branch("fDate", &now);
792  tree->Branch("fCalorimeter", &aCalorimeter);
793  tree->Branch("system", &system);
794  tree->Branch("period", &period);
795  tree->Branch("pass", &pass);
796  tree->Branch("fTrigger", &fTrigger);
797  tree->Branch("run", &RunId, "run/I");
798  tree->Branch("xe", &xe, "xe/D");
799 
800  tree->Branch("Nevent", &Nevent, "Nevent/D");
801  tree->Branch("NMatchClustersP", &NMatchClustersP, "NMatchClustersP/D");
802  tree->Branch("NMatchClustersPRMS", &NMatchClustersPRMS, "NMatchClustersPRMS/D");
803  tree->Branch("CellMean", &CellMean, "CellMean/D");
804  tree->Branch("CellRMS", &CellRMS, "CellRMS/D");
805  tree->Branch("ClusterMean", &ClusterMean, "ClusterMean/D");
806  tree->Branch("ClusterRMS", &ClusterRMS, "ClusterRMS/D");
807  tree->Branch("EtotalMean", &EtotalMean, "EtotalMean/D");
808  tree->Branch("EtotalRMS", &EtotalRMS, "EtotalRMS/D");
809 
810  tree->Branch("CellPerClusterMean", &CellPerClusterMean, "CellPerClusterMean/D");
811  tree->Branch("CellPerClusterRMS", &CellPerClusterRMS, "CellPerClusterRMS/D");
812 
813  tree->Branch("Npi0EMCAL", &Npi0EMCAL, "Npi0EMCAL/D");
814  tree->Branch("Npi0EMCALErr", &Npi0EMCALErr, "Npi0EMCALErr/D");
815  tree->Branch("Npi0DCAL", &Npi0DCAL, "Npi0DCAL/D");
816  tree->Branch("Npi0DCALErr", &Npi0DCALErr, "Npi0DCALErr/D");
817  tree->Branch("MeanPosEMCAL", &MeanPosEMCAL, "MeanPosEMCAL/D");
818  tree->Branch("MeanPosEMCALErr", &MeanPosEMCALErr, "MeanPosEMCALErr/D");
819  tree->Branch("MeanPosDCAL", &MeanPosDCAL, "MeanPosDCAL/D");
820  tree->Branch("MeanPosDCALErr", &MeanPosDCALErr, "MeanPosDCALErr/D");
821  tree->Branch("WidthEMCAL", &WidthEMCAL, "WidthEMCAL/D");
822  tree->Branch("WidthEMCALErr", &WidthEMCALErr, "WidthEMCALErr/D");
823  tree->Branch("WidthDCAL", &WidthDCAL, "WidthDCAL/D");
824  tree->Branch("WidthDCALErr", &WidthDCALErr, "WidthDCALErr/D");
825  tree->Branch("Chi2NdfPi0EMCAL", &Chi2NdfPi0EMCAL, "Chi2NdfPi0EMCAL/D");
826  tree->Branch("Chi2NdfPi0DCAL", &Chi2NdfPi0DCAL, "Chi2NdfPi0DCAL/D");
827  tree->Branch("NggEMCAL", &NggEMCAL, "NggEMCAL/D");
828  tree->Branch("NggEMCALErr", &NggEMCALErr, "NggEMCALErr/D");
829  tree->Branch("NggDCAL", &NggDCAL, "NggDCAL/D");
830  tree->Branch("NggDCALErr", &NggDCALErr, "NggDCALErr/D");
831  tree->Branch("SignifEMCAL", &SignifEMCAL, "SignifEMCAL/D");
832  tree->Branch("SignifDCALErr", &SignifDCALErr, "SignifDCALErr/D");
833  tree->Branch("SignifEMCAL", &SignifEMCAL, "SignifEMCAL/D");
834  tree->Branch("SignifDCALErr", &SignifDCALErr, "SignifDCALErr/D");
835 
836  tree->Branch("nSM", &n, "nSM/I");
837 
838  int nMax = 22;
839  Double_t CellMeanSM[nMax];
840  Double_t CellRMSSM[nMax];
841  Double_t ClusterMeanSM[nMax];
842  Double_t ClusterTotSM[nMax];
843  Double_t ClusterRMSSM[nMax];
844  Double_t EtotalMeanSM[nMax]; // mean total energy deposited per event
845  Double_t EtotalRMSSM[nMax];
846  Double_t CellPerClusterMeanSM[nMax];
847  Double_t CellPerClusterRMSSM[nMax];
848  Double_t ECell1MeanSM[nMax]; // total energy deposited per event without 1-cell clusters
849  Double_t ECell1RMSSM[nMax];
850 
851  Double_t MeanPosSM[nMax];
852  Double_t MeanPosErrSM[nMax];
853  Double_t WidthSM[nMax];
854  Double_t WidthErrSM[nMax];
855  Double_t Npi0SM[nMax];
856  Double_t Npi0ErrSM[nMax];
857 
858  tree->Branch("CellMeanSM", CellMeanSM, TString::Format("CellMeanSM[%i]/D", nMax));
859  tree->Branch("CellRMSSM", CellRMSSM, TString::Format("CellRMSSM[%i]/D", nMax));
860  tree->Branch("ClusterMeanSM", ClusterMeanSM, TString::Format("ClusterMeanSM[%i]/D", nMax));
861  tree->Branch("ClusterTotSM", ClusterTotSM, TString::Format("ClusterTotSM[%i]/D", nMax));
862  tree->Branch("ClusterRMSSM", ClusterRMSSM, TString::Format("ClusterRMSSM[%i]/D", nMax));
863  tree->Branch("EtotalMeanSM", EtotalMeanSM, TString::Format("EtotalMeanSM[%i]/D", nMax));
864  tree->Branch("EtotalRMSSM", EtotalRMSSM, TString::Format("EtotalRMSSM[%i]/D", nMax));
865  tree->Branch("CellPerClusterMeanSM", CellPerClusterMeanSM, TString::Format("CellPerClusterMeanSM[%i]/D", nMax));
866  tree->Branch("CellPerClusterRMSSM", CellPerClusterRMSSM, TString::Format("CellPerClusterRMSSM[%i]/D", nMax));
867  tree->Branch("ECell1MeanSM", ECell1MeanSM, TString::Format("ECell1MeanSM[%i]/D", nMax));
868  tree->Branch("ECell1RMSSM", ECell1RMSSM, TString::Format("ECell1RMSSM[%i]/D", nMax));
869 
870  tree->Branch("MeanPosSM", MeanPosSM, TString::Format("MeanPosSM[%i]/D", nMax));
871  tree->Branch("MeanPosErrSM", MeanPosErrSM, TString::Format("MeanPosErrSM[%i]/D", nMax));
872  tree->Branch("WidthSM", WidthSM, TString::Format("WidthSM[%i]/D", nMax));
873  tree->Branch("WidthErrSM", WidthErrSM, TString::Format("WidthErrSM[%i]/D", nMax));
874  tree->Branch("Npi0SM", Npi0SM, TString::Format("Npi0SM[%i]/D", nMax));
875  tree->Branch("Npi0ErrSM", Npi0ErrSM, TString::Format("Npi0ErrSM[%i]/D", nMax));
876 
877  TF1* fitMass = new TF1("fitMass", pi0massP2, 100, 250, 6);
878  fitMass->SetParName(0, "A");
879  fitMass->SetParName(1, "m_{0}");
880  fitMass->SetParName(2, "sigma");
881  fitMass->SetParName(3, "a_{0}");
882  fitMass->SetParName(4, "a_{1}");
883  fitMass->SetParName(5, "a_{2}");
884  fitMass->SetParLimits(0, 1.e-5, 1.e5);
885  fitMass->SetParLimits(1, 0.11, 0.16);
886  fitMass->SetParLimits(2, 0.001, 0.06);
887 
888  TList* outputList;
889 
890  TH1F* fhNEvents;
891  TH1F* fhE;
892  TH1F* fhNClusters = 0x0;
893  TH1F* fhECharged;
894  TH1F* fhNCells = 0x0;
895 
896  TH1F* NCells[n];
897  TH1F* NClusters[n];
898  TH2F* NCellsPerCluster[n];
899  TH1F* E[n];
900 
901  TH2F* fhIM;
902  TH2F* fhIMDCAL;
903  TH1F* fhMggEMCAL;
904  TH1F* fhMggDCAL;
905  TH2F* IM[n];
906  TH1F* MggSM[n];
907 
908  TPMERegexp r("_\\w+");
909  TIter next(TriggersList);
910  int ret = 0;
911 
912  while(TObject *obj = next()){
913  fTrigger= TString(obj->GetName());
914  TString namefile = QAPATH + period + "_" + pass + fTrigger(r).Data() + "_" + RunId + "_data.txt";
915  ofstream QAData(namefile, ios::app); // write checks at the end
916 
917  NTClusters=0;
918  NTClustersRMS=0;
919  NCClusters=0;
920  NCClustersRMS=0;
921  NMatchClustersP=0;
922  NMatchClustersPRMS=0;
923  CellMean=0;
924  CellRMS=0;
925  ClusterMean=0;
926  ClusterRMS=0;
927  EtotalMean=0;
928  EtotalRMS=0;
929  CellPerClusterMean=0;
930  CellPerClusterRMS=0;
931 
932  // EMCal
933  Npi0EMCAL=0;
934  Npi0EMCALErr=0;
935  MeanPosEMCAL=0;
936  MeanPosEMCALErr=0;
937  WidthEMCAL=0;
938  WidthEMCALErr=0;
939  Chi2NdfPi0EMCAL=0;
940  NggEMCAL=0;
941  NggEMCALErr=0;
942  SignifEMCAL=0;
943  SignifEMCALErr=0;
944 
945  // DCal
946  Npi0DCAL=0;
947  Npi0DCALErr=0;
948  MeanPosDCAL=0;
949  MeanPosDCALErr=0;
950  WidthDCAL=0;
951  WidthDCALErr=0;
952  Chi2NdfPi0DCAL=0;
953  NggDCAL=0;
954  NggDCALErr=0;
955  SignifDCAL=0;
956  SignifDCALErr=0;
957 
958  memset (CellMeanSM, 0, sizeof (Double_t) * nMax);
959  memset (CellRMSSM, 0, sizeof (Double_t) * nMax);
960  memset (ClusterMeanSM, 0, sizeof (Double_t) * nMax);
961  memset (ClusterTotSM, 0, sizeof (Double_t) * nMax);
962  memset (ClusterRMSSM, 0, sizeof (Double_t) * nMax);
963  memset (EtotalMeanSM, 0, sizeof (Double_t) * nMax);
964  memset (EtotalRMSSM, 0, sizeof (Double_t) * nMax);
965  memset (CellPerClusterMeanSM, 0, sizeof (Double_t) * nMax);
966  memset (CellPerClusterRMSSM, 0, sizeof (Double_t) * nMax);
967  memset (ECell1MeanSM, 0, sizeof (Double_t) * nMax);
968  memset (ECell1RMSSM, 0, sizeof (Double_t) * nMax);
969 
970  memset (MeanPosSM, 0, sizeof (Double_t) * nMax);
971  memset (MeanPosErrSM, 0, sizeof (Double_t) * nMax);
972  memset (WidthSM, 0, sizeof (Double_t) * nMax);
973  memset (WidthErrSM, 0, sizeof (Double_t) * nMax);
974  memset (Npi0SM, 0, sizeof (Double_t) * nMax);
975  memset (Npi0ErrSM, 0, sizeof (Double_t) * nMax);
976 
977  TString dirname;
978  if(!fTrigger.Contains("QA")){
979  dirname = "CaloQA_";
980  }
981  dirname += fTrigger;
982 
983  Bool_t dirok = f->cd(dirname);
984  if(!dirok){
985  Error(__FUNCTION__, Form("No input directory %s", dirname.Data()));
986  tree->Fill(); ftree->cd();
987  tree->Write();
988  ret= -1;
989  continue;
990  }
991 
992  outputList = (TList*)gDirectory->Get(dirname);
993  if(!outputList){
994  Error(__FUNCTION__, Form("No input list! %s", dirname.Data()));
995  tree->Fill(); ftree->cd();
996  tree->Write();
997  ret=-2;
998  continue;;
999  }
1000  outputList->SetOwner();
1001 
1002  // number of events
1003  fhNEvents =(TH1F *)outputList->FindObject("hNEvents");
1004  if(!fhNEvents){
1005  Error(__FUNCTION__, Form("NEvent histogram not found for trigger %s", fTrigger.Data()));
1006  tree->Fill(); ftree->cd();
1007  tree->Write(); ret=-3;
1008  continue;
1009  }
1010 
1011  Nevent=fhNEvents->GetEntries();
1012 
1013  if(Nevent==0){
1014  Error(__FUNCTION__, Form("No event in trigger %s", fTrigger.Data()));
1015  tree->Fill(); ftree->cd();
1016  tree->Write();
1017  ret=-4;
1018  continue;
1019  }
1020  if(Nevent<20){
1021  Error(__FUNCTION__, Form("Less than 20 events in trigger %s", fTrigger.Data()));
1022  tree->Fill(); ftree->cd();
1023  tree->Write();
1024  ret=-5;
1025  continue;
1026  }
1027 
1028  TH1F* hE =(TH1F *)outputList->FindObject("EMCAL_hE");
1029  if(!hE || (!hE->GetEntries())){
1030  Error(__FUNCTION__, Form("hE Histogram not found or it is empty for trigger %s, EMCal was not in this run?", fTrigger.Data())); EtotalMean=0;
1031  tree->Fill(); ftree->cd();
1032  tree->Write();
1033  ret=-4;
1034  continue;
1035  }
1036 
1037  // first do clusters trending
1038  fhE = (TH1F *)outputList->FindObject(fCalorimeter+"_hE");
1039  fhECharged = (TH1F *)outputList->FindObject(fCalorimeter+"_hECharged");
1040 
1041  Double_t energy = 0. ;
1042 
1043  for(Int_t ibin = fhE->FindBin(0.6) ; ibin <fhE->FindBin(50.) ; ibin++){
1044  energy+=fhE->GetBinCenter(ibin)*fhE->GetBinContent(ibin);
1045  }
1046 
1047  if(fhE->Integral(fhE->FindBin(0.6), fhE->FindBin(50.))==0){
1048  Error(__FUNCTION__, Form("Not enough events"));
1049  tree->Fill(); ftree->cd();
1050  tree->Write();
1051  ret=-6;
1052  continue;
1053  }
1054 
1055  EtotalMean=energy/fhE->Integral(fhE->FindBin(0.6), fhE->FindBin(50.)) ;
1056  EtotalRMS=fhE->GetMeanError();
1057 
1058  TString nameNCell = Form("%s_hNCells_Mod", fCalorimeter.Data());
1059  TString nameNCluster = Form("%s_hNClusters_Mod", fCalorimeter.Data());
1060  TString nameE = Form("%s_hE_Mod", fCalorimeter.Data());
1061  TH2F* hNCellsMod= (TH2F*)outputList->FindObject(nameNCell.Data());
1062  TH2F* hNClusterMod=(TH2F*)outputList->FindObject(nameNCluster.Data());
1063  TH2F* hEMod=(TH2F*)outputList->FindObject(nameE.Data());
1064 
1065  if(!hNCellsMod || !hNClusterMod || !hEMod){
1066  Error(__FUNCTION__, "A requiered histogram was not found (the imput QAresult.root might be too old)!");
1067  tree->Fill(); ftree->cd();
1068  tree->Write();
1069  ret=-7;
1070  continue;
1071  }
1072 
1073  //::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
1074  // EM addition (Nov. 2016) to observe the PAR Time issue
1075  TString trig_title = "";
1076  if(fTrigger.Contains("EMC"))
1077  trig_title += "EMC";
1078  else
1079  trig_title += "MB";
1080 
1081  // BC/4 = 0
1082  TCanvas* c_TimeSM_BC0 = new TCanvas("TimeSM_BC0", "Cell Time (for BC/4=0) for each SM", 1200, 1200);
1083  c_TimeSM_BC0->SetFillColor(0);
1084  c_TimeSM_BC0->SetBorderSize(0);
1085  c_TimeSM_BC0->SetFrameBorderMode(0);
1086  Bool_t mod3_BC0=0;
1087  if(n%3)
1088  mod3_BC0=1;
1089  c_TimeSM_BC0->Divide(3, n/3+mod3_BC0);
1090 
1091  TH2F* hTimeSM_BC0 = (TH2F*)outputList->FindObject("EMCAL_hTimePerSM_BC0");
1092  TH1F* hTimeSM_BC0_proj[n];
1093  TString TimeSM_BC0_projName;
1094 
1095  if(hTimeSM_BC0){
1096  cout << "---------------------" << endl;
1097  cout << " BC/4=0 | Run " << RunId << endl;
1098  for(Int_t ism = 0 ; ism < n ; ism++){
1099  c_TimeSM_BC0->cd(ism+1);
1100  gPad->SetTopMargin(0.1);
1101  gPad->SetBottomMargin(0.15);
1102  // gPad->SetLogy();
1103 
1104  TimeSM_BC0_projName = Form("SM_%d", ism);
1105  hTimeSM_BC0_proj[ism] = (TH1F *)hTimeSM_BC0->ProjectionX(TimeSM_BC0_projName.Data(), ism+1, ism+1, "") ;
1106  hTimeSM_BC0_proj[ism]->Draw("hist");
1107 
1108  hTimeSM_BC0_proj[ism]->SetTitle(Form("#it{t}_{cell} for super module %i (BC/4=0) %s", ism, trig_title.Data()));
1109  hTimeSM_BC0_proj[ism]->SetTitleSize(0.1);
1110  hTimeSM_BC0_proj[ism]->SetXTitle("#it{t}_{cell} (ns)");
1111  hTimeSM_BC0_proj[ism]->SetYTitle("Nb of entries");
1112  hTimeSM_BC0_proj[ism]->SetStats(kTRUE);
1113 
1114  hTimeSM_BC0_proj[ism]->GetXaxis()->SetRangeUser(400, 900);
1115  hTimeSM_BC0_proj[ism]->GetXaxis()->SetLabelSize(0.05);
1116  hTimeSM_BC0_proj[ism]->GetXaxis()->SetLabelOffset(0.015);
1117  hTimeSM_BC0_proj[ism]->GetXaxis()->SetTitleSize(0.075);
1118  hTimeSM_BC0_proj[ism]->GetXaxis()->SetTitleOffset(0.9);
1119 
1120  hTimeSM_BC0_proj[ism]->GetYaxis()->SetLabelSize(0.05);
1121  hTimeSM_BC0_proj[ism]->GetYaxis()->SetTitleSize(0.06);
1122  hTimeSM_BC0_proj[ism]->GetYaxis()->SetTitleOffset(0.6);
1123  }
1124 
1125  TString outfilename_BC0 = QAPATH + "TimeSM_BC0" + fTrigger(r) + ".pdf";
1126  TString outfilename2_BC0 = QAPATH + "TimeSM_BC0" + fTrigger(r) + ".png";
1127  if(SavePlots==2)
1128  c_TimeSM_BC0->SaveAs(outfilename_BC0);
1129  if(SavePlots)
1130  c_TimeSM_BC0->SaveAs(outfilename2_BC0);
1131  }
1132 
1133  // BC/4 = 1
1134  TCanvas* c_TimeSM_BC1 = new TCanvas("TimeSM_BC1", "Cell Time (for BC/4=1) for each SM", 1200, 1200);
1135  c_TimeSM_BC1->SetFillColor(0);
1136  c_TimeSM_BC1->SetBorderSize(0);
1137  c_TimeSM_BC1->SetFrameBorderMode(0);
1138  Bool_t mod3_BC1=0;
1139  if(n%3)
1140  mod3_BC1=1;
1141  c_TimeSM_BC1->Divide(3, n/3+mod3_BC1);
1142 
1143  TH2F* hTimeSM_BC1 = (TH2F*)outputList->FindObject("EMCAL_hTimePerSM_BC1");
1144  TH1F* hTimeSM_BC1_proj[n];
1145  TString TimeSM_BC1_projName;
1146 
1147  if(hTimeSM_BC1){
1148  cout << "---------------------" << endl;
1149  cout << " BC/4=1 | Run " << RunId << endl;
1150  for(Int_t ism = 0 ; ism < n ; ism++){
1151  c_TimeSM_BC1->cd(ism+1);
1152  gPad->SetTopMargin(0.1);
1153  gPad->SetBottomMargin(0.15);
1154  // gPad->SetLogy();
1155 
1156  TimeSM_BC1_projName = Form("SM_%d", ism);
1157  hTimeSM_BC1_proj[ism] = (TH1F *)hTimeSM_BC1->ProjectionX(TimeSM_BC1_projName.Data(), ism+1, ism+1, "") ;
1158  hTimeSM_BC1_proj[ism]->Draw("hist");
1159 
1160  hTimeSM_BC1_proj[ism]->SetTitle(Form("#it{t}_{cell} for super module %i (BC/4=1) %s", ism, trig_title.Data()));
1161  hTimeSM_BC1_proj[ism]->SetTitleSize(0.1);
1162  hTimeSM_BC1_proj[ism]->SetXTitle("#it{t}_{cell} (ns)");
1163  hTimeSM_BC1_proj[ism]->SetYTitle("Nb of entries");
1164  hTimeSM_BC1_proj[ism]->SetStats(kTRUE);
1165 
1166  hTimeSM_BC1_proj[ism]->GetXaxis()->SetRangeUser(400, 900);
1167  hTimeSM_BC1_proj[ism]->GetXaxis()->SetLabelSize(0.05);
1168  hTimeSM_BC1_proj[ism]->GetXaxis()->SetLabelOffset(0.015);
1169  hTimeSM_BC1_proj[ism]->GetXaxis()->SetTitleSize(0.075);
1170  hTimeSM_BC1_proj[ism]->GetXaxis()->SetTitleOffset(0.9);
1171 
1172  hTimeSM_BC1_proj[ism]->GetYaxis()->SetLabelSize(0.05);
1173  hTimeSM_BC1_proj[ism]->GetYaxis()->SetTitleSize(0.06);
1174  hTimeSM_BC1_proj[ism]->GetYaxis()->SetTitleOffset(0.6);
1175  }
1176 
1177  TString outfilename_BC1 = QAPATH + "TimeSM_BC1" + fTrigger(r) + ".pdf";
1178  TString outfilename2_BC1 = QAPATH + "TimeSM_BC1" + fTrigger(r) + ".png";
1179  if(SavePlots==2)
1180  c_TimeSM_BC1->SaveAs(outfilename_BC1);
1181  if(SavePlots)
1182  c_TimeSM_BC1->SaveAs(outfilename2_BC1);
1183  }
1184 
1185  // BC/4 = 2
1186  TCanvas* c_TimeSM_BC2 = new TCanvas("TimeSM_BC2", "Cell Time (for BC/4=2) for each SM", 1200, 1200);
1187  c_TimeSM_BC2->SetFillColor(0);
1188  c_TimeSM_BC2->SetBorderSize(0);
1189  c_TimeSM_BC2->SetFrameBorderMode(0);
1190  Bool_t mod3_BC2=0;
1191  if(n%3)
1192  mod3_BC2=1;
1193  c_TimeSM_BC2->Divide(3, n/3+mod3_BC2);
1194 
1195  TH2F* hTimeSM_BC2 = (TH2F*)outputList->FindObject("EMCAL_hTimePerSM_BC2");
1196  TH1F* hTimeSM_BC2_proj[n];
1197  TString TimeSM_BC2_projName;
1198 
1199  if(hTimeSM_BC2){
1200  cout << "---------------------" << endl;
1201  cout << " BC/4=2 | Run " << RunId << endl;
1202  for(Int_t ism = 0 ; ism < n ; ism++){
1203  c_TimeSM_BC2->cd(ism+1);
1204  gPad->SetTopMargin(0.1);
1205  gPad->SetBottomMargin(0.15);
1206  // gPad->SetLogy();
1207 
1208  TimeSM_BC2_projName = Form("SM_%d", ism);
1209  hTimeSM_BC2_proj[ism] = (TH1F *)hTimeSM_BC2->ProjectionX(TimeSM_BC2_projName.Data(), ism+1, ism+1, "") ;
1210  hTimeSM_BC2_proj[ism]->Draw("hist");
1211 
1212  hTimeSM_BC2_proj[ism]->SetTitle(Form("#it{t}_{cell} for super module %i (BC/4=2) %s", ism, trig_title.Data()));
1213  hTimeSM_BC2_proj[ism]->SetTitleSize(0.1);
1214  hTimeSM_BC2_proj[ism]->SetXTitle("#it{t}_{cell} (ns)");
1215  hTimeSM_BC2_proj[ism]->SetYTitle("Nb of entries");
1216  hTimeSM_BC2_proj[ism]->SetStats(kTRUE);
1217 
1218  hTimeSM_BC2_proj[ism]->GetXaxis()->SetRangeUser(400, 900);
1219  hTimeSM_BC2_proj[ism]->GetXaxis()->SetLabelSize(0.05);
1220  hTimeSM_BC2_proj[ism]->GetXaxis()->SetLabelOffset(0.015);
1221  hTimeSM_BC2_proj[ism]->GetXaxis()->SetTitleSize(0.075);
1222  hTimeSM_BC2_proj[ism]->GetXaxis()->SetTitleOffset(0.9);
1223 
1224  hTimeSM_BC2_proj[ism]->GetYaxis()->SetLabelSize(0.05);
1225  hTimeSM_BC2_proj[ism]->GetYaxis()->SetTitleSize(0.06);
1226  hTimeSM_BC2_proj[ism]->GetYaxis()->SetTitleOffset(0.6);
1227  }
1228 
1229  TString outfilename_BC2 = QAPATH + "TimeSM_BC2" + fTrigger(r) + ".pdf";
1230  TString outfilename2_BC2 = QAPATH + "TimeSM_BC2" + fTrigger(r) + ".png";
1231  if(SavePlots==2)
1232  c_TimeSM_BC2->SaveAs(outfilename_BC2);
1233  if(SavePlots)
1234  c_TimeSM_BC2->SaveAs(outfilename2_BC2);
1235  }
1236 
1237  // BC/4 = 3
1238  TCanvas* c_TimeSM_BC3 = new TCanvas("TimeSM_BC3", "Cell Time (for BC/4=3) for each SM", 1200, 1200);
1239  c_TimeSM_BC3->SetFillColor(0);
1240  c_TimeSM_BC3->SetBorderSize(0);
1241  c_TimeSM_BC3->SetFrameBorderMode(0);
1242  Bool_t mod3_BC3=0;
1243  if(n%3)
1244  mod3_BC3=1;
1245  c_TimeSM_BC3->Divide(3, n/3+mod3_BC3);
1246 
1247  TH2F* hTimeSM_BC3 = (TH2F*)outputList->FindObject("EMCAL_hTimePerSM_BC3");
1248  TH1F* hTimeSM_BC3_proj[n];
1249  TString TimeSM_BC3_projName;
1250 
1251  if(hTimeSM_BC3){
1252  cout << "---------------------" << endl;
1253  cout << " BC/4=3 | Run " << RunId << endl;
1254  for(Int_t ism = 0 ; ism < n ; ism++){
1255  c_TimeSM_BC3->cd(ism+1);
1256  gPad->SetTopMargin(0.1);
1257  gPad->SetBottomMargin(0.15);
1258  // gPad->SetLogy();
1259 
1260  TimeSM_BC3_projName = Form("SM_%d", ism);
1261  hTimeSM_BC3_proj[ism] = (TH1F *)hTimeSM_BC3->ProjectionX(TimeSM_BC3_projName.Data(), ism+1, ism+1, "") ;
1262  hTimeSM_BC3_proj[ism]->Draw("hist");
1263 
1264  hTimeSM_BC3_proj[ism]->SetTitle(Form("#it{t}_{cell} for super module %i (BC/4=3) %s", ism, trig_title.Data()));
1265  hTimeSM_BC3_proj[ism]->SetTitleSize(0.1);
1266  hTimeSM_BC3_proj[ism]->SetXTitle("#it{t}_{cell} (ns)");
1267  hTimeSM_BC3_proj[ism]->SetYTitle("Nb of entries");
1268  hTimeSM_BC3_proj[ism]->SetStats(kTRUE);
1269 
1270  hTimeSM_BC3_proj[ism]->GetXaxis()->SetRangeUser(400, 900);
1271  hTimeSM_BC3_proj[ism]->GetXaxis()->SetLabelSize(0.05);
1272  hTimeSM_BC3_proj[ism]->GetXaxis()->SetLabelOffset(0.015);
1273  hTimeSM_BC3_proj[ism]->GetXaxis()->SetTitleSize(0.075);
1274  hTimeSM_BC3_proj[ism]->GetXaxis()->SetTitleOffset(0.9);
1275 
1276  hTimeSM_BC3_proj[ism]->GetYaxis()->SetLabelSize(0.05);
1277  hTimeSM_BC3_proj[ism]->GetYaxis()->SetTitleSize(0.06);
1278  hTimeSM_BC3_proj[ism]->GetYaxis()->SetTitleOffset(0.6);
1279  }
1280 
1281  TString outfilename_BC3 = QAPATH + "TimeSM_BC3" + fTrigger(r) + ".pdf";
1282  TString outfilename2_BC3 = QAPATH + "TimeSM_BC3" + fTrigger(r) + ".png";
1283  if(SavePlots==2)
1284  c_TimeSM_BC3->SaveAs(outfilename_BC3);
1285  if(SavePlots)
1286  c_TimeSM_BC3->SaveAs(outfilename2_BC3);
1287  }
1288 
1289  //::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
1290 
1291  TCanvas* c1 = new TCanvas("Pi0InvMassSM", "Pi0 Invariant Mass for each SM", 600, 600);
1292  gStyle->SetOptStat(1); // MG modif
1293  gStyle->SetOptFit(1); // MG modif
1294  c1->SetFillColor(0);
1295  c1->SetBorderSize(0);
1296  c1->SetFrameBorderMode(0);
1297  // c1->SetOptStat(1); // MG modif
1298 
1299  Bool_t mod3=0;
1300  if(n%3)
1301  mod3=1;
1302  c1->Divide(3, n/3+mod3);
1303 
1304  // per sm trending
1305  TString nameNCellPerCluster;
1306  for(Int_t ism = 0 ; ism < n ; ism++){
1307  cout << "#########################"<< endl;
1308  cout << " Super Module " << ism << " Run " << RunId << endl;
1309 
1310  // first do clusters trending
1311  nameNCellPerCluster = Form("%s_hNCellsPerCluster_Mod%d", fCalorimeter.Data(), ism);
1312  NCellsPerCluster[ism] = (TH2F*)outputList->FindObject(nameNCellPerCluster.Data());
1313 
1314  if(! (TH2F*)outputList->FindObject(nameNCellPerCluster.Data()) ){
1315  Error(__FUNCTION__, Form("NCellsPerCluster histogram not found for super module %i", ism));
1316  ret=-8;
1317  continue;
1318  }
1319 
1320  NCellsPerCluster[ism] = (TH2F*)outputList->FindObject(nameNCellPerCluster.Data());
1321  NCells[ism] = (TH1F*)hNCellsMod->ProjectionX(Form("NCells%d", ism), ism+1, ism+1, "");
1322  NClusters[ism] = (TH1F*)hNClusterMod->ProjectionX(Form("NClusters%d", ism), ism+1, ism+1, "");
1323  E[ism] = (TH1F*)hEMod->ProjectionX(Form("E%d", ism), ism+1, ism+1, "");
1324 
1325  CellMeanSM[ism]=NCells[ism]->GetMean();
1326  CellRMSSM[ism]=NCells[ism]->GetMeanError();
1327 
1328  Int_t binmax;
1329  binmax = NClusters[ism]->GetNbinsX();
1330  for(Int_t ibin = 1 ;ibin <=binmax ;ibin++){
1331  ClusterTotSM[ism]+=(NClusters[ism]->GetBinCenter(ibin))*(NClusters[ism]->GetBinContent(ibin));
1332  }
1333  if(NClusters[ism]->GetEntries()>0)
1334  ClusterTotSM[ism]= ClusterTotSM[ism]/NClusters[ism]->GetEntries();
1335 
1336  // ClusterTotSM[ism]=NClusters[ism]->GetMean();
1337  ClusterMeanSM[ism]=NClusters[ism]->GetMean();
1338  ClusterRMSSM[ism]=NClusters[ism]->GetMeanError();
1339  CellPerClusterMeanSM[ism]=NCellsPerCluster[ism]->GetMean(2);
1340  CellPerClusterRMSSM[ism]=NCellsPerCluster[ism]->GetMeanError(2);
1341 
1342  ECell1MeanSM[ism] =NCellsPerCluster[ism]->ProjectionX("", 2, 50, "")->Integral(5, 50)/(Nevent);
1343  ECell1RMSSM[ism] =NCellsPerCluster[ism]->ProjectionX("", 2, 50, "")->GetMeanError();
1344 
1345  Double_t energySM = 0. ;
1346  for(Int_t ibin = E[ism]->FindBin(0.6) ; ibin <E[ism]->FindBin(50.) ; ibin++){
1347  energySM+=E[ism]->GetBinCenter(ibin)*(E[ism]->GetBinContent(ibin));
1348  }
1349 
1350  if(E[ism]->Integral(E[ism]->FindBin(0.6), E[ism]->FindBin(50.))==0){
1351  Error(__FUNCTION__, Form("Energy: Not enough events/SM"));
1352  continue;
1353  }
1354 
1355  EtotalMeanSM[ism]=energySM/(E[ism]->Integral(E[ism]->FindBin(0.6), E[ism]->FindBin(50.)));
1356  EtotalRMSSM[ism]=E[ism]->GetMeanError();
1357 
1358  if(ism==0){
1359  fhNCells = (TH1F*)NCells[ism]->Clone("NCells");
1360  fhNClusters = (TH1F*)NClusters[ism]->Clone("NClusters");
1361  }
1362  else{
1363  fhNCells->Add(NCells[ism], 1);
1364  fhNClusters->Add(NClusters[ism], 1);
1365  }
1366 
1367  // Pi0
1368  c1->cd(ism+1);
1369  TString namePair = Form("%s_hIM_Mod%d", fCalorimeter.Data(), ism);
1370  IM[ism] = (TH2F*)outputList->FindObject(namePair.Data());
1371  IM[ism]->Sumw2();
1372 
1373  TString projname = Form("SM_%d", ism);
1374  MggSM[ism] = (TH1F *)IM[ism]->ProjectionY(projname.Data(), 2, 20, "") ; // MG modif
1375 
1376  if(MggSM[ism]->GetEntries()>100){
1377  fitMass->SetParameter(0, MggSM[ism]->GetBinContent(MggSM[ism]->GetMaximumBin()));
1378  fitMass->SetParameter(1, mPDG);
1379  fitMass->SetParameter(2, 15.);
1380  fitMass->SetParameter(3, 0.);
1381  fitMass->SetParameter(4, MggSM[ism]->GetBinContent(MggSM[ism]->FindBin(0.11)));
1382  fitMass->SetParameter(5, MggSM[ism]->GetBinContent(MggSM[ism]->FindBin(0.20)));
1383 
1384  if(MggSM[ism]->GetEntries()<1000)
1385  MggSM[ism]->Rebin(5);
1386  else
1387  MggSM[ism]->Rebin();
1388 
1389  MggSM[ism]->Fit("fitMass", "WL R +", "", 0.1, 0.30);
1390  MggSM[ism]->SetTitle(Form("Pi0 Mass for super module %i", ism));
1391  MggSM[ism]->SetTitleSize(0.1);
1392  MggSM[ism]->SetXTitle("Pi0 Mass");
1393  MggSM[ism]->SetYTitle("Nb of entries");
1394  MggSM[ism]->GetXaxis()->SetLabelSize(0.05);
1395  MggSM[ism]->GetXaxis()->SetTitleSize(0.07);
1396  MggSM[ism]->GetXaxis()->SetTitleOffset(0.68);
1397  MggSM[ism]->GetYaxis()->SetLabelSize(0.05);
1398  MggSM[ism]->GetYaxis()->SetTitleSize(0.06);
1399  MggSM[ism]->GetYaxis()->SetTitleOffset(0.78);
1400 
1401  MeanPosSM[ism] = MggSM[ism]->GetFunction("fitMass")->GetParameter(1)*1000;
1402  MeanPosErrSM[ism] = MggSM[ism]->GetFunction("fitMass")->GetParError(1)*1000;
1403  WidthSM[ism] = MggSM[ism]->GetFunction("fitMass")->GetParameter(2)*1000;
1404  WidthErrSM[ism] = MggSM[ism]->GetFunction("fitMass")->GetParError(2)*1000;
1405  Npi0SM[ism] = MggSM[ism]->GetFunction("fitMass")->GetParameter(0)*(MggSM[ism]->GetFunction("fitMass")->GetParameter(2))*TMath::Sqrt(2*TMath::Pi())/(Nevent*MggSM[ism]->GetBinWidth(1));
1406  Npi0ErrSM[ism] = TMath::Sqrt((MggSM[ism]->GetFunction("fitMass")->GetParError(0)/MggSM[ism]->GetFunction("fitMass")->GetParameter(0))*(MggSM[ism]->GetFunction("fitMass")->GetParError(0)/MggSM[ism]->GetFunction("fitMass")->GetParameter(0))
1407  +(MggSM[ism]->GetFunction("fitMass")->GetParError(2)/MggSM[ism]->GetFunction("fitMass")->GetParameter(2))*(MggSM[ism]->GetFunction("fitMass")->GetParError(2)/MggSM[ism]->GetFunction("fitMass")->GetParameter(2)));
1408  Npi0ErrSM[ism] = 0.;
1409 
1410  } // end if enough events for Pi0 fit and trending
1411  else
1412  Info(__FUNCTION__, Form("Not enough events for Pi0 fit and trending for super module %i", ism));
1413  } // per SM loop
1414 
1415  // Now Pi0 global trending in EMCAL
1416  TCanvas* c2 = new TCanvas("Pi0InvMassEMCAL", "Pi0 Invariant Mass in EMCal", 600, 600);
1417  c2->SetFillColor(0);
1418  c2->SetBorderSize(0);
1419  c2->SetFrameBorderMode(0);
1420  // c2->SetOptStat(1); // MG modif
1421 
1422  fhIM = (TH2F *)outputList->FindObject(fCalorimeter+"_hIM");
1423  fhIM->Sumw2();
1424  fhMggEMCAL = (TH1F *)fhIM->ProjectionY("MggEMCAL", 2, 20, "") ; // to modify projection range
1425 
1426  if(fhMggEMCAL->GetEntries()==0){
1427  Error(__FUNCTION__, "The Pi0 histogram in EMCal is empty !"); tree->Fill();
1428  ret=-8;
1429  continue;
1430  }
1431 
1432  fitMass->SetParameter(0, 4500);
1433  fitMass->SetParameter(1, mPDG);
1434  fitMass->SetParameter(2, 0.01);
1435  fitMass->SetParameter(3, 0.);
1436  fitMass->SetParameter(4, fhMggEMCAL->GetBinContent(fhMggEMCAL->FindBin(0.11)));
1437  fitMass->SetParameter(5, fhMggEMCAL->GetBinContent(fhMggEMCAL->FindBin(0.20)));
1438 
1439  if(fhMggEMCAL->GetEntries()<5000)
1440  fhMggEMCAL->Rebin(5); // MG modif
1441  else
1442  fhMggEMCAL->Rebin();
1443 
1444  fhMggEMCAL->Fit("fitMass", "L R +", "", 0.05, 0.20);
1445 
1446  fhMggEMCAL->SetTitle("Pi0 Mass in EMCal");
1447  fhMggEMCAL->SetTitleSize(0.1);
1448  fhMggEMCAL->SetXTitle("Pi0 Mass");
1449  fhMggEMCAL->SetYTitle("Nb of entries");
1450  fhMggEMCAL->GetXaxis()->SetLabelSize(0.03);
1451  fhMggEMCAL->GetXaxis()->SetTitleSize(0.03);
1452  fhMggEMCAL->GetXaxis()->SetTitleOffset(1.3);
1453  fhMggEMCAL->GetYaxis()->SetLabelSize(0.03);
1454  fhMggEMCAL->GetYaxis()->SetTitleSize(0.03);
1455  fhMggEMCAL->GetYaxis()->SetTitleOffset(1.3);
1456 
1457  MeanPosEMCAL = fhMggEMCAL->GetFunction("fitMass")->GetParameter(1)*1000;
1458  MeanPosEMCALErr = fhMggEMCAL->GetFunction("fitMass")->GetParError(1)*1000;
1459  WidthEMCAL = fhMggEMCAL->GetFunction("fitMass")->GetParameter(2)*1000;
1460  WidthEMCALErr = fhMggEMCAL->GetFunction("fitMass")->GetParError(2)*1000;
1461  Chi2NdfPi0EMCAL = fhMggEMCAL->GetFunction("fitMass")->GetChisquare()/fhMggEMCAL->GetFunction("fitMass")->GetNDF();
1462  Npi0EMCAL = fhMggEMCAL->GetFunction("fitMass")->GetParameter(0)*fhMggEMCAL->GetFunction("fitMass")->GetParameter(2)*TMath::Sqrt(2*TMath::Pi())/(Nevent*fhMggEMCAL->GetBinWidth(10));
1463  Npi0EMCALErr = TMath::Sqrt((fhMggEMCAL->GetFunction("fitMass")->GetParError(0)/fhMggEMCAL->GetFunction("fitMass")->GetParameter(0))*(fhMggEMCAL->GetFunction("fitMass")->GetParError(0)/fhMggEMCAL->GetFunction("fitMass")->GetParameter(0))+(WidthEMCALErr/WidthEMCAL)*(WidthEMCALErr/WidthEMCAL));
1464  Npi0EMCALErr = 0.;
1465  NggEMCAL = fhMggEMCAL->GetFunction("fitMass")->Integral(0.11, 0.16)/(Nevent*fhMggEMCAL->GetBinWidth(10));
1466  NggEMCALErr = fhMggEMCAL->GetFunction("fitMass")->IntegralError(0.11, 0.16)/(fhMggEMCAL->Integral()*fhMggEMCAL->GetBinWidth(10));
1467  SignifEMCAL = Npi0EMCAL/NggEMCAL;
1468  SignifEMCALErr = TMath::Sqrt((Npi0EMCALErr/Npi0EMCAL*(Npi0EMCALErr/Npi0EMCAL)+(NggEMCALErr/NggEMCAL*(NggEMCALErr/NggEMCAL))));
1469  SignifEMCALErr = SignifEMCAL*SignifEMCALErr;
1470 
1471  TCanvas* c3 = new TCanvas("Pi0InvMassDCAL", "Pi0 Invariant Mass in DCal", 600, 600);
1472  c3->SetFillColor(0);
1473  c3->SetBorderSize(0);
1474  c3->SetFrameBorderMode(0);
1475 
1476  fhIMDCAL = (TH2F *)outputList->FindObject(fCalorimeter+"_hIMDCAL");
1477  if(fhIMDCAL){
1478  // Now Pi0 global trending in DCal
1479  // TCanvas* c3 = new TCanvas("Pi0InvMassDCAL", "Pi0 Invariant Mass in DCal", 600, 600);
1480  // c3->SetFillColor(0);
1481  // c3->SetBorderSize(0);
1482  // c3->SetFrameBorderMode(0);
1483  // c2->SetOptStat(1); // MG modif
1484 
1485  // fhIMDCAL = (TH2F *)outputList->FindObject(fCalorimeter+"_hIMDCAL");
1486  // if(!fhIMDCAL) continue;
1487  fhIMDCAL->Sumw2();
1488 
1489  fhMggDCAL = (TH1F *)fhIMDCAL->ProjectionY("MggDCAL", 2, 20, "") ; // to modify projection range
1490  if(fhMggDCAL->GetEntries()==0){
1491  Error(__FUNCTION__, "The Pi0 histogram in DCal is empty !"); tree->Fill();
1492  ret=-8;
1493  continue;
1494  }
1495 
1496  fitMass->SetParameter(0, 4500);
1497  fitMass->SetParameter(1, mPDG);
1498  fitMass->SetParameter(2, 0.01);
1499  fitMass->SetParameter(3, 0.);
1500  fitMass->SetParameter(4, fhMggDCAL->GetBinContent(fhMggDCAL->FindBin(0.11)));
1501  fitMass->SetParameter(5, fhMggDCAL->GetBinContent(fhMggDCAL->FindBin(0.20)));
1502 
1503  if(fhMggDCAL->GetEntries()<5000)
1504  fhMggDCAL->Rebin(5); // MG modif
1505  else
1506  fhMggDCAL->Rebin();
1507 
1508  fhMggDCAL->Fit("fitMass", "L R +", "", 0.05, 0.20);
1509 
1510  fhMggDCAL->SetTitle("Pi0 Mass in DCal");
1511  fhMggDCAL->SetTitleSize(0.1);
1512  fhMggDCAL->SetXTitle("Pi0 Mass");
1513  fhMggDCAL->SetYTitle("Nb of entries");
1514  fhMggDCAL->GetXaxis()->SetLabelSize(0.03);
1515  fhMggDCAL->GetXaxis()->SetTitleSize(0.03);
1516  fhMggDCAL->GetXaxis()->SetTitleOffset(1.3);
1517  fhMggDCAL->GetYaxis()->SetLabelSize(0.03);
1518  fhMggDCAL->GetYaxis()->SetTitleSize(0.03);
1519  fhMggDCAL->GetYaxis()->SetTitleOffset(1.3);
1520 
1521  MeanPosDCAL = fhMggDCAL->GetFunction("fitMass")->GetParameter(1)*1000;
1522  MeanPosDCALErr = fhMggDCAL->GetFunction("fitMass")->GetParError(1)*1000;
1523  WidthDCAL = fhMggDCAL->GetFunction("fitMass")->GetParameter(2)*1000;
1524  WidthDCALErr = fhMggDCAL->GetFunction("fitMass")->GetParError(2)*1000;
1525  Chi2NdfPi0DCAL = fhMggDCAL->GetFunction("fitMass")->GetChisquare()/fhMggDCAL->GetFunction("fitMass")->GetNDF();
1526  Npi0DCAL = fhMggDCAL->GetFunction("fitMass")->GetParameter(0)*fhMggDCAL->GetFunction("fitMass")->GetParameter(2)*TMath::Sqrt(2*TMath::Pi())/(Nevent*fhMggDCAL->GetBinWidth(10));
1527  Npi0DCALErr = TMath::Sqrt((fhMggDCAL->GetFunction("fitMass")->GetParError(0)/fhMggDCAL->GetFunction("fitMass")->GetParameter(0))*(fhMggDCAL->GetFunction("fitMass")->GetParError(0)/fhMggDCAL->GetFunction("fitMass")->GetParameter(0))+(WidthDCALErr/WidthDCAL)*(WidthDCALErr/WidthDCAL));
1528  Npi0DCALErr = 0.;
1529  NggDCAL = fhMggDCAL->GetFunction("fitMass")->Integral(0.11, 0.16)/(Nevent*fhMggDCAL->GetBinWidth(10));
1530  NggDCALErr = fhMggDCAL->GetFunction("fitMass")->IntegralError(0.11, 0.16)/(fhMggDCAL->Integral()*fhMggDCAL->GetBinWidth(10));
1531  SignifDCAL = Npi0DCAL/NggDCAL;
1532  SignifDCALErr = TMath::Sqrt((Npi0DCALErr/Npi0DCAL*(Npi0DCALErr/Npi0DCAL)+(NggDCALErr/NggDCAL*(NggDCALErr/NggDCAL))));
1533  SignifDCALErr = SignifDCAL*SignifDCALErr;
1534  }
1535 
1536  cout<<"******************"<<endl;
1537  // end of global trending
1538 
1539  // cout<<"********************************** Total number of clusters **********************************"<<endl;
1540  NTClusters=fhE->IntegralAndError(fhE->GetXaxis()->FindBin(0.), fhE->GetXaxis()->FindBin(200.), NTClustersRMS);
1541  NCClusters=fhECharged->IntegralAndError(fhECharged->GetXaxis()->FindBin(0.), fhECharged->GetXaxis()->FindBin(50.), NCClustersRMS);
1542 
1543  if(NTClusters!=0)
1544  NMatchClustersP = NCClusters/NTClusters;
1545  else
1546  NMatchClustersP=0;
1547  // cout<<" here is searched value "<<NMatchClustersP<<" "<<endl;
1548 
1549  if(NCClusters!=0)
1550  NMatchClustersPRMS = NMatchClustersP* TMath::Sqrt(TMath::Power(NCClustersRMS/NCClusters, 2)+TMath::Power(NTClustersRMS/NTClusters, 2));
1551  else
1552  NMatchClustersPRMS=1;
1553 
1554  ClusterMean=fhNClusters->GetMean();
1555  ClusterRMS=fhNClusters->GetMeanError();
1556  CellMean=fhNCells->GetMean();
1557  CellRMS=fhNCells->GetMeanError();
1558 
1559  tree->Fill();
1560 
1561  TString outfilename = QAPATH + "Pi0InvMassEMCAL" + fTrigger(r) + ".pdf";
1562  TString outfilename2 = QAPATH + "Pi0InvMassEMCAL" + fTrigger(r) + ".png";
1563  if(SavePlots==2)
1564  c2->SaveAs(outfilename);
1565  if(SavePlots)
1566  c2->SaveAs(outfilename2);
1567 
1568  outfilename = QAPATH + "Pi0InvMassSM" + fTrigger(r) + ".pdf";
1569  outfilename2 = QAPATH + "Pi0InvMassSM" + fTrigger(r) + ".png";
1570  if(SavePlots==2)
1571  c1->SaveAs(outfilename);
1572  if(SavePlots)
1573  c1->SaveAs(outfilename2);
1574  // cout<<"********************************** "<<n<<" **********************************"<<endl; // n = number of SMs
1575 
1576  if(fhIMDCAL){
1577  TString outfilename4 = QAPATH + "Pi0InvMassDCAL" + fTrigger(r) + ".pdf";
1578  TString outfilename5= QAPATH + "Pi0InvMassDCAL" + fTrigger(r) + ".png";
1579  if(SavePlots==2)
1580  c3->SaveAs(outfilename4);
1581  if(SavePlots)
1582  c3->SaveAs(outfilename5);
1583  }
1584 
1585  fout->cd();
1586  fout->Cd(Form("%s/%s/%ld/%s/%s", period.Data(), pass.Data(), RunId, "RunLevelQA", fTrigger.Data()));
1587 
1588  // Cell time distribution for each BC and SM
1589  if(hTimeSM_BC0) c_TimeSM_BC0->Write();
1590  if(hTimeSM_BC1) c_TimeSM_BC1->Write();
1591  if(hTimeSM_BC2) c_TimeSM_BC2->Write();
1592  if(hTimeSM_BC3) c_TimeSM_BC3->Write();
1593 
1594  c2->Write();
1595  c1->Write();
1596  if(fhIMDCAL)
1597  c3->Write();
1598 
1599  delete c1;
1600  delete c2;
1601  if(fhIMDCAL) delete c3;
1602  if(hTimeSM_BC0) delete c_TimeSM_BC0;
1603  if(hTimeSM_BC1) delete c_TimeSM_BC1;
1604  if(hTimeSM_BC2) delete c_TimeSM_BC2;
1605  if(hTimeSM_BC3) delete c_TimeSM_BC3;
1606  if(outputList) outputList->Delete();
1607 
1608  QAData << RunId << " " << Nevent
1609  << "\n";
1610 
1611  QAData.close();
1612  }
1613 
1614  ftree->cd();
1615  tree->Write();
1616  ftree->Close();
1617 
1618  return ret;
1619 
1620 }
1621 
1622 //-------------------------------------------------------------------------
1623 TH2F* FormatRunHisto(TH2F* aHisto, const char* title, const char* YTitle){
1624 
1625  if(!aHisto){
1626  Error(__FUNCTION__, Form("The histogram with title \"%s\" was not found!", title));
1627  return new TH2F();
1628  }
1629 
1630  aHisto->SetStats(kFALSE);
1631  aHisto->SetTitle(title);
1632  aHisto->SetStats(kFALSE);
1633  aHisto->SetYTitle(YTitle);
1634  aHisto->GetYaxis()->SetTitleOffset(1.2);
1635  aHisto->GetYaxis()->SetLabelSize(0.03);
1636  aHisto->GetZaxis()->SetLabelSize(0.02);
1637 
1638  return aHisto;
1639 
1640 }
1641 
1642 //--------------------------------------------------------------------------------------------------------------
1643 TH2F* HistoPerMod(TH2F* hTmpPerMod, const char* title){
1644 
1645  if(!hTmpPerMod){
1646  Error(__FUNCTION__, Form("The histogram with title \"%s\" was not found!", title));
1647  return new TH2F();
1648  }
1649 
1650  hTmpPerMod->SetStats(kFALSE);
1651  hTmpPerMod->SetTitle(title);
1652  hTmpPerMod->SetTitleSize(0.1);
1653  hTmpPerMod->GetXaxis()->SetTitleOffset(1.1);
1654  hTmpPerMod->GetXaxis()->SetTitleSize(0.05);
1655  hTmpPerMod->GetXaxis()->SetLabelSize(0.06);
1656  hTmpPerMod->GetYaxis()->SetTitleOffset(1.1);
1657  hTmpPerMod->GetYaxis()->SetTitleSize(0.05);
1658  hTmpPerMod->GetYaxis()->SetLabelSize(0.06);
1659  hTmpPerMod->GetZaxis()->SetLabelSize(0.04);
1660 
1661  return hTmpPerMod;
1662 
1663 }
1664 
1665 //---------------------------------------------------------------------------------------------------
1666 TH2* AutoZoom(TH2* H, Option_t* aType, Int_t EntryMin){
1667 
1668  Int_t shiftX = (Int_t)(H->GetNbinsX()/30.);
1669  Int_t shiftY = (Int_t)(H->GetNbinsY()/30.);
1670 
1671  TString opt = aType;
1672  opt.ToLower();
1673 
1674  int minX = 0;
1675  int maxX = H->GetNbinsX();
1676  int New_minX = minX;
1677  int New_maxX = maxX;
1678 
1679  int minY = 0;
1680  int maxY = H->GetNbinsY();
1681  int New_minY = minY;
1682  int New_maxY = maxY;
1683 
1684  if(opt.Contains("all")) opt = TString("minx, maxx, miny, maxy");
1685 
1686  if(opt.Contains("maxx")){
1687  for(New_maxX = maxX;New_maxX >=minX; New_maxX--){
1688  Stat_t c = 0;
1689  for(int i_y = maxY; i_y >= minY;i_y--){
1690  c = H->GetBinContent(New_maxX, i_y);
1691  if(c>EntryMin)
1692  break;
1693  }
1694  if(c>EntryMin)
1695  break;
1696  }
1697  }
1698 
1699  if(opt.Contains("maxy")){
1700  for(New_maxY = maxY;New_maxY >=minY;New_maxY--){
1701  Stat_t c = 0;
1702  for(int i_x=maxX; i_x>=minX;i_x--){
1703  c = H->GetBinContent(i_x, New_maxY );
1704  if(c>EntryMin)
1705  break;
1706  }
1707  if(c>EntryMin)
1708  break;
1709  }
1710  }
1711 
1712  if(opt.Contains("minx")){
1713  for(New_minX = minX;New_minX <=maxX; New_minX++){
1714  Stat_t c = 0;
1715  for(int i_y = minY; i_y <= maxY;i_y++){
1716  c = H->GetBinContent(New_minX, i_y);
1717  if(c>EntryMin)
1718  break;
1719  }
1720  if(c>EntryMin)
1721  break;
1722  }
1723  }
1724 
1725  if(opt.Contains("miny")){
1726  for(New_minY = minY;New_minY <=maxY;New_minY++){
1727  Stat_t c = 0;
1728  for(int i_x=minX; i_x<=maxX;i_x++){
1729  c = H->GetBinContent(i_x, New_minY );
1730  if(c>EntryMin)
1731  break;
1732  }
1733  if(c>EntryMin)
1734  break;
1735  }
1736  }
1737 
1738  if(New_maxX!=-1 && New_maxY!=-1)
1739  H->GetXaxis()->SetRange(New_minX - shiftX, New_maxX + shiftX);
1740  if(New_maxX!=-1 && New_maxY!=-1)
1741  H->GetYaxis()->SetRange(New_minY - shiftY, New_maxY + shiftY);
1742 
1743  return H;
1744 
1745 }
1746 
1747 //----------------------------------------------------------------------------------------------------
1748 int FindNumberOfSM(TFile* f, TString fTrigger, TString period){
1749 
1750  TString direct;
1751  if(!fTrigger.Contains("QA")){
1752  direct = "CaloQA_";
1753  }
1754  direct += fTrigger;
1755 
1756  Int_t nSMt=-1;
1757  Int_t year = 2000 + TString(period(3, 2)).Atoi();
1758  if ( year == 2010 ) { nSMt=4; }
1759  else if( year == 2011 || year == 2012 ){ nSMt=10; }
1760  else if( year == 2013 || year == 2014 ){ nSMt=12; }
1761  else { nSMt=20; }
1762 
1763  TList *outputList;
1764  Bool_t dirok = f->cd(direct);
1765 
1766  // ////
1767  // cout << "==============================================================" << endl;
1768  // cout << "direct = " << direct << endl;
1769  // cout << "dirok = " << dirok << endl;
1770  // cout << "==============================================================" << endl;
1771  // ////
1772 
1773  if(!dirok)
1774  Error(__FUNCTION__, Form("No input directory %s, the number SMs will be returned based on the year!", direct.Data()));
1775  else
1776  outputList = (TList*)gDirectory->Get(direct);
1777 
1778  if(!outputList)
1779  Error(__FUNCTION__, "No input list, the number SMs will be returned based on the year! ");
1780  else{
1781  outputList->SetOwner();
1782  TH2F* hNSM =(TH2F *)outputList->FindObject("EMCAL_hE_Mod");
1783  if(!hNSM || (!hNSM->GetEntries()))
1784  Error(__FUNCTION__, "hNSM Histogram not found or it is empty, the number SMs will be returned based on the year!");
1785  else
1786  nSMt = hNSM->GetYaxis()->GetBinUpEdge(hNSM->FindLastBinAbove(0, 2));
1787  }
1788 
1789  if(outputList)
1790  outputList->Delete();
1791 
1792  return nSMt;
1793 
1794 }
FormatRunHisto
TH2F * FormatRunHisto(TH2F *aHisto, const char *title, const char *YTitle="")
Definition: CreateEMCALRunQA.C:1623
filename
const char * filename
Definition: TestFCM.C:1
Double_t
double Double_t
Definition: External.C:58
QAPATH
TString QAPATH
Definition: CreateEMCALRunQA.C:47
TH2F
Definition: External.C:236
title
const char * title
Definition: MakeQAPdf.C:27
energy
energy
Definition: HFPtSpectrum.C:44
gSystem
TSystem * gSystem
DrawOccupancy
Int_t DrawOccupancy(Long_t run, TString period, TString pass, TString fTrigger, TString system, TFile *f, TFile *fout, AliEMCALGeometry *geom, Int_t SavePlots)
Definition: CreateEMCALRunQA.C:186
c
TCanvas * c
Definition: TestFitELoss.C:172
FindNumberOfSM
int FindNumberOfSM(TFile *f, TString fTrigger, TString period)
Definition: CreateEMCALRunQA.C:1748
pi0massP1
Double_t pi0massP1(Double_t *x, Double_t *par)
Definition: CreateEMCALRunQA.C:82
pi0massP2
Double_t pi0massP2(Double_t *x, Double_t *par)
Definition: CreateEMCALRunQA.C:66
DrawRun
Int_t DrawRun(Long_t run, TString period, TString pass, TString fTrigger, TFile *f, TFile *fout, Int_t SavePlots, Int_t nSM, Bool_t kFilter)
Definition: CreateEMCALRunQA.C:373
QAPATHF
TString QAPATHF
Definition: CreateEMCALRunQA.C:48
Int_t
int Int_t
Definition: External.C:63
TrendingEMCALTree
Int_t TrendingEMCALTree(Long_t RunId, TString fCalorimeter, TString system, TString period, TString pass, const int n, TList *TriggersList, TFile *f, TFile *fout, Int_t SavePlots)
Definition: CreateEMCALRunQA.C:733
TH2D
Definition: External.C:228
AutoZoom
TH2 * AutoZoom(TH2 *H, Option_t *aType="all", Int_t EntryMin=0)
Definition: CreateEMCALRunQA.C:1666
set_plot_style
void set_plot_style()
Definition: CreateEMCALRunQA.C:51
HistoPerMod
TH2F * HistoPerMod(TH2F *name, const char *title)
Definition: CreateEMCALRunQA.C:1643
Data
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
Definition: charmCutsOptimization.C:371
TH2
Definition: External.C:220
fitE
Double_t fitE(Double_t *x, Double_t *par)
Definition: CreateEMCALRunQA.C:92
Option_t
const char Option_t
Definition: External.C:48
Bool_t
bool Bool_t
Definition: External.C:53
fout
TFile * fout
input train file
Definition: DrawAnaCaloTrackQA.C:90
CreateEMCALRunQA
int CreateEMCALRunQA(const char *filename, TString RunId, TString period, TString pass, Int_t SavePlots=0, Bool_t filter=0, TString fTrigger="", TString system="", TString fCalorimeter="EMCAL")
Definition: CreateEMCALRunQA.C:103