AliPhysics  32e057f (32e057f)
DrawAnaCaloTrackQA.C
Go to the documentation of this file.
1 
40 //---------------------------------------------------------
41 // Set includes and declare methods for compilation
42 
43 #if !defined(__CINT__) || defined(__MAKECINT__)
44 
45 #include <TFile.h>
46 #include <TDirectoryFile.h>
47 #include <TList.h>
48 #include <TString.h>
49 #include <TROOT.h>
50 #include <TStyle.h>
51 #include <TH1F.h>
52 #include <TH2F.h>
53 #include <TCanvas.h>
54 #include <TPad.h>
55 #include <TLegend.h>
56 #include <TObject.h>
57 #include <TAxis.h>
58 #include <TGaxis.h>
59 #include <TLine.h>
60 #include <TF1.h>
61 #include <TMath.h>
62 #include <TLatex.h>
63 
64 #endif
65 
66 void ProcessTrigger(TString trigName = "default",
67  Bool_t checkList = kTRUE);
68 void CaloQA (Int_t icalo);
69 void TrackQA ();
70 void Pi0QA (Int_t icalo);
71 void IsolQA (Int_t icalo);
72 void CorrelQA (Int_t icalo);
73 void MCQA (Int_t icalo);
74 void ScaleAxis (TAxis *a, Double_t scale);
75 void ScaleXaxis(TH1 *h, Double_t scale);
76 void MyLatexMakeUp(TLatex *currentLatex, Int_t textfont, Double_t textsize, Int_t textcolor);
77 
78 TObject * GetHisto (TString histoName);
79 void SaveHisto (TObject* histo, Bool_t tag = kTRUE);
80 void SaveCanvas(TCanvas* canvas );
81 Bool_t GetList (TString trigName );
82 //
83 //---------------------------------------------------------
84 
85 //-----------------------
86 // Some global variables
87 TDirectoryFile *dir = 0;
88 TList *histArr = 0;
89 TFile *file = 0;
90 TFile *fout = 0;
92 TString format = "eps";
94 TLatex *text[5];
96 Bool_t addOkFlag = kTRUE;
97 
99 Int_t color[]={kBlack,kRed,kOrange+1,kYellow+1,kGreen+2,kBlue,kCyan+1,kViolet,kMagenta+2,kGray,kCyan-2,kViolet-2};
100 //
101 //-----------------------
102 
119 //_______________________________________________________________________
121 (
122  TString listName = "Pi0IM_GammaTrackCorr_EMCAL",
123  TString fileName = "AnalysisResults.root",
124  Int_t exportTo = 1,
125  TString fileFormat = "eps",
126  TString outFileName = "CaloTrackCorrQA_output",
127  Bool_t okFlag = kTRUE
128 )
129 {
130  format = fileFormat;
131  exportToFile = exportTo;
132  addOkFlag = okFlag;
133 
134  printf("Open <%s>; Get Trigger List : <%s>; Export option <%d>; format %s; outputFileName %s.root, ok flag %d\n",
135  fileName.Data(),listName.Data(),exportToFile, format.Data(),outFileName.Data(),addOkFlag);
136 
137  // Get file and list container, global variables
138  //
139  file = new TFile(fileName,"read");
140  if ( !file )
141  {
142  printf("File not found, do nothing\n");
143  return;
144  }
145 
146  dir = (TDirectoryFile*) file->Get(listName);
147  if ( !dir )
148  {
149  printf("DirectoryFile not found, do nothing\n");
150  return;
151  }
152 
153  //---------------
154  // output file with plots
155  if(exportToFile == 1)
156  {
157  fout = TFile::Open(Form("%s.root",outFileName.Data()),"UPDATE");
158  if(!fout)
159  fout = new TFile(Form("%s.root",outFileName.Data()),"RECREATE");
160 
161  //fout->ls();
162 
163  TDirectoryFile *cdd = (TDirectoryFile*)fout->Get("GA");
164  if(!cdd)
165  {
166  printf("Warning: GA <dir> doesn't exist, creating a new one");
167  cdd = (TDirectoryFile*)fout->mkdir("GA");
168  }
169  cdd->cd();
170  cdd->ls();
171  }
172  //---------------
173 
174  // Initialize OK messages
175  //
176  if(addOkFlag)
177  {
178  text[0] = new TLatex(0.35,0.4,"#color[3]{OK}"); MyLatexMakeUp(text[0],42,0.1,1);
179  text[1] = new TLatex(0.30,0.4,"#color[2]{NOT OK}"); MyLatexMakeUp(text[1],42,0.1,2);
180  text[2] = new TLatex(0.20,0.4,"#color[6]{Likely OK}"); MyLatexMakeUp(text[2],42,0.1,3);
181  text[3] = new TLatex(0.20,0.4,"#color[4]{Expert plot}"); MyLatexMakeUp(text[3],42,0.1,4);
182  text[4] = new TLatex(0.20,0.4,"#color[6]{Low stat}"); MyLatexMakeUp(text[4],42,0.1,3);
183  }
184 
185  // Process each of the triggers
186  //
187  ProcessTrigger("default" );
188  ProcessTrigger("EMCAL_L0");
189  ProcessTrigger("EMCAL_L1");
190  ProcessTrigger("EMCAL_L2");
191  ProcessTrigger("DCAL_L0" );
192  ProcessTrigger("DCAL_L1" );
193  ProcessTrigger("DCAL_L2" );
194 
195  if(exportToFile == 1)
196  {
197  fout->cd();
198  fout->Close();
199  }
200 
201  file->Close();
202 }
203 
221 //_______________________________________________________________________
222 void ProcessTrigger( TString trigName, Bool_t checkList)
223 {
224  // Access the list of histograms, global variables
225  //
226  if(checkList)
227  {
228  Bool_t ok = GetList(trigName);
229  printf("\t -- Process trigger %s, ok %d\n",trigName.Data(), ok);
230 
231  if ( !ok ) return;
232  }
233 
234  // Check the number of events, if no trigger, skip plots
235  //
236  nEvents = ((TH1D*) GetHisto("hNEvents"))->GetBinContent(1);
237  printf("\t \t n Events %2.3e\n",nEvents);
238 
239  if( nEvents < 1 )
240  {
241  printf("Skip list %s, no events\n",trigName.Data());
242  return;
243  }
244 
245  // Style options
246  //
247  gStyle->SetOptTitle(1);
248  gStyle->SetOptStat(0);
249  gStyle->SetOptFit(000000);
250  gStyle->SetPadRightMargin(0.15);
251  //gStyle->SetPadTopMargin(0.02);
252  //gStyle->SetPadLeftMargin(0.15);
253  gStyle->SetTitleFontSize(0.05);
254 
255  // Calorimeter, DCal or EMCal, to be checked, depends on trigger
256  //
257  Int_t nCalo = 2;
258  Int_t calo = 0;
259  if (trigName.Contains("EMCAL")) { calo = 0 ; nCalo = 1 ; }
260  else if(trigName.Contains("DCAL" )) { calo = 1 ; nCalo = 2 ; }
261 
262  TString caloString [] = {"EMCAL","DCAL"};
263 
264  histoTag = trigName;
265 
267  // Plotting
269 
270  // Plot basic Track QA
271  TrackQA();
272 
273  for(Int_t icalo = calo; icalo < nCalo; icalo++)
274  {
275  if(trigName.Contains("default")) histoTag=Form("%s_%s",caloString[icalo].Data(),trigName.Data());
276 
277  // Plot basic QA
278  CaloQA(icalo);
279 
280  // Plot basic Pi0 QA
281  Pi0QA(icalo);
282 
283  gStyle->SetPadRightMargin(0.02);
284  // Plot basic isolation QA
285  IsolQA(icalo);
286 
287  // Plot basic correlation QA
288  CorrelQA(icalo);
289 
290  // MC basic QA plots, cluster origins (only if it run on MC)
291  MCQA(icalo);
292 
293  // Re-set default setting
294  gStyle->SetPadRightMargin(0.15);
295  }
296 }
297 
306 //______________________________________
307 void CaloQA(Int_t icalo)
308 {
309  Int_t ok = 0;
310 
311  //-----------------------------
312  // Cluster spectra and track match residuals
313  //
314  TCanvas * ccalo = new TCanvas(Form("%s_CaloHisto_SpectraTM" ,histoTag.Data()),
315  Form("Cluster spectra and track match residuals for %s",histoTag.Data()),
316  1000,1000);
317  ccalo->Divide(2,2);
318 
319  ccalo->cd(1);
320  gPad->SetLogy();
321  gPad->SetLogx();
322 
323  TH1F* hClusterEnergy = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPt_Cut_6_Fidutial",icalo));
324  if(!hClusterEnergy) return;
325  hClusterEnergy->SetYTitle("Entries");
326  hClusterEnergy->SetTitle("Cluster-cell energy spectra");
327  hClusterEnergy->SetTitleOffset(1.5,"Y");
328  hClusterEnergy->Sumw2();
329  hClusterEnergy->SetMarkerColor(1);
330  hClusterEnergy->SetMarkerStyle(20);
331  hClusterEnergy->SetAxisRange(0.,50.,"X");
332  hClusterEnergy->Draw();
333 
334  TLegend l(0.15,0.15,0.3,0.3);
335  l.SetTextSize(0.04);
336  l.AddEntry(hClusterEnergy,"Good Cluster","P");
337  l.SetBorderSize(0);
338  l.SetFillColor(0);
339 
340  TH2F* h2CellAmplitude = (TH2F*) GetHisto("QA_Cell_hAmp_Mod");
341  TH1F* hCellAmplitude = 0;
342  if(h2CellAmplitude)
343  {
344  if(histoTag.Contains("default"))
345  {
346  if ( icalo == 0 ) hCellAmplitude = (TH1F*) h2CellAmplitude->ProjectionX(Form("%s_hCellAmp",histoTag.Data()), 1,12);
347  else hCellAmplitude = (TH1F*) h2CellAmplitude->ProjectionX(Form("%s_hCellAmp",histoTag.Data()),12,20);
348  }
349  else hCellAmplitude = (TH1F*) h2CellAmplitude->ProjectionX(Form("%s_hCellAmp",histoTag.Data()),0,100);
350 
351  hCellAmplitude->Sumw2();
352  hCellAmplitude->SetMarkerColor(4);
353  hCellAmplitude->SetMarkerStyle(25);
354  hCellAmplitude->Draw("same");
355  l.AddEntry(hCellAmplitude,"Cell","P");
356  }
357 
358  l.Draw();
359 
360  // ok message
361  // check for sudden spikes in the spectrum, due to bad channels
362  if(addOkFlag)
363  {
364  Float_t minClusterE = 1;
365  if ( histoTag.Contains("L0") ) minClusterE = 4;
366  else if ( histoTag.Contains("L2") ) minClusterE = 10;
367  else if ( histoTag.Contains("L1") ) minClusterE = 14;
368 
369  ok=0;
370 
371  if ( hClusterEnergy->GetEntries() < 1000 ) ok = 4;
372  else
373  {
374  for(Int_t ibin = 1; ibin < hClusterEnergy->GetNbinsX(); ibin++)
375  {
376  if ( hClusterEnergy->GetBinCenter(ibin) < minClusterE ) continue;
377 
378  if ( hClusterEnergy->GetBinContent(ibin) < 100 ) continue;
379 
380  if ( hClusterEnergy->GetBinContent(ibin)*1.4 < hClusterEnergy->GetBinContent(ibin+1) )
381  {
382  ok=2;
383  //printf("ibin %d, E %2.2f, %2.1f < 1.5* %2.1f\n",ibin, hClusterEnergy->GetBinCenter(ibin), hClusterEnergy->GetBinContent(ibin), hClusterEnergy->GetBinContent(ibin+1));
384  break;
385  }
386  }
387  }
388  text[ok]->Draw();
389  }
390 
391  ccalo->cd(2);
392  //gPad->SetLogy();
393  gPad->SetLogx();
394 
395  //TH1F* hRaw = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPt_Cut_6_Fidutial",icalo));
396  //TH1F* hCorr = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPt_Cut_4_NCells" ,icalo));
397  TH1F* hTM = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPt_Cut_7_Matching",icalo));
398  TH1F* hShSh = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPt_Cut_9_PID" ,icalo));
399 
400  //hRaw->Sumw2();
401 
402 // hCorr->SetTitle("Ratio after cluster cuts application");
403 // hCorr->SetYTitle("Selected clusters / Raw clusters");
404 // hCorr->SetTitleOffset(1.5,"Y");
405 // hCorr->Sumw2();
406 // hCorr->SetMarkerColor(1);
407 // hCorr->SetMarkerStyle(20);
408 // hCorr->Divide(hRaw);
409 // hCorr->SetAxisRange(0.,30.,"X");
410 // hCorr->SetMaximum(1.1);
411 // hCorr->SetMinimum(0);
412 // hCorr->Draw();
413 
414  hTM ->SetTitle("Ratio after cluster cuts application");
415  hTM ->SetYTitle("Selected clusters / Good clusters");
416  hTM ->SetTitleOffset(1.5,"Y");
417  hTM ->Sumw2();
418  hTM ->SetAxisRange(0.,50.,"X");
419  hTM ->SetMarkerColor(2);
420  hTM ->SetMarkerStyle(21);
421  hTM ->SetMaximum(1.1);
422  hTM ->SetMinimum(0);
423  hTM ->Divide(hClusterEnergy);
424  hTM ->Draw("");
425 
426  hShSh->Sumw2();
427  hShSh->SetMarkerColor(4);
428  hShSh->SetMarkerStyle(22);
429  hShSh->Divide(hClusterEnergy);
430  hShSh->Draw("same");
431 
432  TLegend l2(0.15,0.15,0.3,0.3);
433  l2.SetTextSize(0.04);
434  //l2.AddEntry(hCorr,"No Exotics + non lin.","P");
435  l2.AddEntry(hTM, "+ Track matching","P");
436  l2.AddEntry(hShSh,"+ #lambda^{2}_{0} < 0.4","P");
437  l2.SetBorderSize(0);
438  l2.SetFillColor(0);
439  l2.Draw();
440 
441  // ok message
442  // Check that the ratio neutral clusters over input cluster and neutral and photonic over input
443  // is what is expected. Currently mild guesses
444  if(addOkFlag)
445  {
446  ok=0;
447  for(Int_t ibin = 1; ibin < hTM->GetNbinsX(); ibin++)
448  {
449  if(hTM->GetBinCenter(ibin) < 0.5) continue;
450  if(hTM->GetBinCenter(ibin) > 20 ) break;
451 
452  Float_t conTM = hTM ->GetBinContent(ibin);
453  Float_t conSh = hShSh->GetBinContent(ibin);
454  Float_t errTM = hTM ->GetBinError(ibin);
455  Float_t errSh = hShSh->GetBinError(ibin);
456 
457  if(errTM > conTM || errSh > conSh) continue;
458 
459  //printf("bin %d, center %2.2f, cont TM %f, Sh %f\n",ibin,hTM->GetBinCenter(ibin),conTM,conTM);
460 
461  if(conTM > 0.999 || conSh > 0.999 || conTM < 0.5 || conSh < 0.2)
462  {
463  ok=1;
464  //printf("bin %d, center %2.2f, cont TM %f, Sh %f\n",ibin,hTM->GetBinCenter(ibin),conTM,conTM);
465  break;
466  }
467  }
468 
469  text[ok]->Draw();
470  }
471 
472  // Plot track-matching residuals
473  TH1F* hTrackMatchResEtaNeg = NULL;
474  TH1F* hTrackMatchResEtaPos = NULL;
475  TH1F* hTrackMatchResPhiNeg = NULL;
476  TH1F* hTrackMatchResPhiPos = NULL;
477 
478  // first test did not have this histogram, add protection
479  TH2F* hTrackMatchResEtaPhi = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hTrackMatchedDEtaDPhiPosNoCut",icalo));
480  if(hTrackMatchResEtaPhi)
481  {
482  hTrackMatchResEtaPhi ->Add( (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hTrackMatchedDEtaDPhiNegNoCut",icalo) ));
483 
484  ccalo->cd(3);
485  gPad->SetLogz();
486 
487  hTrackMatchResEtaPhi->SetAxisRange(-0.025,0.025,"X");
488  hTrackMatchResEtaPhi->SetAxisRange(-0.025,0.025,"Y");
489  hTrackMatchResEtaPhi->SetTitleOffset(1.5,"Y");
490  hTrackMatchResEtaPhi->SetTitleOffset(1.5,"Z");
491  hTrackMatchResEtaPhi->SetTitle("Track-cluster residual #Delta#varphi vs #Delta#eta, #it{E}>0.5 GeV");
492  hTrackMatchResEtaPhi->SetXTitle("#Delta #eta");
493  hTrackMatchResEtaPhi->SetYTitle("#Delta #varphi");
494  hTrackMatchResEtaPhi->SetZTitle("Entries");
495  hTrackMatchResEtaPhi->Draw("colz");
496 
497  if(addOkFlag) text[3]->Draw();
498 
499  ccalo->cd(4);
500  //gPad->SetLogy();
501  TGaxis::SetMaxDigits(3);
502 
503  TH2F* h2TrackMatchResEtaNeg = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hTrackMatchedDEtaNegNoCut",icalo));
504  TH2F* h2TrackMatchResEtaPos = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hTrackMatchedDEtaPosNoCut",icalo));
505  TH2F* h2TrackMatchResPhiNeg = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hTrackMatchedDPhiNegNoCut",icalo));
506  TH2F* h2TrackMatchResPhiPos = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hTrackMatchedDPhiPosNoCut",icalo));
507 
508  Float_t binMin = hClusterEnergy->FindBin(0.5);
509  hTrackMatchResEtaNeg = (TH1F*) h2TrackMatchResEtaNeg->ProjectionY(Form("%s_hTrackMatchProjClusEnEtaNeg",histoTag.Data()),binMin, 1000);
510  hTrackMatchResEtaPos = (TH1F*) h2TrackMatchResEtaPos->ProjectionY(Form("%s_hTrackMatchProjClusEnEtaPos",histoTag.Data()),binMin, 1000);
511  hTrackMatchResPhiNeg = (TH1F*) h2TrackMatchResPhiNeg->ProjectionY(Form("%s_hTrackMatchProjClusEnPhiNeg",histoTag.Data()),binMin, 1000);
512  hTrackMatchResPhiPos = (TH1F*) h2TrackMatchResPhiPos->ProjectionY(Form("%s_hTrackMatchProjClusEnPhiPos",histoTag.Data()),binMin, 1000);
513 
514  hTrackMatchResEtaNeg->SetXTitle("#Delta #eta, #Delta #varphi");
515  hTrackMatchResEtaNeg->SetYTitle("Entries");
516  hTrackMatchResEtaNeg->SetTitle("Track-cluster residuals, #it{E} > 1 GeV");
517  hTrackMatchResEtaNeg->SetTitleOffset(1.5,"Y");
518  hTrackMatchResEtaNeg->SetAxisRange(-0.025,0.025,"X");
519  hTrackMatchResEtaNeg->Sumw2();
520  hTrackMatchResEtaNeg->SetMarkerStyle(25);
521  hTrackMatchResEtaNeg->SetMarkerColor(2);
522  hTrackMatchResEtaNeg->Draw("");
523 
524  hTrackMatchResEtaPos->Sumw2();
525  hTrackMatchResEtaPos->SetMarkerStyle(25);
526  hTrackMatchResEtaPos->SetMarkerColor(4);
527  hTrackMatchResEtaPos->Draw("same");
528 
529  hTrackMatchResPhiNeg->Sumw2();
530  hTrackMatchResPhiNeg->SetMarkerStyle(24);
531  hTrackMatchResPhiNeg->SetMarkerColor(2);
532  hTrackMatchResPhiNeg->Draw("same");
533 
534  hTrackMatchResPhiPos->Sumw2();
535  hTrackMatchResPhiPos->SetMarkerStyle(24);
536  hTrackMatchResPhiPos->SetMarkerColor(4);
537  hTrackMatchResPhiPos->Draw("same");
538 
539  TLine l0(0,hTrackMatchResEtaNeg->GetMinimum(),0,hTrackMatchResEtaNeg->GetMaximum()*1.);
540  l0.Draw("same");
541 
542  TLegend l3(0.55,0.7,0.83,0.85);
543  l3.SetTextSize(0.04);
544  l3.AddEntry(hTrackMatchResEtaNeg,"#Delta #eta, Negative","P");
545  l3.AddEntry(hTrackMatchResEtaPos,"#Delta #eta, Positive","P");
546  l3.AddEntry(hTrackMatchResPhiNeg,"#Delta #varphi, Negative","P");
547  l3.AddEntry(hTrackMatchResPhiPos,"#Delta #varphi, Positive","P");
548  l3.SetBorderSize(0);
549  l3.SetFillColor(0);
550  l3.Draw();
551 
552  // ok message
553  // check that the residuals maxima are displaced from 0 within reasonable limits.
554  // quite open.
555  if(addOkFlag)
556  {
557  ok=0;
558 
559  if ( hTrackMatchResPhiPos->GetEntries() < 10 || hTrackMatchResPhiNeg->GetEntries() < 10 ||
560  hTrackMatchResEtaPos->GetEntries() < 10 || hTrackMatchResEtaNeg->GetEntries() < 10 ) ok = 1;
561 
562  // printf("Mean (%2.2f, %2.2f, %2.2f, %2.2f), RMS (%2.2f, %2.2f, %2.2f, %2.2f)\n",
563  // hTrackMatchResPhiPos->GetMean()*1000, hTrackMatchResEtaPos->GetMean()*1000,
564  // hTrackMatchResPhiNeg->GetMean()*1000, hTrackMatchResEtaNeg->GetMean()*1000,
565  // hTrackMatchResPhiPos->GetRMS ()*1000, hTrackMatchResEtaPos->GetRMS ()*1000,
566  // hTrackMatchResPhiNeg->GetRMS ()*1000, hTrackMatchResEtaNeg->GetRMS ()*1000
567  // );
568 
569  // just a guess, not based on anything
570  if( hTrackMatchResPhiPos->GetMean()*1000 < -2 || hTrackMatchResEtaPos->GetMean()*1000 < -2 ||
571  hTrackMatchResPhiNeg->GetMean()*1000 > 2 || hTrackMatchResEtaNeg->GetMean()*1000 > 2) ok = 2;
572 
573  text[ok]->Draw();
574  }
575  }
576 
577  ccalo->Print(Form("%s_CaloHisto_ClusterSpectraAndTrackResiduals.%s",histoTag.Data(),format.Data()));
578 
579  // cleanup or save
580  //
581  if(exportToFile!=1)
582  {
583  delete hCellAmplitude;
584  delete hTrackMatchResEtaNeg;
585  delete hTrackMatchResEtaPos;
586  delete hTrackMatchResPhiNeg;
587  delete hTrackMatchResPhiPos;
588 
589  delete ccalo;
590  }
591  else
592  {
593  SaveHisto(hCellAmplitude ,kFALSE);
594  SaveHisto(hTrackMatchResEtaNeg,kFALSE);
595  SaveHisto(hTrackMatchResEtaPos,kFALSE);
596  SaveHisto(hTrackMatchResPhiNeg,kFALSE);
597  SaveHisto(hTrackMatchResPhiPos,kFALSE);
598 
599  SaveCanvas(ccalo);
600  }
601 
602  //-----------------------------
603  // Cell and cluster hit maps
604  //
605  TCanvas * ccalo2 = new TCanvas(Form("%s_CaloHisto_CellClusterHit" ,histoTag.Data()),
606  Form("Cell and cluster hit maps for %s",histoTag.Data()),
607  1000,1000);
608  ccalo2->Divide(2,2);
609 
610  ccalo2->cd(1);
611  gPad->SetLogz();
612 
613  TH2F* hCellActivity = (TH2F*) GetHisto("QA_Cell_hGridCells");
614  if(hCellActivity)
615  {
616  if(icalo == 0) hCellActivity->SetAxisRange( 0,127,"Y");
617  else hCellActivity->SetAxisRange(128,220,"Y");
618  hCellActivity->SetTitle("Hits per cell (#it{E} > 0.2 GeV)");
619  hCellActivity->SetTitleOffset(1.5,"Y");
620  hCellActivity->SetZTitle("Entries");
621  hCellActivity->SetTitleOffset(1.5,"Z");
622  hCellActivity->Draw("colz");
623 
624  // ok message
625  // Check that the dead regions of the EMCal are not a too large fraction
626  if(addOkFlag)
627  {
628  //Float_t meanCellHits = 0;
629  Int_t nHitCell = 0;
630  //Int_t nHotCell = 0;
631  //Int_t nColdCell = 0;
632  ok=0;
633  for(Int_t ibin = 2; ibin < hCellActivity->GetNbinsX(); ibin++)
634  {
635  for(Int_t jbin = 2; jbin < hCellActivity->GetNbinsY(); jbin++)
636  {
637  if(histoTag.Contains("default"))
638  {
639  if ( icalo == 0 && jbin > 129 ) continue ; // EMC
640  if ( icalo == 1 && jbin < 130 ) continue ; // DMC
641  }
642 
643  //meanCellHits += hCellActivity->GetBinContent(ibin,jbin);
644  if ( hCellActivity->GetBinContent(ibin,jbin) > 2 ) nHitCell++;
645  }
646  }
647 
648  //if ( nHitCell > 0 ) meanCellHits/=nHitCell;
649 
650  Float_t totalCells = 10*24*48+2*8*48; // EMC
651  if(icalo==1) totalCells = 6*24*32+2*8*48; // DMC
652 
653  Float_t active = nHitCell / totalCells ;
654 
655  if ( active > 0.90) ok = 0;
656  else if ( active > 0.7 ) ok = 2;
657  else ok = 1;
658 
659  text[ok]->Draw();
660 
661  // for(Int_t ibin = 1; ibin < hCellActivity->GetNbinsX(); ibin++)
662  // {
663  // for(Int_t jbin = 1; jbin < hCellActivity->GetNbinsY(); jbin++)
664  // {
665  // if(icalo == 0 && (jbin > 127)) continue ; // EMC
666  // if(icalo == 1 && (jbin <= 127)) continue ; // DMC
667  //
668  // if ( hCellActivity->GetBinContent(ibin,jbin) > meanCellHits*10 ) nHotCell++;
669  // if ( hCellActivity->GetBinContent(ibin,jbin) < meanCellHits*0.1 ) nColdCell++;
670  // }
671  // }
672  //
673  // printf("cell hits %d, total %d, active %2.2f, mean %f, cold %d, hot %d\n",nHitCell,totalCells,active,meanCellHits,nHotCell,nColdCell);
674  }
675  }
676 
677  ccalo2->cd(2);
678 
679  TH2F* hCellActivityE = (TH2F*) GetHisto("QA_Cell_hGridCellsE");
680  if(hCellActivityE)
681  {
682  if(icalo == 0) hCellActivityE->SetAxisRange( 0,127,"Y");
683  else hCellActivityE->SetAxisRange(128,220,"Y");
684 
685  hCellActivityE->SetTitle("Mean energy per cell (#it{E} > 0.2 GeV)");
686 
687  if(icalo != 1 && !histoTag.Contains("default")) // ratio already done for calo=0
688  hCellActivityE->Divide(hCellActivity);
689 
690  hCellActivityE->SetTitleOffset(1.5,"Y");
691  hCellActivityE->SetZTitle("#Sigma #it{E}_{cell} / Entries_{per cell}");
692  hCellActivityE->SetTitleOffset(1.5,"Z");
693 
694  hCellActivityE->Draw("colz");
695 
696  // ok message
697  // Check that the dead regions of the EMCal are not a too large fraction
698  if(addOkFlag)
699  {
700  //Float_t meanCellHits = 0;
701  Int_t nHitCell = 0;
702  //Int_t nHotCell = 0;
703  //Int_t nColdCell = 0;
704  ok=0;
705  for(Int_t ibin = 2; ibin < hCellActivityE->GetNbinsX(); ibin++)
706  {
707  //printf("col %d ",ibin);
708  for(Int_t jbin = 2; jbin < hCellActivityE->GetNbinsY(); jbin++)
709  {
710  if(histoTag.Contains("default"))
711  {
712  if ( icalo == 0 && jbin > 129 ) continue ; // EMC
713  if ( icalo == 1 && jbin < 130 ) continue ; // DMC
714  }
715 
716  //meanCellHits += hCellActivityE->GetBinContent(ibin,jbin);
717  if ( hCellActivityE->GetBinContent(ibin,jbin) > 0 )
718  {
719  //printf("1");
720  nHitCell++;
721  }
722  //else printf("0");
723  }
724  //printf("\n");
725  }
726 
727  //if ( nHitCell > 0 ) meanCellHits/=nHitCell;
728 
729  Float_t totalCells = 10*24*48+2*8*48; // EMC
730  if(icalo==1) totalCells = 6*24*32+2*8*48; // DMC
731 
732  Float_t active = nHitCell / totalCells ;
733 
734  if ( active > 0.90) ok = 0;
735  else if ( active > 0.7 ) ok = 2;
736  else ok = 1;
737 
738  text[ok]->Draw();
739  //printf("cell hits %d, total %2.0f, active %2.2f\n",nHitCell,totalCells,active);
740 
741  // for(Int_t ibin = 1; ibin < hCellActivityE->GetNbinsX(); ibin++)
742  // {
743  // for(Int_t jbin = 1; jbin < hCellActivityE->GetNbinsY(); jbin++)
744  // {
745  // if(icalo == 0 && (jbin > 127)) continue ; // EMC
746  // if(icalo == 1 && (jbin <= 127)) continue ; // DMC
747  //
748  // if ( hCellActivityE->GetBinContent(ibin,jbin) > meanCellHits*10 ) nHotCell++;
749  // if ( hCellActivityE->GetBinContent(ibin,jbin) < meanCellHits*0.1 ) nColdCell++;
750  // }
751  // }
752  //
753  // printf("cellE hits %d, total %d, active %2.2f, mean %f, cold %d, hot %d\n",nHitCell,totalCells,active,meanCellHits,nHotCell,nColdCell);
754  }
755  }
756 
757  ccalo2->cd(3);
758  gPad->SetLogz();
759 
760  TH2F* hClusterActivity = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hEBin0_Cluster_ColRow_PID",icalo));
761 
762  if(histoTag.Contains("default")) hClusterActivity->SetTitle("Clusters per col-row 0.5<#it{E}<3 GeV");
763  else if(histoTag.Contains("L0")) hClusterActivity->SetTitle("Clusters per col-row 2<#it{E}<5 GeV");
764  else hClusterActivity->SetTitle("Clusters per col-row 5<#it{E}<12 GeV");
765 
766  hClusterActivity->SetTitleOffset(1.5,"Y");
767  hClusterActivity->SetZTitle("Entries");
768  hClusterActivity->SetTitleOffset(1.5,"Z");
769 
770  hClusterActivity->Draw("colz");
771 
772  // ok message
773  // Check that the dead regions of the EMCal are not a too large fraction
774  if(addOkFlag)
775  {
776  //Float_t meanCellHits = 0;
777  Int_t nHitCell = 0;
778  //Int_t nHotCell = 0;
779  //Int_t nColdCell = 0;
780  ok=0;
781  for(Int_t ibin = 1; ibin < hClusterActivity->GetNbinsX(); ibin++)
782  {
783  for(Int_t jbin = 1; jbin < hClusterActivity->GetNbinsY(); jbin++)
784  {
785  //meanCellHits += hClusterActivity->GetBinContent(ibin,jbin);
786  if ( hClusterActivity->GetBinContent(ibin,jbin) > 2 ) nHitCell++;
787  }
788  }
789 
790  //if ( nHitCell > 0 ) meanCellHits/=nHitCell;
791 
792  Float_t totalCells = 10*24*48+2*8*48 -48*12-10*24-4*8; // EMC minus close to border
793  if(icalo==1) totalCells = 6*24*32+2*8*48 -32*8-48*4-6*24-8*2; // DMC minus close to border
794 
795  Float_t active = nHitCell / totalCells ;
796 
797  if ( active > 0.90) ok = 0;
798  else if ( active > 0.7 ) ok = 2;
799  else ok = 1;
800 
801  if(ok==1 && nEvents < 1000000) ok = 2;
802 
803  text[ok]->Draw();
804  //printf("cellE hits %d, total %2.0f, active %2.2f\n",nHitCell,totalCells,active);
805 
806  //printf("cluster cell hits %d, total %d, active %2.2f\n",nHitCell,totalCells,active);
807 
808  // for(Int_t ibin = 1; ibin < hClusterActivity->GetNbinsX(); ibin++)
809  // {
810  // for(Int_t jbin = 1; jbin < hClusterActivity->GetNbinsY(); jbin++)
811  // {
812  // if ( hClusterActivity->GetBinContent(ibin,jbin) > meanCellHits*10 ) nHotCell++;
813  // if ( hClusterActivity->GetBinContent(ibin,jbin) < meanCellHits*0.1 ) nColdCell++;
814  // }
815  // }
816  //
817  // printf("cluster cell hits %d, total %d, active %2.2f, mean %f, cold %d, hot %d\n",nHitCell,totalCells,active,meanCellHits,nHotCell,nColdCell);
818  }
819 
820  ccalo2->cd(4);
821  gPad->SetLogz();
822 
823  TH2F* hClusterActivity2 = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hEBin1_Cluster_ColRow_PID",icalo));
824 
825  if(histoTag.Contains("default")) hClusterActivity2->SetTitle("Clusters per col-row #it{E} > 3 GeV");
826  else if(histoTag.Contains("L0")) hClusterActivity2->SetTitle("Clusters per col-row #it{E} > 5 GeV");
827  else hClusterActivity2->SetTitle("Clusters per col-row #it{E} > 12 GeV");
828 
829  hClusterActivity2->SetTitleOffset(1.5,"Y");
830  hClusterActivity2->SetZTitle("Entries");
831  hClusterActivity2->SetTitleOffset(1.5,"Z");
832 
833  hClusterActivity2->Draw("colz");
834 
835  // ok message
836  // Check that the dead regions of the EMCal are not a too large fraction
837  if(addOkFlag)
838  {
839  //Float_t meanCellHits = 0;
840  Int_t nHitCell = 0;
841  //Int_t nHotCell = 0;
842  //Int_t nColdCell = 0;
843  ok=0;
844  for(Int_t ibin = 1; ibin < hClusterActivity2->GetNbinsX(); ibin++)
845  {
846  for(Int_t jbin = 1; jbin < hClusterActivity2->GetNbinsY(); jbin++)
847  {
848  //meanCellHits += hClusterActivity2->GetBinContent(ibin,jbin);
849  if ( hClusterActivity2->GetBinContent(ibin,jbin) > 0 ) nHitCell++;
850  }
851  }
852 
853  //if ( nHitCell > 0 ) meanCellHits/=nHitCell;
854 
855  Float_t totalCells = 10*24*48+2*8*48 -48*12-10*24-4*8; // EMC minus close to border
856  if(icalo==1) totalCells = 6*24*32+2*8*48 -32*8-48*4-6*24-8*2; // DMC minus close to border
857 
858  Float_t active = nHitCell / totalCells ;
859 
860  if ( active > 0.80) ok = 0;
861  else if ( active > 0.50) ok = 2;
862  else ok = 1;
863 
864  if(ok==1 && nEvents < 1000000) ok = 2;
865 
866  text[ok]->Draw();
867  //printf("cluster2 cell hits %d, total %d, active %2.2f\n",nHitCell,totalCells,active);
868 
869 // for(Int_t ibin = 1; ibin < hClusterActivity2->GetNbinsX(); ibin++)
870 // {
871 // for(Int_t jbin = 1; jbin < hClusterActivity2->GetNbinsY(); jbin++)
872 // {
873 // if ( hClusterActivity2->GetBinContent(ibin,jbin) > meanCellHits*10 ) nHotCell++;
874 // if ( hClusterActivity2->GetBinContent(ibin,jbin) < meanCellHits*0.1 ) nColdCell++;
875 // }
876 // }
877 //
878 // printf("cluster2 cell hits %d, total %d, active %2.2f, mean %f, cold %d, hot %d\n",nHitCell,totalCells,active,meanCellHits,nHotCell,nColdCell);
879  }
880 
881  ccalo2->Print(Form("%s_CaloHisto_CellClusterHit.%s",histoTag.Data(),format.Data()));
882 
883  // cleanup or save
884  //
885  if(exportToFile!=1) delete ccalo2;
886  else SaveCanvas(ccalo2);
887 
888  //-----------------------------
889  // Cluster time, shape, ncells
890  //
891  TCanvas * ccalo3 = new TCanvas(Form("%s_CaloHisto_ClusterTimeShape" ,histoTag.Data()),
892  Form("Cluster time, shape, ncells for %s",histoTag.Data()),
893  1000,1000);
894  ccalo3->Divide(2,2);
895 
896  ccalo3->cd(1);
897 
898  gPad->SetLogz();
899 
900  TH2F* hClusterTime = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hTimePt",icalo));
901  hClusterTime->SetTitle("Cluster #it{E} vs #it{time}");
902  hClusterTime->SetYTitle("#it{time} (ns)");
903  //hClusterTime->SetAxisRange(300.,900.,"Y");
904  hClusterTime->SetAxisRange(0.,30.,"X");
905  hClusterTime->SetTitleOffset(1.5,"Y");
906  hClusterTime->SetZTitle("Entries");
907  hClusterTime->SetTitleOffset(1.5,"Z");
908  hClusterTime->Draw("colz");
909 
910  if(addOkFlag) text[3]->Draw();
911 
912  ccalo3->cd(2);
913 
914  gPad->SetLogz();
915 
916  TH2F* hClusterL0 = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hLam0E",icalo));
917  hClusterL0->SetTitle("Cluster #sigma_{long}");
918  hClusterL0->SetAxisRange(0.,30.,"X");
919  hClusterL0->SetTitleOffset(1.5,"Y");
920  hClusterL0->SetZTitle("Entries");
921  hClusterL0->SetTitleOffset(1.5,"Z");
922 
923  hClusterL0->Draw("colz");
924  if(addOkFlag) text[3]->Draw();
925 
926  ccalo3->cd(3);
927 
928  gPad->SetLogz();
929 
930  TH2F* hClusterNCell = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hNCellsE",icalo));
931  hClusterNCell->SetTitle("Number of cells in cluster");
932  hClusterNCell->SetAxisRange(0.,30.,"X");
933  hClusterNCell->SetTitleOffset(1.5,"Y");
934  hClusterNCell->SetZTitle("Entries");
935  hClusterNCell->SetTitleOffset(1.5,"Z");
936 
937  hClusterNCell->Draw("colz");
938  if(addOkFlag) text[3]->Draw();
939 
940  ccalo3->cd(4);
941 
942  gPad->SetLogz();
943 
944  TH2F* hClusterECell = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hCellsE",icalo));
945  hClusterECell->SetTitle("cells in cluster #it{E} vs cluster #it{E}");
946  hClusterECell->SetAxisRange(0.,30.,"X");
947  hClusterECell->SetAxisRange(0.,20.,"Y");
948  hClusterECell->SetTitleOffset(1.5,"Y");
949  hClusterECell->SetZTitle("Entries");
950  hClusterECell->SetTitleOffset(1.5,"Z");
951 
952  hClusterECell->Draw("colz");
953  if(addOkFlag) text[3]->Draw();
954 
955  ccalo3->Print(Form("%s_CaloHisto_TimeShapeNCells.%s",histoTag.Data(),format.Data()));
956 
957  // cleanup or save
958  //
959  if(exportToFile!=1) delete ccalo3;
960  else SaveCanvas(ccalo3);
961 
962 }
963 
970 //______________________________________
971 void TrackQA()
972 {
973  Int_t ok = 0;
974 
975  TCanvas * ctrack = new TCanvas(Form("%s_TrackHisto" ,histoTag.Data()),
976  Form("Hybrid tracks for %s",histoTag.Data()),
977  1000,1000);
978  ctrack->Divide(2,2);
979 
980  ctrack->cd(1);
981  //gPad->SetLogz();
982  TH2F * hTrackEtaPhi = (TH2F*) GetHisto("AnaHadrons_hEtaPhiNegative");
983  if(!hTrackEtaPhi) return;
984  hTrackEtaPhi ->Add( (TH2F*) GetHisto("AnaHadrons_hEtaPhiPositive"));
985  hTrackEtaPhi ->SetAxisRange(-0.9,0.9,"X");
986  hTrackEtaPhi ->SetTitleOffset(1.5,"Y");
987  hTrackEtaPhi ->SetTitle("Hybrid tracks #eta vs #varphi #it{p}_{T} > 0.2 GeV/#it{c}");
988  hTrackEtaPhi->SetZTitle("Entries");
989  hTrackEtaPhi->SetTitleOffset(1.5,"Z");
990 
991  hTrackEtaPhi ->Draw("colz");
992 
993  // ok message
994  if(addOkFlag) text[3]->Draw();
995 
996  ctrack->cd(2);
997  //gPad->SetLogy();
998  TH2F * hTrackEtaPhiSPD = (TH2F*) GetHisto("AnaHadrons_hEtaPhiSPDRefitPt02");
999  TH2F * hTrackEtaPhiNoSPD = (TH2F*) GetHisto("AnaHadrons_hEtaPhiNoSPDRefitPt02");
1000 
1001  TH1F* hPhiSPD = (TH1F*)hTrackEtaPhiSPD ->ProjectionY(Form("%s_hTrackPhiSPD" ,histoTag.Data()),0,1000);
1002  TH1F* hPhiNoSPD = (TH1F*)hTrackEtaPhiNoSPD->ProjectionY(Form("%s_hTrackPhiNoSPD",histoTag.Data()),0,1000);
1003  TH1F* hPhi = (TH1F*)hPhiSPD->Clone( Form("%s_hTrackPhi" ,histoTag.Data()));
1004  hPhi->Add(hPhiNoSPD);
1005 
1006  hPhi ->SetTitle("Hybrid track type #varphi, 0.2<#it{p}_{T}<2 GeV/#it{c}");
1007  hPhi ->SetLineColor(1);
1008  hPhiSPD ->SetLineColor(2);
1009  hPhiNoSPD->SetLineColor(4);
1010 
1011  hPhi ->SetMinimum(1);
1012  hPhi ->SetMaximum(hPhi->GetMaximum()*1.3);
1013  hPhi ->SetTitleOffset(1.5,"Y");
1014  hPhi ->SetYTitle("Entries");
1015 
1016  TGaxis::SetMaxDigits(3);
1017 
1018  hPhi ->Draw("H");
1019  hPhiSPD ->Draw("Hsame");
1020  hPhiNoSPD->Draw("Hsame");
1021 
1022  TLegend l(0.2,0.75,0.4,0.89);
1023  l.SetTextSize(0.04);
1024  l.AddEntry(hPhi,"Sum","L");
1025  l.AddEntry(hPhiSPD ,"SPD+Refit","L");
1026  l.AddEntry(hPhiNoSPD,"No SPD+Refit","L");
1027  l.SetBorderSize(0);
1028  l.SetFillColor(0);
1029  l.Draw();
1030 
1031  // ok message
1032  // check that the tracks phi distribution is uniform (only for Min Bias) within a window.
1033  // quite open
1034  if(addOkFlag)
1035  {
1036  ok=0;
1037 
1038  if ( hPhi->GetEntries() < 1000 ) ok = 4;
1039  else if ( !histoTag.Contains("default")) ok = 3;
1040  else
1041  {
1042  hPhi->Fit("pol0","QRL","",0, 7);
1043  Float_t fitParam = hPhi->GetFunction("pol0")->GetParameter(0);
1044  //printf("fit param %f\n",fitParam);
1045  for(Int_t ibin = 1; ibin < hPhi->GetNbinsX(); ibin++)
1046  {
1047  if(fitParam < 100)
1048  {
1049  ok=1;
1050  break;
1051  }
1052 
1053  if ( hPhi->GetBinContent(ibin) < 100 ) continue;
1054 
1055  Float_t ratToFit = hPhi->GetBinContent(ibin)/fitParam;
1056  //printf("\t bin %d, rat to fit %2.2f\n",ibin,ratToFit);
1057  if ( ratToFit < 0.85 || ratToFit > 1.15 ) // Just guessing, not sure it is ok.
1058  {
1059  ok=1;
1060  //printf("ibin %d, ratio %f\n",ibin, ratToFit);
1061  break;
1062  }
1063  else if ( (ratToFit < 0.98 && ratToFit >= 0.85) || (ratToFit <= 1.15 && ratToFit > 1.02) )
1064  {
1065  ok=2;
1066  }
1067  }
1068  }
1069  text[ok]->Draw();
1070  }
1071 
1072  ctrack->cd(3);
1073  gPad->SetLogy();
1074 
1075  TH1F* hTOF = (TH1F*) GetHisto("AnaHadrons_hTOFSignalPtCut");
1076  hTOF->SetYTitle("Entries");
1077  hTOF->SetTitleOffset(1.5,"Y");
1078 
1079  hTOF->Draw("");
1080 
1081  // ok message
1082  // check that most of the tracks have TOF at BC0
1083  if(addOkFlag)
1084  {
1085  ok=0;
1086 
1087  if ( hTOF->GetEntries() < 1000 ) ok = 4;
1088  else
1089  {
1090  Float_t intBC0 = hTOF->Integral(hTOF->FindBin(0),hTOF->FindBin(25));
1091  Float_t intAll = hTOF->Integral();
1092 
1093  Float_t rat = 0;
1094  if ( intAll > 0 ) rat = intBC0 / intAll;
1095  //printf("ratio TOF %f\n",rat);
1096 
1097  // Just guessing
1098  if (rat < 0.7) ok = 1;
1099  else if(rat < 0.9) ok = 2;
1100  }
1101 
1102  text[ok]->Draw();
1103  }
1104 
1105  ctrack->cd(4);
1106  gPad->SetLogy();
1107  gPad->SetLogx();
1108 
1109  TH1F* hPt = (TH1F*) GetHisto("AnaHadrons_hPt");
1110  TH1F* hPtSPD = (TH1F*) GetHisto("AnaHadrons_hPtSPDRefit");
1111  TH1F* hPtNoSPD = (TH1F*) GetHisto("AnaHadrons_hPtNoSPDRefit");
1112  hPt ->SetLineColor(1);
1113  hPtSPD ->SetLineColor(2);
1114  hPtNoSPD->SetLineColor(4);
1115 
1116  hPt ->SetTitle("Hybrid track type #it{p}_{T}");
1117  hPt ->SetYTitle("Entries");
1118  hPt ->SetTitleOffset(1.5,"Y");
1119 
1120  hPt ->Draw("");
1121  hPtSPD ->Draw("same");
1122  hPtNoSPD->Draw("same");
1123 
1124  // ok message
1125  // check for sudden spikes in the spectrum
1126  if(addOkFlag)
1127  {
1128  ok=0;
1129 
1130  if ( hPt->GetEntries() < 1000 ) ok = 4;
1131  else
1132  {
1133  for(Int_t ibin = 1; ibin < hPt->GetNbinsX(); ibin++)
1134  {
1135  if ( hPt->GetBinContent(ibin) < 100 ) continue;
1136 
1137  if ( hPt->GetBinContent(ibin)*1.4 < hPt->GetBinContent(ibin+1) )
1138  {
1139  ok=1;
1140  //printf("ibin %d, E %2.2f, %2.1f < 1.5* %2.1f\n",ibin, hClusterEnergy->GetBinCenter(ibin), hClusterEnergy->GetBinContent(ibin), hClusterEnergy->GetBinContent(ibin+1));
1141  break;
1142  }
1143  }
1144  }
1145 
1146  text[ok]->Draw();
1147  }
1148 
1149 
1150 // ctrack->cd(3);
1151 // gPad->SetLogz();
1152 //
1153 // TH2F* hPtDCAxy = (TH2F*) GetHisto("AnaHadrons_hPtDCAxy");
1154 // hPtDCAxy->SetAxisRange(-1,1,"Y");
1155 // hPtDCAxy->SetAxisRange(0,30,"X");
1156 // hPtDCAxy->Draw("colz");
1157 //
1158 // ctrack->cd(4);
1159 // gPad->SetLogz();
1160 //
1161 // TH2F* hPtDCAz = (TH2F*) GetHisto("AnaHadrons_hPtDCAz");
1162 // hPtDCAz->SetAxisRange(-1,1,"Y");
1163 // hPtDCAz->SetAxisRange(0,30,"X");
1164 // hPtDCAz->Draw("colz");
1165 
1166  ctrack->Print(Form("%s_TrackHisto.%s",histoTag.Data(),format.Data()));
1167 
1168  // cleanup or save
1169  //
1170  if(exportToFile!=1)
1171  {
1172  delete hPhi;
1173  delete hPhiSPD;
1174  delete hPhiNoSPD;
1175 
1176  delete ctrack;
1177  }
1178  else
1179  {
1180  SaveHisto(hPhi ,kFALSE);
1181  SaveHisto(hPhiSPD ,kFALSE);
1182  SaveHisto(hPhiNoSPD,kFALSE);
1183 
1184  SaveCanvas(ctrack);
1185  }
1186 }
1187 
1196 //_____________________________
1197 void Pi0QA(Int_t icalo)
1198 {
1199  Int_t ok = 0;
1200 
1201  TCanvas * cpi0 = new TCanvas(Form("%s_InvariantMassHisto" ,histoTag.Data()),
1202  Form("Neutral mesons inv. mass for %s",histoTag.Data()),
1203  1000,1000);
1204  cpi0->Divide(2,2);
1205 
1206  TH2F* hMassE[10];
1207  TH2F* hMixMassE[10];
1208  for(Int_t icen = 0; icen < 10; icen++)
1209  {
1210  hMassE [icen] = (TH2F*) GetHisto(Form("AnaPi0_Calo%d_hRe_cen%d_pidbit0_asy0_dist1",icalo,icen));
1211  hMixMassE[icen] = (TH2F*) GetHisto(Form("AnaPi0_Calo%d_hMi_cen%d_pidbit0_asy0_dist1",icalo,icen));
1212  }
1213  if(!hMassE[0]) return;
1214 
1215  // 2D Invariant mass vs E, in PbPb from 60 to 100 %, all in pp
1216  cpi0->cd(1);
1217  gPad->SetLogz();
1218  TH2F* h2DMass;
1219 
1220  if(hMassE[1]) // Plot centrality from 60 to 100%
1221  {
1222  h2DMass = (TH2F*) hMassE[6]->Clone(Form("%s_h2DMass",histoTag.Data()));
1223  for(Int_t icen = 7; icen < 10; icen++) h2DMass->Add(hMassE[icen]);
1224  h2DMass->SetTitle("Inv. mass vs #it{p}_{T,pair}, Cen: 60-100%");
1225  }
1226  else
1227  {
1228  h2DMass = (TH2F*) hMassE[0]->Clone(Form("%s_h2DMass",histoTag.Data()));
1229  h2DMass->SetTitle("Inv. mass vs #it{p}_{T,pair}");
1230  }
1231 
1232  h2DMass->SetTitleOffset(1.5,"Y");
1233  h2DMass->SetAxisRange(0.0,0.7,"Y");
1234  h2DMass->SetAxisRange(0,30,"X");
1235  h2DMass->Draw("colz");
1236 
1237  // ok message
1238  if(addOkFlag) text[3]->Draw();
1239 
1240  // Pi0 Invariant mass projection, in PbPb 6 centrality bins from 0 to 50%, all in pp
1241  cpi0->cd(2);
1242  TH1F* hMass [10];
1243  TH1F* hMix [10];
1244  TH1F* hMassEta[10];
1245  TH1F* hMassPi0[10];
1246 
1247  //Init to 0
1248  for(Int_t icen=0; icen<10; icen++ )
1249  {
1250  hMass [icen] = 0;
1251  hMix [icen] = 0;
1252  hMassEta[icen] = 0;
1253  hMassPi0[icen] = 0;
1254  }
1255 
1256  TH1F * hX = (TH1F*) hMassE[0]->ProjectionX(Form("%s_hEPairCen0",histoTag.Data()),0,10000);
1257  Int_t binmin = hX->FindBin(2); // Project histo from 2 GeV pairs
1258  Int_t binmax = hX->FindBin(10); // Project histo up to 10 GeV pairs
1259  if(histoTag.Contains("L0"))
1260  {
1261  binmin = hX->FindBin(5); // Project histo from 5 GeV pairs
1262  binmax = hX->FindBin(10); // Project histo up to 10 GeV pairs
1263  }
1264  else if(histoTag.Contains("L2"))
1265  {
1266  binmin = hX->FindBin(8); // Project histo from 8 GeV pairs
1267  binmax = hX->FindBin(12); // Project histo up to 12 GeV pairs
1268  }
1269  else if(histoTag.Contains("L1"))
1270  {
1271  binmin = hX->FindBin(10); // Project histo from 10 GeV pairs
1272  binmax = hX->FindBin(15); // Project histo up to 15 GeV pairs
1273  }
1274 
1275  Float_t maxPi0 = 0;
1276  Float_t maxEta = 0;
1277  Float_t minPi0 = 1e6;
1278  Float_t minEta = 1e6;
1279  for(Int_t icen = 0; icen < 6; icen++)
1280  {
1281  if(!hMassE[icen]) continue;
1282 
1283  hMass[icen] = (TH1F*) hMassE [icen]->ProjectionY(Form("%s_hMassCen%d",histoTag.Data(),icen),binmin,binmax);
1284  hMix [icen] = (TH1F*) hMixMassE[icen]->ProjectionY(Form("%s_hMixCen%d" ,histoTag.Data(),icen),binmin,binmax);
1285  hMass[icen]->Sumw2();
1286  hMix [icen]->Sumw2();
1287 
1288  hMassPi0[icen] = (TH1F*) hMass[icen]->Clone(Form("%s_hMassPi0Cen%d",histoTag.Data(),icen));
1289  hMassEta[icen] = (TH1F*) hMass[icen]->Clone(Form("%s_hMassEtaCen%d",histoTag.Data(),icen));
1290 
1291  hMassPi0[icen]->Divide(hMix[icen]);
1292  hMassPi0[icen]->Fit("pol0","Q","",0.25,0.35);
1293  Float_t scale = 1;
1294  if(hMassPi0[icen]->GetFunction("pol0")) scale = hMassPi0[icen]->GetFunction("pol0")->GetParameter(0);
1295  //printf("Scale factor %f for cen %d\n",scale,icen);
1296  hMassPi0[icen]->Scale(1./scale);
1297  hMassPi0[icen]->SetMarkerStyle(24);
1298  hMassPi0[icen]->SetMarkerColor(color[icen]);
1299  hMassPi0[icen]->SetLineColor(color[icen]);
1300  hMassPi0[icen]->SetAxisRange(0.04,0.24);
1301  //hMassPi0[icen]->SetMarkerSize(0.5);
1302 
1303  hMassEta[icen]->Rebin(4);
1304  hMix [icen]->Rebin(4);
1305  hMassEta[icen]->Divide(hMix[icen]);
1306  hMassEta[icen]->SetMarkerStyle(25);
1307  hMassEta[icen]->SetMarkerColor(color[icen]);
1308  hMassEta[icen]->SetLineColor(color[icen]);
1309  hMassEta[icen]->SetAxisRange(0.4,0.9);
1310  //hMassEta[icen]->SetMarkerSize(0.5);
1311  hMassEta[icen]->Scale(1./scale);
1312 
1313  if(maxEta < hMassEta[icen]->GetMaximum()) maxEta = hMassEta[icen]->GetMaximum();
1314  if(maxPi0 < hMassPi0[icen]->GetMaximum()) maxPi0 = hMassPi0[icen]->GetMaximum();
1315 
1316  if(minEta > hMassEta[icen]->GetMinimum()) minEta = hMassEta[icen]->GetMinimum();
1317  if(minPi0 > hMassPi0[icen]->GetMinimum()) minPi0 = hMassPi0[icen]->GetMinimum();
1318  }
1319 
1320  //gPad->SetLogy();
1321  //gPad->SetGridy();
1322  hMassPi0[0]->SetMinimum(minPi0);
1323  hMassPi0[0]->SetTitleOffset(1.5,"Y");
1324  hMassPi0[0]->SetYTitle("Real / Mixed");
1325  hMassPi0[0]->SetTitle("#pi^{0} peak, 2 < #it{E}_{pair}< 10 GeV");
1326  if (histoTag.Contains("L0")) hMassPi0[0]->SetTitle("#pi^{0} peak, 5 < #it{E}_{pair}< 10 GeV");
1327  else if(histoTag.Contains("L2")) hMassPi0[0]->SetTitle("#pi^{0} peak, 8 < #it{E}_{pair}< 12 GeV");
1328  else if(histoTag.Contains("L1")) hMassPi0[0]->SetTitle("#pi^{0} peak, 10 < #it{E}_{pair}< 15 GeV");
1329 
1330  hMassPi0[0]->Draw();
1331 
1332  // ok message
1333  // Check that at low or high mass with respect the pi0 peak there is not
1334  // too much activity with respect the peak, sign of bad channels
1335  if(addOkFlag)
1336  {
1337  ok=0;
1338 
1339  if ( hMassPi0[0]->GetEntries() < 1000 ) ok = 4;
1340  else
1341  {
1342  Float_t intLowMass = hMassPi0[0]->Integral(hMassPi0[0]->FindBin(0.04),hMassPi0[0]->FindBin(0.08));
1343  Float_t intPi0Mass = hMassPi0[0]->Integral(hMassPi0[0]->FindBin(0.10),hMassPi0[0]->FindBin(0.14));
1344  Float_t intHigMass = hMassPi0[0]->Integral(hMassPi0[0]->FindBin(0.18),hMassPi0[0]->FindBin(0.22));
1345 
1346  Float_t ratL = 0;
1347  if ( intPi0Mass > 0 ) ratL = intLowMass / intPi0Mass;
1348  Float_t ratH = 0;
1349  if ( intPi0Mass > 0 ) ratH = intHigMass / intPi0Mass;
1350  //printf("ratio InvMass Low %f High %f\n",ratL, ratH);
1351 
1352  // Just guessing
1353  if (ratL > 0.80 || ratH > 0.80 ) ok = 1;
1354  else if(ratL > 0.25 || ratH > 0.25 ) ok = 2;
1355 
1356 // hMassPi0[0]->Fit("gaus","RL","",0.11,0.14);
1357 // Float_t massFit = hMassPi0[0]->GetFunction("gaus")->GetParameter(1);
1358 // printf("mass fit = %f\n",massFit);
1359 //
1360 // if(massFit < 0.12 || massFit > 0.16) ok = 1;
1361  }
1362 
1363  text[ok]->Draw();
1364  }
1365 
1366  if(hMass[1]) // PbPb
1367  {
1368  hMassPi0[0]->SetMaximum(maxPi0*1.2);
1369  hMassPi0[5]->Draw("Hsame");
1370  hMassPi0[4]->Draw("Hsame");
1371  hMassPi0[3]->Draw("Hsame");
1372  hMassPi0[2]->Draw("Hsame");
1373  hMassPi0[1]->Draw("Hsame");
1374  hMassPi0[0]->Draw("Hsame");
1375  //hMass[6]->Draw("Hsame");
1376  //hMass[7]->Draw("same");
1377  //hMass[8]->Draw("same");
1378  //hMass[9]->Draw("same");
1379 
1380  TLegend l(0.12,0.6,0.4,0.85);
1381  l.SetTextSize(0.04);
1382  l.AddEntry(hMassPi0[0],"0-10%","P");
1383  l.AddEntry(hMassPi0[1],"10-20%","P");
1384  l.AddEntry(hMassPi0[2],"20-30%","P");
1385  l.AddEntry(hMassPi0[3],"30-40%","P");
1386  l.AddEntry(hMassPi0[4],"40-70%","P");
1387  l.AddEntry(hMassPi0[5],"50-60%","P");
1388  l.SetBorderSize(0);
1389  l.SetFillColor(0);
1390  l.Draw();
1391  }
1392 
1393  TLine l1(0.04,1,0.24,1);
1394  l1.Draw("same");
1395 
1396  // Pi0 invariant mass per EMCal super module
1397  cpi0->cd(3);
1398 
1399  TH1F* hSM [20];
1400  TH1F* hMixSM[20];
1401 
1402  //Init to 0
1403  for(Int_t ism = 0; ism < 20; ism++)
1404  {
1405  hSM [ism] = 0;
1406  hMixSM[ism] = 0;
1407  }
1408 
1409  binmin = hX->FindBin(4); // Project histo from 3 GeV pairs
1410  binmax = hX->FindBin(20); // Project histo up to 20 GeV pairs
1411  Float_t maxSM = 0;
1412 
1413  Int_t first = 0;
1414  if(icalo==1) first = 12;
1415 
1416  for(Int_t ism = first; ism < 20; ism++)
1417  {
1418  TH2F* hTmpSM = (TH2F*) GetHisto(Form("AnaPi0_Calo%d_hReMod_%d",icalo,ism));
1419  if(!hTmpSM) continue;
1420 
1421  hSM[ism] = (TH1F*) hTmpSM->ProjectionY(Form("%s_hMassSM%d",histoTag.Data(),ism),binmin,binmax);
1422  hSM[ism]->Sumw2();
1423  hSM[ism]->SetMarkerStyle(26);
1424  hSM[ism]->Rebin(2);
1425  //hSM[ism]->Scale(1./hSM[ism]->Integral(0,10000));
1426  hSM[ism]->SetMarkerColor(color[ism-first]);
1427  hSM[ism]->SetLineColor(color[ism-first]);
1428  //hSM[ism]->SetMarkerSize(0.5);
1429 
1430  TH2F* hTmpMixSM = (TH2F*) GetHisto(Form("AnaPi0_hMiMod_%d",ism));
1431  if(hTmpMixSM)
1432  {
1433  hMixSM[ism] = (TH1F*) hTmpMixSM->ProjectionY(Form("%s_hMassMixSM%d",histoTag.Data(),ism),binmin,binmax);
1434  hMixSM[ism]->Sumw2();
1435  hMixSM[ism]->Rebin(2);
1436  hSM[ism]->Divide(hMixSM[ism]);
1437  hSM[ism]->Fit("pol0","Q","",0.25,0.35);
1438  Float_t scale = 1;
1439  if(hSM[ism]->GetFunction("pol0")) scale = hSM[ism]->GetFunction("pol0")->GetParameter(0);
1440  //printf("Scale factor %f for cen %d\n",scale,icen);
1441  hSM[ism]->Scale(1./scale);
1442  }
1443 
1444  if(maxSM < hSM[ism]->GetMaximum()) maxSM = hSM[ism]->GetMaximum();
1445  }
1446 
1447  hSM[first]->SetTitle("#pi^{0} peak in SM, 4 < #it{E}_{pair}< 10 GeV");
1448  hSM[first]->SetTitleOffset(1.5,"Y");
1449  hSM[first]->SetAxisRange(0.04,0.24);
1450  hSM[first]->SetMaximum(maxSM*1.2);
1451  hSM[first]->SetMinimum(0.8);
1452  hSM[first]->SetYTitle("Real / Mixed");
1453 
1454  hSM[first]->Draw("H");
1455 
1456  // ok message
1457  // Check per SM that at low or high mass with respect the pi0 peak there is not
1458  // too much activity with respect the peak, sign of bad channels
1459  Bool_t okSM[20] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0};
1460  Bool_t okcheck = 0;
1461 
1462  if(addOkFlag)
1463  {
1464  if ( hMassPi0[0]->GetEntries() < 1000 ) text[4]->Draw();
1465  else
1466  {
1467  for(Int_t ism = first+1; ism < 20; ism++)
1468  {
1469  if(!hSM[ism]) continue;
1470 
1471  Float_t intLowMass = hSM[ism]->Integral(hSM[ism]->FindBin(0.04),hSM[ism]->FindBin(0.08));
1472  Float_t intPi0Mass = hSM[ism]->Integral(hSM[ism]->FindBin(0.10),hSM[ism]->FindBin(0.14));
1473  Float_t intHigMass = hSM[ism]->Integral(hSM[ism]->FindBin(0.18),hSM[ism]->FindBin(0.22));
1474 
1475  Float_t ratL = 0;
1476  if ( intPi0Mass > 0 ) ratL = intLowMass / intPi0Mass;
1477  Float_t ratH = 0;
1478  if ( intPi0Mass > 0 ) ratH = intHigMass / intPi0Mass;
1479  //printf("ratio InvMass Low %f High %f, ism %d\n",ratL, ratH,ism);
1480 
1481  // Just guessing
1482  if (ratL > 0.80 || ratH > 0.80 ) { okSM[ism] = 1; okcheck=1;}
1483  else if(ratL > 0.2 || ratH > 0.2 ) { okSM[ism] = 2; okcheck=1;}
1484  }
1485  }
1486  //printf("ok check %d\n",okcheck);
1487  if ( !okcheck ) {text[ok]->Draw();}
1488  }
1489 
1490  TLegend lsm(0.12,0.5,0.35,0.85);
1491  lsm.SetTextSize(0.04);
1492  if( addOkFlag && okcheck && okSM[first]) lsm.AddEntry(hSM[first],Form("#color[6]{Mod %d; not ok?}",first),"P");
1493  else lsm.AddEntry(hSM[first],Form("Mod %d",first),"P");
1494 
1495  for(Int_t ism = first+1; ism < 20; ism++)
1496  {
1497  if(!hSM[ism]) continue;
1498 
1499  hSM[ism]->Draw("Hsame");
1500 
1501  if( addOkFlag && okcheck && okSM[ism]) lsm.AddEntry(hSM[ism],Form("#color[6]{Mod %d; not ok?}",ism),"P");
1502  else lsm.AddEntry(hSM[ism],Form("Mod %d",ism),"P");
1503  }
1504 
1505  lsm.SetBorderSize(0);
1506  lsm.SetFillColor(0);
1507  lsm.Draw();
1508 
1509  l1.Draw("same");
1510 
1511  // Pi0 Invariant mass projection, in PbPb 6 centrality bins from 0 to 50%, all in pp
1512  cpi0->cd(4);
1513 
1514  //gPad->SetLogy();
1515  //gPad->SetGridy();
1516  hMassEta[0]->SetMinimum(minEta);
1517  hMassEta[0]->SetTitleOffset(1.5,"Y");
1518  hMassEta[0]->SetYTitle("Real / Mixed");
1519  hMassEta[0]->SetTitle("#eta peak, 2 <#it{E}_{pair}< 10 GeV");
1520  if (histoTag.Contains("L0")) hMassEta[0]->SetTitle("#eta peak, 5 < #it{E}_{pair}< 10 GeV");
1521  else if(histoTag.Contains("L2")) hMassEta[0]->SetTitle("#eta peak, 8 < #it{E}_{pair}< 12 GeV");
1522  else if(histoTag.Contains("L1")) hMassEta[0]->SetTitle("#eta peak, 10 < #it{E}_{pair}< 15 GeV");
1523 
1524  hMassEta[0]->Draw("H");
1525 
1526  // ok message
1527  // Check that at low or high mass with respect the eta peak there is not
1528  // too much activity with respect the peak, sign of bad channels
1529  if(addOkFlag)
1530  {
1531  ok=0;
1532 
1533  if ( hMassEta[0]->GetEntries() < 1000 ) ok = 4;
1534  else
1535  {
1536  Float_t intLowMass = hMassEta[0]->Integral(hMassEta[0]->FindBin(0.40),hMassEta[0]->FindBin(0.47));
1537  Float_t intEtaMass = hMassEta[0]->Integral(hMassEta[0]->FindBin(0.48),hMassEta[0]->FindBin(0.55));
1538  Float_t intHigMass = hMassEta[0]->Integral(hMassEta[0]->FindBin(0.60),hMassEta[0]->FindBin(0.67));
1539 
1540  Float_t ratL = 0;
1541  if ( intEtaMass > 0 ) ratL = intLowMass / intEtaMass;
1542  Float_t ratH = 0;
1543  if ( intEtaMass > 0 ) ratH = intHigMass / intEtaMass;
1544  //printf("ratio InvMass Low %f High %f\n",ratL, ratH);
1545 
1546  // Just guessing
1547  if (ratL > 1.50 || ratH > 1.50 ) ok = 1;
1548  else if(ratL > 1.10 || ratH > 1.10 ) ok = 2;
1549 
1550 // hMassEta[0]->Fit("gaus","RL","",0.48,0.55);
1551 // Float_t massFit = hMassEta[0]->GetFunction("gaus")->GetParameter(1);
1552 // printf("mass fit = %f\n",massFit);
1553 //
1554 // if(massFit < 0.48 || massFit > 0.55) ok = 1;
1555  }
1556 
1557  text[ok]->Draw();
1558  }
1559 
1560  if(hMass[1]) // PbPb
1561  {
1562  hMassEta[0]->SetMaximum(maxEta*1.2);
1563  hMassEta[5]->Draw("Hsame");
1564  hMassEta[4]->Draw("Hsame");
1565  hMassEta[3]->Draw("Hsame");
1566  hMassEta[2]->Draw("Hsame");
1567  hMassEta[1]->Draw("Hsame");
1568  hMassEta[0]->Draw("Hsame");
1569 
1570  TLegend l2(0.12,0.6,0.4,0.85);
1571  l2.SetTextSize(0.04);
1572  l2.AddEntry(hMassEta[0],"0-10%","P");
1573  l2.AddEntry(hMassEta[1],"10-20%","P");
1574  l2.AddEntry(hMassEta[2],"20-30%","P");
1575  l2.AddEntry(hMassEta[3],"30-40%","P");
1576  l2.AddEntry(hMassEta[4],"40-70%","P");
1577  l2.AddEntry(hMassEta[5],"50-60%","P");
1578  l2.SetBorderSize(0);
1579  l2.SetFillColor(0);
1580  l2.Draw();
1581  }
1582 
1583  cpi0->Print(Form("%s_Pi0Histo.%s",histoTag.Data(),format.Data()));
1584 
1585  // cleanup or save
1586  //
1587  if(exportToFile!=1)
1588  {
1589  delete h2DMass;
1590  delete hX;
1591 
1592  for(Int_t icen=0; icen<10; icen++ )
1593  {
1594  if ( hMass [icen] ) delete hMass [icen];
1595  if ( hMix [icen] ) delete hMix [icen];
1596  if ( hMassPi0[icen] ) delete hMassPi0[icen];
1597  if ( hMassEta[icen] ) delete hMassEta[icen];
1598  }
1599 
1600  for(Int_t ism = first; ism < 20; ism++)
1601  {
1602  if ( hSM [ism] ) delete hSM [ism];
1603  if ( hMixSM[ism] ) delete hMixSM[ism];
1604  }
1605 
1606  delete cpi0;
1607  }
1608  else
1609  {
1610  SaveHisto(h2DMass,kFALSE);
1611 
1612  for(Int_t icen=0; icen<10; icen++ )
1613  {
1614  SaveHisto(hMass [icen],kFALSE);
1615  SaveHisto(hMix [icen],kFALSE);
1616  SaveHisto(hMassPi0[icen],kFALSE);
1617  SaveHisto(hMassEta[icen],kFALSE);
1618  }
1619 
1620  for(Int_t ism = first; ism < 20; ism++)
1621  {
1622  SaveHisto(hSM [ism],kFALSE);
1623  SaveHisto(hMixSM[ism],kFALSE);
1624  }
1625 
1626  SaveCanvas(cpi0);
1627  }
1628 }
1629 
1637 //__________________________________________________
1638 void IsolQA(Int_t icalo)
1639 {
1640  Int_t ok = 0;
1641 
1642  TCanvas * cIsolation = new TCanvas(Form("%s_IsolationHisto" ,histoTag.Data()),
1643  Form("Isolation cone for %s",histoTag.Data()),
1644  1000,1000);
1645  cIsolation->Divide(2,2);
1646 
1647  Float_t minClusterE = 5;
1648  if ( histoTag.Contains("L0") ) minClusterE = 5;
1649  else if ( histoTag.Contains("L2") ) minClusterE = 10;
1650  else if ( histoTag.Contains("L1") ) minClusterE = 12;
1651 
1652  TH1F * hIsolated = (TH1F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hPt" ,icalo));
1653  TH1F * hNotIsolated = (TH1F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hPtNoIso",icalo));
1654 
1655  if(!hIsolated) return;
1656 
1657  Int_t minClusterEBin = hIsolated->FindBin(minClusterE);
1658  Float_t nTrig = hIsolated->Integral(minClusterEBin,100000)+hNotIsolated->Integral(minClusterE,100000);
1659 
1660  if ( nTrig <=0 ) return ;
1661 
1662  //
1663  // Candidates Pt
1664  //
1665  cIsolation->cd(1);
1666  gPad->SetLogy();
1667  gPad->SetLogx();
1668 
1669  hIsolated ->Sumw2();
1670  hNotIsolated->Sumw2();
1671  hIsolated ->SetMarkerColor(4);
1672  hNotIsolated->SetMarkerColor(2);
1673  hIsolated ->SetLineColor (4);
1674  hNotIsolated->SetLineColor (2);
1675  hIsolated ->SetMarkerStyle(24);
1676  hNotIsolated->SetMarkerStyle(20);
1677 
1678  hNotIsolated->SetTitle("(non) isolated cluster spectra, #it{R}=0.4, #Sigma #it{p}_{T}<2 GeV/#it{c}");
1679  hNotIsolated->SetYTitle("Entries");
1680 
1681  hNotIsolated->SetAxisRange(5,100,"X");
1682 
1683  hNotIsolated->Draw();
1684  hIsolated ->Draw("same");
1685 
1686  TLegend lI(0.5,0.7,0.88,0.88);
1687  lI.SetTextSize(0.04);
1688  lI.SetBorderSize(0);
1689  lI.SetFillColor(0);
1690  lI.AddEntry(hIsolated ,"Isolated candidates","P");
1691  lI.AddEntry(hNotIsolated,"NOT Isolated candidates","P");
1692  lI.Draw("same");
1693 
1694  // ok message
1695  // Check for spikes on the candidate cluster (not) isolated spectrum, sign of bad channels.
1696  // Also the isolated spectrum yield should be lower than the non isolated.
1697  if(addOkFlag)
1698  {
1699  Float_t minClusterE = 5;
1700  if ( histoTag.Contains("L2") ) minClusterE = 10;
1701  else if ( histoTag.Contains("L1") ) minClusterE = 14;
1702 
1703  ok=0;
1704 
1705  if ( hNotIsolated->GetEntries() < 1000 ) ok = 4;
1706  else
1707  {
1708  for(Int_t ibin = 1; ibin < hNotIsolated->GetNbinsX(); ibin++)
1709  {
1710  if ( hNotIsolated->GetBinCenter(ibin) < minClusterE ) continue;
1711 
1712  if ( hNotIsolated->GetBinContent(ibin) < 100 ) continue;
1713 
1714  if ( hNotIsolated->GetBinContent(ibin)*1.4 < hNotIsolated->GetBinContent(ibin+1) ||
1715  hIsolated ->GetBinContent(ibin)*1.4 < hIsolated ->GetBinContent(ibin+1) ||
1716  hNotIsolated->GetBinContent(ibin) < hIsolated ->GetBinContent(ibin)
1717  )
1718  {
1719  ok=1;
1720 // printf("ibin %d, E %2.2f, Iso: %2.1f (+1 bin) %2.1f; Not Iso %2.1f (+1 bin) %2.1f; \n",
1721 // ibin, hIsolated->GetBinCenter(ibin),
1722 // hIsolated->GetBinContent(ibin), hIsolated->GetBinContent(ibin+1),
1723 // hNotIsolated->GetBinContent(ibin), hNotIsolated->GetBinContent(ibin+1));
1724  break;
1725  }
1726  }
1727  }
1728  text[ok]->Draw();
1729  }
1730 
1731  //
1732  // Pt in cone
1733  //
1734  cIsolation->cd(2);
1735  gPad->SetLogy();
1736 
1737  TLegend l(0.55,0.55,0.88,0.88);
1738  l.SetTextSize(0.04);
1739  l.SetBorderSize(0);
1740  l.SetFillColor(0);
1741 
1742  TH2F* h2PtInCone = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hPtInCone" ,icalo));
1743  TH2F* h2PtInConeCluster = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hPtClusterInCone" ,icalo));
1744  TH2F* h2PtInConeTrack = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hPtTrackInCone" ,icalo));
1745  TH2F* h2PtInConeTrackPerp = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hPtInPerpCone" ,icalo));
1746  TH2F* h2PtInEtaBandTrack = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hEtaBandTrackPt" ,icalo));
1747  TH2F* h2PtInEtaBandCluster = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hEtaBandClusterPt",icalo));
1748 
1749  TH1F* hPtInCone = (TH1F*) h2PtInCone ->ProjectionY(Form("%s_hPtInCone_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1750  TH1F* hPtInConeCluster = (TH1F*) h2PtInConeCluster ->ProjectionY(Form("%s_hPtInConeCluster_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1751  TH1F* hPtInConeTrack = (TH1F*) h2PtInConeTrack ->ProjectionY(Form("%s_hPtInConeTrack_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1752  TH1F* hPtInConeTrackPerp = (TH1F*) h2PtInConeTrackPerp ->ProjectionY(Form("%s_hPtInConePerp_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1753  TH1F* hPtInEtaBandTrack = (TH1F*) h2PtInEtaBandTrack ->ProjectionY(Form("%s_hPtInConeEtaBandTrack_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1754  TH1F* hPtInEtaBandCluster = (TH1F*) h2PtInEtaBandCluster->ProjectionY(Form("%s_hPtInConeEtaBandCluster_TrigEnMin%2.0fGeV",histoTag.Data(),minClusterE),minClusterEBin,10000);
1755 
1756  hPtInCone ->Sumw2();
1757  hPtInConeCluster ->Sumw2();
1758  hPtInConeTrack ->Sumw2();
1759  hPtInConeTrackPerp ->Sumw2();
1760  hPtInEtaBandCluster->Sumw2();
1761  hPtInEtaBandTrack ->Sumw2();
1762 
1763  Int_t rb = 1;
1764  hPtInCone ->Rebin(rb);
1765  hPtInConeCluster ->Rebin(rb);
1766  hPtInConeTrack ->Rebin(rb);
1767  hPtInConeTrackPerp ->Rebin(rb);
1768  hPtInEtaBandCluster->Rebin(rb);
1769  hPtInEtaBandTrack ->Rebin(rb);
1770 
1771  hPtInCone ->Scale(1./nTrig);
1772  hPtInConeCluster ->Scale(1./nTrig);
1773  hPtInConeTrack ->Scale(1./nTrig);
1774  hPtInConeTrackPerp ->Scale(1./nTrig);
1775  hPtInEtaBandCluster->Scale(1./nTrig);
1776  hPtInEtaBandTrack ->Scale(1./nTrig);
1777 
1778  hPtInCone ->SetAxisRange(0,20);
1779  hPtInConeCluster ->SetAxisRange(0,20);
1780  hPtInConeTrack ->SetAxisRange(0,20);
1781  hPtInConeTrackPerp ->SetAxisRange(0,20);
1782  hPtInEtaBandCluster->SetAxisRange(0,20);
1783  hPtInEtaBandTrack ->SetAxisRange(0,20);
1784 
1785  hPtInCone ->SetMarkerStyle(24);
1786  hPtInConeCluster ->SetMarkerStyle(20);
1787  hPtInConeTrack ->SetMarkerStyle(20);
1788  hPtInConeTrackPerp ->SetMarkerStyle(27);
1789  hPtInEtaBandCluster->SetMarkerStyle(21);
1790  hPtInEtaBandTrack ->SetMarkerStyle(21);
1791 
1792  hPtInCone ->SetMarkerColor(1);
1793  hPtInConeCluster ->SetMarkerColor(2);
1794  hPtInConeTrack ->SetMarkerColor(4);
1795  hPtInConeTrackPerp ->SetMarkerColor(4);
1796  hPtInEtaBandCluster->SetMarkerColor(2);
1797  hPtInEtaBandTrack ->SetMarkerColor(4);
1798 
1799  hPtInCone ->SetLineColor(1);
1800  hPtInConeCluster ->SetLineColor(2);
1801  hPtInConeTrack ->SetLineColor(4);
1802  hPtInConeTrackPerp ->SetLineColor(4);
1803  hPtInEtaBandCluster->SetLineColor(2);
1804  hPtInEtaBandTrack ->SetLineColor(4);
1805 
1806  hPtInCone->SetTitleOffset(1.5,"Y");
1807  hPtInCone->SetYTitle("Entries / #it{N}_{candidates}");
1808  hPtInCone->SetTitle(Form("Track/cluster spectra in cone p_{T,cand}>%2.0f GeV/#it{c}, #it{R}=0.4",minClusterE));
1809 
1810  Float_t max = hPtInCone->GetMaximum();
1811  if(max < hPtInConeTrack ->GetMaximum()) max = hPtInConeTrack ->GetMaximum();
1812  if(max < hPtInConeCluster ->GetMaximum()) max = hPtInConeCluster ->GetMaximum();
1813  if(max < hPtInConeTrackPerp ->GetMaximum()) max = hPtInConeTrackPerp ->GetMaximum();
1814  if(max < hPtInEtaBandCluster->GetMaximum()) max = hPtInEtaBandCluster->GetMaximum();
1815  if(max < hPtInEtaBandTrack ->GetMaximum()) max = hPtInEtaBandTrack ->GetMaximum();
1816  hPtInCone->SetMaximum(max*2);
1817  hPtInCone->SetMinimum(1e-4);
1818 
1819  hPtInCone ->Draw("");
1820  hPtInConeCluster ->Draw("same");
1821  hPtInConeTrack ->Draw("same");
1822  hPtInConeTrackPerp ->Draw("same");
1823  hPtInEtaBandCluster->Draw("same");
1824  hPtInEtaBandTrack ->Draw("same");
1825 
1826  l.AddEntry(hPtInCone ,"Tracks+Clusters","P");
1827  l.AddEntry(hPtInConeCluster ,"Clusters inside cone","P");
1828  l.AddEntry(hPtInConeTrack ,"Tracks inside cone","P");
1829  l.AddEntry(hPtInConeTrackPerp ,"Tracks inside #perp cones","P");
1830  l.AddEntry(hPtInEtaBandTrack ,"Tracks #eta band","P");
1831  l.AddEntry(hPtInEtaBandCluster,"Clusters #eta band","P");
1832 
1833  l.Draw("same");
1834 
1835  // ok message
1836  if(addOkFlag)
1837  {
1838  ok=0;
1839 
1840  if ( nTrig < 10000 ) ok = 4;
1841  else
1842  {
1843  ok=3;
1844 
1845 // Float_t minPt = 0.5;
1846 //
1847 // for(Int_t ibin = 1; ibin < hPtInCone->GetNbinsX(); ibin++)
1848 // {
1849 // if ( hPtInCone->GetBinCenter(ibin) < minPt ) continue;
1850 //
1851 // if ( hPtInCone->GetBinContent(ibin) < 1e-4 ) continue;
1852 //
1853 // if (
1854 // hPtInCone ->GetBinContent(ibin)*2 < hPtInCone ->GetBinContent(ibin+1) ||
1855 // hPtInConeCluster ->GetBinContent(ibin) > hPtInConeTrack ->GetBinContent(ibin)
1856 // )
1857 // {
1858 // ok=1;
1859 //
1865 // break;
1866 // }
1867 // }
1868  }
1869  text[ok]->Draw();
1870  }
1871 
1872  //
1873  // Sum Pt in cone
1874  //
1875  cIsolation->cd(3);
1876  gPad->SetLogy();
1877  gPad->SetLogx();
1878 
1879  TLegend l2(0.55,0.55,0.88,0.88);
1880  l2.SetTextSize(0.04);
1881  l2.SetBorderSize(0);
1882  l2.SetFillColor(0);
1883 
1884  TH2F* h2SumPtCone = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hConePtSum" ,icalo));
1885  TH2F* h2SumPtConeCluster = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hConePtSumCluster" ,icalo));
1886  TH2F* h2SumPtConeTrack = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hConePtSumTrack" ,icalo));
1887  TH2F* h2SumPtConeTrackPerp = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hPerpConePtSum" ,icalo));
1888  TH2F* h2SumPtEtaBandTrack = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hConeSumPtEtaUENormTrack" ,icalo));
1889  TH2F* h2SumPtEtaBandCluster = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hConeSumPtEtaUENormCluster",icalo));
1890  TH2F* h2SumPtConeSubTrack = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hConeSumPtEtaUESubTrack" ,icalo));
1891  TH2F* h2SumPtConeSubCluster = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hConeSumPtEtaUESubCluster" ,icalo));
1892  TH2F* h2SumPtConeSub = (TH2F*) GetHisto(Form("AnaIsolPhoton_Calo%d_hConeSumPtEtaUESub" ,icalo));
1893 
1894  TH1F* hSumPtCone = (TH1F*) h2SumPtCone ->ProjectionY(Form("%s_hSumPtCone_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1895  TH1F* hSumPtConeCluster = (TH1F*) h2SumPtConeCluster ->ProjectionY(Form("%s_hSumPtConeCluster_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1896  TH1F* hSumPtConeTrack = (TH1F*) h2SumPtConeTrack ->ProjectionY(Form("%s_hSumPtConeTrack_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1897  TH1F* hSumPtConeTrackPerp = (TH1F*) h2SumPtConeTrackPerp ->ProjectionY(Form("%s_hSumPtConePerp_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1898  TH1F* hSumPtEtaBandTrack = (TH1F*) h2SumPtEtaBandTrack ->ProjectionY(Form("%s_hSumPtConeEtaBandTrack_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1899  TH1F* hSumPtEtaBandCluster = (TH1F*) h2SumPtEtaBandCluster->ProjectionY(Form("%s_hSumPtConeEtaBandCluster_TrigEnMin%2.0fGeV",histoTag.Data(),minClusterE),minClusterEBin,10000);
1900  TH1F* hSumPtConeSub = (TH1F*) h2SumPtConeSub ->ProjectionY(Form("%s_hSumPtConeSub_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1901  TH1F* hSumPtConeSubCluster = (TH1F*) h2SumPtConeSubCluster->ProjectionY(Form("%s_hSumPtConeSubCluster_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1902  TH1F* hSumPtConeSubTrack = (TH1F*) h2SumPtConeSubTrack ->ProjectionY(Form("%s_hSumPtConeSubTrack_TrigEnMin%2.0fGeV" ,histoTag.Data(),minClusterE),minClusterEBin,10000);
1903 
1904  hSumPtCone ->Sumw2();
1905  hSumPtConeCluster ->Sumw2();
1906  hSumPtConeTrack ->Sumw2();
1907  hSumPtConeSub ->Sumw2();
1908  hSumPtConeSubCluster->Sumw2();
1909  hSumPtConeSubTrack ->Sumw2();
1910  hSumPtConeTrackPerp ->Sumw2();
1911  hSumPtEtaBandCluster->Sumw2();
1912  hSumPtEtaBandTrack ->Sumw2();
1913 
1914  rb = 1;
1915  hSumPtCone ->Rebin(rb);
1916  hSumPtConeCluster ->Rebin(rb);
1917  hSumPtConeTrack ->Rebin(rb);
1918  hSumPtConeSub ->Rebin(rb);
1919  hSumPtConeSubCluster->Rebin(rb);
1920  hSumPtConeSubTrack ->Rebin(rb);
1921  hSumPtConeTrackPerp ->Rebin(rb);
1922  hSumPtEtaBandCluster->Rebin(rb);
1923  hSumPtEtaBandTrack ->Rebin(rb);
1924 
1925  hSumPtCone ->Scale(1./nTrig);
1926  hSumPtConeCluster ->Scale(1./nTrig);
1927  hSumPtConeTrack ->Scale(1./nTrig);
1928  hSumPtConeSub ->Scale(1./nTrig);
1929  hSumPtConeSubCluster->Scale(1./nTrig);
1930  hSumPtConeSubTrack ->Scale(1./nTrig);
1931  hSumPtConeTrackPerp ->Scale(1./nTrig);
1932  hSumPtEtaBandCluster->Scale(1./nTrig);
1933  hSumPtEtaBandTrack ->Scale(1./nTrig);
1934 
1935  hSumPtCone ->SetAxisRange(0,500);
1936  hSumPtConeCluster ->SetAxisRange(0,500);
1937  hSumPtConeTrack ->SetAxisRange(0,500);
1938  hSumPtConeSub ->SetAxisRange(-5,500);
1939  hSumPtConeSubCluster->SetAxisRange(-5,500);
1940  hSumPtConeSubTrack ->SetAxisRange(-5,500);
1941  hSumPtConeTrackPerp ->SetAxisRange(0,500);
1942  hSumPtEtaBandCluster->SetAxisRange(0,500);
1943  hSumPtEtaBandTrack ->SetAxisRange(0,500);
1944 
1945  hSumPtCone ->SetMarkerStyle(24);
1946  hSumPtConeCluster ->SetMarkerStyle(20);
1947  hSumPtConeTrack ->SetMarkerStyle(20);
1948  hSumPtConeSub ->SetMarkerStyle(25);
1949  hSumPtConeSubCluster->SetMarkerStyle(25);
1950  hSumPtConeSubTrack ->SetMarkerStyle(25);
1951  hSumPtConeTrackPerp ->SetMarkerStyle(27);
1952  hSumPtEtaBandCluster->SetMarkerStyle(21);
1953  hSumPtEtaBandTrack ->SetMarkerStyle(21);
1954 
1955  hSumPtCone ->SetMarkerColor(1);
1956  hSumPtConeCluster ->SetMarkerColor(2);
1957  hSumPtConeTrack ->SetMarkerColor(4);
1958  hSumPtConeSub ->SetMarkerColor(1);
1959  hSumPtConeSubCluster->SetMarkerColor(2);
1960  hSumPtConeSubTrack ->SetMarkerColor(4);
1961  hSumPtConeTrackPerp ->SetMarkerColor(4);
1962  hSumPtEtaBandCluster->SetMarkerColor(2);
1963  hSumPtEtaBandTrack ->SetMarkerColor(4);
1964 
1965  hSumPtCone ->SetLineColor(1);
1966  hSumPtConeCluster ->SetLineColor(2);
1967  hSumPtConeTrack ->SetLineColor(4);
1968  hSumPtConeSub ->SetLineColor(1);
1969  hSumPtConeSubCluster->SetLineColor(2);
1970  hSumPtConeSubTrack ->SetLineColor(4);
1971  hSumPtConeTrackPerp ->SetLineColor(4);
1972  hSumPtEtaBandCluster->SetLineColor(2);
1973  hSumPtEtaBandTrack ->SetLineColor(4);
1974 
1975  hSumPtCone->SetTitleOffset(1.5,"Y");
1976  hSumPtCone->SetYTitle("Entries / #it{N}_{candidates}");
1977  hSumPtCone->SetTitle(Form("Track/cluster #Sigma #it{p}_{T}, p_{T,cand}>%2.0f GeV/#it{c}, #it{R}=0.4, ",minClusterE));
1978 
1979  max = hSumPtCone->GetMaximum();
1980  if(max < hSumPtConeTrack ->GetMaximum()) max = hSumPtConeTrack ->GetMaximum();
1981  if(max < hSumPtConeCluster ->GetMaximum()) max = hSumPtConeCluster ->GetMaximum();
1982  if(max < hSumPtConeSub ->GetMaximum()) max = hSumPtConeSub ->GetMaximum();
1983  if(max < hSumPtConeSubTrack ->GetMaximum()) max = hSumPtConeSubTrack ->GetMaximum();
1984  if(max < hSumPtConeSubCluster->GetMaximum()) max = hSumPtConeSubCluster->GetMaximum();
1985  if(max < hSumPtConeTrackPerp ->GetMaximum()) max = hSumPtConeTrackPerp ->GetMaximum();
1986  if(max < hSumPtEtaBandCluster->GetMaximum()) max = hSumPtEtaBandCluster->GetMaximum();
1987  if(max < hSumPtEtaBandTrack ->GetMaximum()) max = hSumPtEtaBandTrack ->GetMaximum();
1988  hSumPtCone->SetMaximum(max*2);
1989 
1990  hSumPtCone ->Draw("");
1991  hSumPtConeCluster ->Draw("same");
1992  hSumPtConeTrack ->Draw("same");
1993 // hSumPtConeSub ->Draw("same");
1994 // hSumPtConeSubCluster->Draw("same");
1995 // hSumPtConeSubTrack ->Draw("same");
1996  hSumPtConeTrackPerp ->Draw("same");
1997  hSumPtEtaBandCluster->Draw("same");
1998  hSumPtEtaBandTrack ->Draw("same");
1999 
2000  l2.AddEntry(hSumPtCone ,"Tracks+Clusters","P");
2001  l2.AddEntry(hSumPtConeCluster ,"Clusters inside cone","P");
2002  l2.AddEntry(hSumPtConeTrack ,"Tracks inside cone","P");
2003  l2.AddEntry(hSumPtConeTrackPerp ,"Tracks inside #perp cones","P");
2004  l2.AddEntry(hSumPtEtaBandTrack ,"Tracks #eta band","P");
2005  l2.AddEntry(hSumPtEtaBandCluster,"Clusters #eta band","P");
2006 // l2.AddEntry(hSumPtConeSub ,"Tracks+Clusters-#eta band","P");
2007 // l2.AddEntry(hSumPtConeSubCluster,"Clusters inside cone-#eta band","P");
2008 // l2.AddEntry(hSumPtConeSubTrack ,"Tracks inside cone-#eta band","P");
2009 
2010  l2.Draw("same");
2011 
2012  // ok message
2013  if(addOkFlag)
2014  {
2015  ok=0;
2016 
2017  if ( nTrig < 10000 ) ok = 4;
2018  else
2019  {
2020  ok=3;
2021 // Float_t minPt = 3;
2022 // Float_t maxPt = 30;
2023 //
2024 // for(Int_t ibin = 1; ibin < hSumPtCone->GetNbinsX(); ibin++)
2025 // {
2026 // if ( hSumPtCone->GetBinCenter(ibin) < minPt ) continue;
2027 // if ( hSumPtCone->GetBinCenter(ibin) > maxPt ) continue;
2028 //
2029 // if ( hSumPtCone->GetBinContent(ibin) < 1e-4 ) continue;
2030 //
2031 // if (
2032 // hSumPtCone ->GetBinContent(ibin)*2 < hSumPtCone ->GetBinContent(ibin+1) ||
2033 // hSumPtConeCluster ->GetBinContent(ibin) > hSumPtConeTrack ->GetBinContent(ibin)
2034 // )
2035 // {
2036 // ok=1;
2042 //
2043 // break;
2044 // }
2045 // }
2046  }
2047  text[ok]->Draw();
2048  }
2049 
2050  //
2051  // Sum Pt in cone, UE subtracted
2052  //
2053  cIsolation->cd(4);
2054  gPad->SetLogy();
2055 
2056  TLegend l3(0.4,0.75,0.8,0.88);
2057  l3.SetTextSize(0.04);
2058  l3.SetBorderSize(0);
2059  l3.SetFillColor(0);
2060 
2061  hSumPtConeSub->SetTitle(Form("Track/Cluster #Sigma #it{p}_{T}-#Sigma #eta band, p_{T,cand}>%2.0f GeV/#it{c}, #it{R}=0.4",minClusterE));
2062  hSumPtConeSub->SetYTitle("Entries / #it{N}_{candidates}");
2063  hSumPtConeSub->SetMaximum(max*2);
2064 
2065  hSumPtConeSub ->Draw("");
2066  hSumPtConeSubCluster->Draw("same");
2067  hSumPtConeSubTrack ->Draw("same");
2068 
2069  l3.AddEntry(hSumPtConeSub ,"Tracks+Clusters-#eta band","P");
2070  l3.AddEntry(hSumPtConeSubCluster,"Clusters inside cone-#eta band","P");
2071  l3.AddEntry(hSumPtConeSubTrack ,"Tracks inside cone-#eta band","P");
2072 
2073  l3.Draw("same");
2074 
2075 
2076  // ok message
2077  if(addOkFlag)
2078  {
2079  ok=0;
2080 
2081  if ( nTrig < 10000 ) ok = 4;
2082  else
2083  {
2084  ok = 3;
2085 // Float_t minPt = 3;
2086 // Float_t maxPt = 30;
2087 //
2088 // for(Int_t ibin = 1; ibin < hSumPtCone->GetNbinsX(); ibin++)
2089 // {
2090 // if ( hSumPtConeSub->GetBinCenter(ibin) < minPt ) continue;
2091 // if ( hSumPtConeSub->GetBinCenter(ibin) > maxPt ) continue;
2092 //
2093 // if ( hSumPtConeSub->GetBinContent(ibin) < 1e-4 ) continue;
2094 //
2095 // if (
2096 // hSumPtConeSub ->GetBinContent(ibin)*2 < hSumPtConeSub ->GetBinContent(ibin+1) //||
2097 // //hSumPtConeSubCluster->GetBinContent(ibin) > hSumPtConeSubTrack ->GetBinContent(ibin)
2098 // )
2099 // {
2100 // ok=1;
2101 // break;
2102 // }
2103 // }
2104  }
2105  text[ok]->Draw();
2106  }
2107 
2108 
2109  cIsolation->Print(Form("%s_IsolationHisto.%s",histoTag.Data(),format.Data()));
2110 
2111  // cleanup or save
2112  //
2113  if(exportToFile!=1)
2114  {
2115  delete hPtInCone ;
2116  delete hPtInConeCluster ;
2117  delete hPtInConeTrack ;
2118  delete hPtInConeTrackPerp ;
2119  delete hPtInEtaBandTrack ;
2120  delete hPtInEtaBandCluster ;
2121 
2122  delete hSumPtCone ;
2123  delete hSumPtConeCluster ;
2124  delete hSumPtConeTrack ;
2125  delete hSumPtConeTrackPerp ;
2126  delete hSumPtEtaBandTrack ;
2127  delete hSumPtEtaBandCluster ;
2128 
2129  delete hSumPtConeSub ;
2130  delete hSumPtConeSubCluster ;
2131  delete hSumPtConeSubTrack ;
2132 
2133  delete cIsolation ;
2134  }
2135  else
2136  {
2137  SaveHisto(hPtInCone ,kFALSE);
2138  SaveHisto(hPtInConeCluster ,kFALSE);
2139  SaveHisto(hPtInConeTrack ,kFALSE);
2140  SaveHisto(hPtInConeTrackPerp ,kFALSE);
2141  SaveHisto(hPtInEtaBandTrack ,kFALSE);
2142  SaveHisto(hPtInEtaBandCluster ,kFALSE);
2143 
2144  SaveHisto(hSumPtCone ,kFALSE);
2145  SaveHisto(hSumPtConeCluster ,kFALSE);
2146  SaveHisto(hSumPtConeTrack ,kFALSE);
2147  SaveHisto(hSumPtConeTrackPerp ,kFALSE);
2148  SaveHisto(hSumPtEtaBandTrack ,kFALSE);
2149  SaveHisto(hSumPtEtaBandCluster,kFALSE);
2150 
2151  SaveHisto(hSumPtConeSub ,kFALSE);
2152  SaveHisto(hSumPtConeSubCluster,kFALSE);
2153  SaveHisto(hSumPtConeSubTrack ,kFALSE);
2154 
2155  SaveCanvas(cIsolation);
2156  }
2157 }
2158 
2165 //__________________________________________________
2166 void CorrelQA(Int_t icalo)
2167 {
2168  Int_t ok = 0;
2169 
2170  TCanvas * cCorrelation = new TCanvas(Form("%s_CorrelationHisto" ,histoTag.Data()),
2171  Form("Trigger cluster - associated track correlation for %s",histoTag.Data()),
2172  1000,500);
2173  cCorrelation->Divide(2,1);
2174 
2175  Float_t minClusterE = 5;
2176  if ( histoTag.Contains("L0") ) minClusterE = 8;
2177  else if ( histoTag.Contains("L2") ) minClusterE = 10;
2178  else if ( histoTag.Contains("L1") ) minClusterE = 12;
2179 
2180  Float_t assocBins[] = {0.5,2.,5.,10.,20.};
2181  Int_t nAssocBins = 4;
2182 
2183  TH1F * hTrigger = (TH1F*) GetHisto(Form("AnaPhotonHadronCorr_Calo%d_hPtTrigger",icalo));
2184 
2185  if(!hTrigger) return;
2186 
2187  Int_t minClusterEBin = hTrigger->FindBin(minClusterE);
2188  Float_t nTrig = hTrigger->Integral(minClusterEBin,100000);
2189 
2190  if ( nTrig <=0 ) return ;
2191 
2192  //Azimuthal correlation
2193  cCorrelation->cd(1);
2194  gPad->SetLogy();
2195  TH1F* hDeltaPhi[4];
2196  for(Int_t i = 0; i < 4; i++) hDeltaPhi[i] = 0;
2197 
2198  TLegend l(0.35,0.6,0.83,0.85);
2199  l.SetHeader(Form("p_{T,T} > %2.1f GeV/c",minClusterE));
2200  l.SetTextSize(0.04);
2201  l.SetBorderSize(0);
2202  l.SetFillColor(0);
2203 
2204  for(Int_t ibin = 0; ibin < nAssocBins; ibin++ )
2205  {
2206  TH2F* hDeltaPhiE =
2207  (TH2F*) GetHisto(Form("AnaPhotonHadronCorr_Calo%d_hDeltaPhiPtAssocPt%2.1f_%2.1f",icalo,assocBins[ibin],assocBins[ibin+1]));
2208  hDeltaPhi[ibin] =
2209  (TH1F*) hDeltaPhiE->ProjectionY(Form("%s_hDeltaPhi_TrackMinPt%2.1fGeV_TrigEnMin%2.0f",histoTag.Data(),assocBins[ibin],minClusterE),minClusterEBin,10000);
2210  hDeltaPhi[ibin]->Sumw2();
2211  hDeltaPhi[ibin]->Rebin(2);
2212  hDeltaPhi[ibin]->Scale(1./nTrig);
2213 
2214  hDeltaPhi[ibin]->Fit("pol0","Q","",1,2);
2215 
2216  Float_t scale = 1;
2217  if(hDeltaPhi[ibin]->GetFunction("pol0"))
2218  {
2219  scale = hDeltaPhi[ibin]->GetFunction("pol0")->GetParameter(0);
2220  hDeltaPhi[ibin]->GetFunction("pol0")->SetRange(6,7); // move from plot
2221  }
2222  hDeltaPhi[ibin]->Scale(1./scale);
2223  //printf("ibin %d, scale %f\n",ibin,scale);
2224 
2225  hDeltaPhi[ibin]->SetAxisRange(-1.6,4.7);
2226 
2227  hDeltaPhi[ibin]->SetMarkerStyle(24);
2228  hDeltaPhi[ibin]->SetMarkerColor(color[ibin]);
2229  hDeltaPhi[ibin]->SetLineColor(color[ibin]);
2230  hDeltaPhi[ibin]->SetTitleOffset(1.5,"Y");
2231  hDeltaPhi[ibin]->SetYTitle("#it{N}_{pairs} / #it{N}_{trig} / ZYAM");
2232  hDeltaPhi[ibin]->SetTitle("#gamma (#lambda_{0}^{2} < 0.4, neutral cluster) trigger");
2233 
2234  if(!addOkFlag) l.AddEntry(hDeltaPhi[ibin],Form("%2.1f<#it{p}_{T,A}< %2.0f GeV/c",assocBins[ibin],assocBins[ibin+1]),"P");
2235  }
2236 
2237  hDeltaPhi[2]->SetMaximum(hDeltaPhi[2]->GetMaximum()*10);
2238  hDeltaPhi[2]->SetMinimum(0.8);
2239 
2240  hDeltaPhi[2]->Draw("H");
2241  hDeltaPhi[1]->Draw("Hsame");
2242  hDeltaPhi[3]->Draw("Hsame");
2243  hDeltaPhi[0]->Draw("Hsame");
2244 
2245  // ok message
2246  // Check that in the near region there is more yield than in the away region and more than in the zyam region
2247  // per each of the associated pT bins
2248  if(addOkFlag)
2249  {
2250  Bool_t okcheck = 0;
2251  if ( nTrig < 1000 ) text[4]->Draw();
2252  else
2253  {
2254  for(Int_t ibin = 0; ibin < nAssocBins; ibin++ )
2255  {
2256  if(!hDeltaPhi[ibin]) continue;
2257 
2258  Float_t intNear = hDeltaPhi[ibin]->Integral(hDeltaPhi[ibin]->FindBin(-0.4),hDeltaPhi[ibin]->FindBin(0.4));
2259  Float_t intZyam = hDeltaPhi[ibin]->Integral(hDeltaPhi[ibin]->FindBin( 1.1),hDeltaPhi[ibin]->FindBin(1.9));
2260  Float_t intAway = hDeltaPhi[ibin]->Integral(hDeltaPhi[ibin]->FindBin( 2.8),hDeltaPhi[ibin]->FindBin(3.6));
2261 
2262  if(intNear > intZyam && intAway > intZyam && intNear > intAway )
2263  {
2264  l.AddEntry(hDeltaPhi[ibin],Form("%2.1f<#it{p}_{T,A}< %2.0f GeV/c",assocBins[ibin],assocBins[ibin+1]),"P");
2265  }
2266  else
2267  {
2268  //printf("bin %d, near %2.2e, zyam %2.2e, away %2.2e\n",ibin,intNear,intZyam,intAway);
2269  l.AddEntry(hDeltaPhi[ibin],Form("#color[6]{%2.1f<#it{p}_{T,A}< %2.0f GeV/c, not ok?}",assocBins[ibin],assocBins[ibin+1]),"P");
2270  okcheck=1;
2271  }
2272  }
2273  }
2274  //printf("ok check %d\n",okcheck);
2275  if ( !okcheck ) {text[ok]->Draw();}
2276  }
2277 
2278 
2279  l.Draw("same");
2280 
2281  // xE correlation
2282  cCorrelation->cd(2);
2283  gPad->SetLogy();
2284 
2285  TLegend l2(0.35,0.6,0.83,0.85);
2286  l2.SetHeader(Form("p_{T,T} > %2.0f GeV/c",minClusterE));
2287  l2.SetTextSize(0.04);
2288  l2.SetBorderSize(0);
2289  l2.SetFillColor(0);
2290 
2291  TH2F* hEXE = (TH2F*) GetHisto(Form("AnaPhotonHadronCorr_Calo%d_hXECharged" ,icalo));
2292  TH2F* hEXEUE = (TH2F*) GetHisto(Form("AnaPhotonHadronCorr_Calo%d_hXEUeCharged",icalo));
2293 
2294  TH1F* hXE = (TH1F*) hEXE->ProjectionY(Form("%s_hXE_TrigEnMin%2.0fGeV",histoTag.Data(),minClusterE),minClusterEBin,10000);
2295  hXE->Sumw2();
2296  hXE->Rebin(2);
2297  hXE->Scale(1./nTrig);
2298  hXE->SetAxisRange(0,1);
2299  hXE->SetMarkerStyle(24);
2300  hXE->SetMarkerColor(1);
2301  hXE->SetLineColor(1);
2302  hXE->SetTitleOffset(1.5,"Y");
2303  hXE->SetYTitle("#it{N}_{pairs} / #it{N}_{trig}");
2304  hXE->SetTitle("#gamma (#lambda_{0}^{2} < 0.4, neutral cluster) trigger");
2305  l2.AddEntry(hXE,"raw x_{E}","P");
2306  hXE->Draw();
2307 
2308  TH1F* hXEUE = (TH1F*) hEXEUE->ProjectionY(Form("%s_hXEUE_TrigEnMin%2.0fGeV",histoTag.Data(),minClusterE),minClusterEBin,10000);
2309  hXEUE->Sumw2();
2310  hXEUE->Rebin(2);
2311  hXEUE->Scale(1./nTrig);
2312  hXEUE->SetAxisRange(0,1);
2313  hXEUE->SetMarkerStyle(25);
2314  hXEUE->SetMarkerColor(2);
2315  hXEUE->SetLineColor(2);
2316  l2.AddEntry(hXEUE,"raw Und. Event x_{E}","P");
2317  hXEUE->Draw("same");
2318 
2319  l2.Draw("same");
2320 
2321  // ok message
2322  // Check that there are no spikes and that the UE yield is lower than the correlated yield.
2323  if(addOkFlag)
2324  {
2325  ok=0;
2326 
2327  if ( nTrig < 1000 ) ok = 4;
2328  else
2329  {
2330  //ok=3;
2331 
2332  Float_t minXE = 0.1;
2333  Float_t maxXE = 0.5;
2334 
2335  for(Int_t ibin = 1; ibin < hXE->GetNbinsX(); ibin++)
2336  {
2337  if ( hXE->GetBinCenter(ibin) < minXE ) continue;
2338  if ( hXE->GetBinCenter(ibin) > maxXE ) continue;
2339 
2340  if ( hXE->GetBinContent(ibin) < 1e-3 ) continue;
2341 
2342  if (
2343  hXE ->GetBinContent(ibin)*2 < hXE ->GetBinContent(ibin+1) ||
2344  //hXEUE->GetBinContent(ibin)*2 < hXEUE->GetBinContent(ibin+1) ||
2345  hXE ->GetBinContent(ibin) < hXEUE->GetBinContent(ibin)
2346  )
2347  {
2348  ok=1;
2349 
2350  printf("ibin %d, XE bin %2.2f, XE: %2.1e (+1 bin) %2.1e; XEUE %2.1e (+1 bin) %2.1e;\n",
2351  ibin, hXE->GetBinCenter(ibin),
2352  hXE->GetBinContent(ibin)*2, hXE->GetBinContent(ibin+1),
2353  hXEUE->GetBinContent(ibin)*2, hXEUE->GetBinContent(ibin+1));
2354  break;
2355  }
2356  }
2357  }
2358  text[ok]->Draw();
2359  }
2360 
2361 
2362  cCorrelation->Print(Form("%s_CorrelationHisto.%s",histoTag.Data(),format.Data()));
2363 
2364  // cleanup or save
2365  //
2366  if(exportToFile!=1)
2367  {
2368  for(Int_t i = 0; i < 4; i++) delete hDeltaPhi[i];
2369 
2370  delete hXE ;
2371  delete hXEUE;
2372 
2373  delete cCorrelation;
2374  }
2375  else
2376  {
2377  for(Int_t i = 0; i < 4; i++) SaveHisto(hDeltaPhi[i],kFALSE);
2378 
2379  SaveHisto(hXE ,kFALSE) ;
2380  SaveHisto(hXEUE,kFALSE);
2381 
2382  SaveCanvas(cCorrelation);
2383  }
2384 
2385 }
2386 
2395 //________________________________________________________
2396 void MCQA(Int_t icalo)
2397 {
2398  TH1F* hClusterPho = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPt_MCPhoton",icalo));
2399  TH1F* hClusterPi0 = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPt_MCPi0" ,icalo));
2400  TH1F* hClusterEta = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPt_MCEta" ,icalo));
2401  TH1F* hClusterPhoPi0 = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPt_MCPhotonPi0Decay",icalo));
2402  TH1F* hClusterPhoEta = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPt_MCPhotonEtaDecay",icalo));
2403 
2404  if(!hClusterPho) return;
2405 
2406  TH1F* hPrimPho = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPtPrim_MCPhoton" ,icalo));
2407  TH1F* hPrimPhoPi0 = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPtPrim_MCPhotonPi0Decay",icalo));
2408  TH1F* hPrimPhoEta = (TH1F*) GetHisto(Form("AnaPhoton_Calo%d_hPtPrim_MCPhotonEtaDecay",icalo));
2409  TH1F* hPrimPi0 = (TH1F*) GetHisto(Form("AnaPi0_Calo%d_hPrimPi0Pt",icalo));
2410  TH1F* hPrimEta = (TH1F*) GetHisto(Form("AnaPi0_Calo%d_hPrimEtaPt",icalo));
2411 
2412  TCanvas * cmc = new TCanvas(Form("%s_MCHisto" ,histoTag.Data()),
2413  Form("Cluster MC origin for %s",histoTag.Data()),
2414  1000,1000);
2415  cmc->Divide(2,2);
2416 
2417  cmc->cd(1);
2418  gPad->SetLogy();
2419 
2420  hClusterPho->SetTitle("Cluster origin spectra, primary spectra in Calo acceptance");
2421  hClusterPho->Sumw2();
2422  hClusterPho->SetMarkerColor(1);
2423  hClusterPho->SetMarkerStyle(20);
2424  hClusterPho->SetAxisRange(0.,50.,"X");
2425  //hClusterPho->SetXTitle("E_{rec,gen} (GeV)");
2426  hClusterPho->SetXTitle("#it{E}_{rec}, #it{p}_{T,gen} (GeV)");
2427  hClusterPho->SetYTitle("Entries");
2428  hClusterPho->Draw("");
2429 
2430  hClusterPhoPi0->Sumw2();
2431  hClusterPhoPi0->SetMarkerColor(4);
2432  hClusterPhoPi0->SetMarkerStyle(20);
2433  hClusterPhoPi0->Draw("same");
2434 
2435  hClusterPhoEta->Sumw2();
2436  hClusterPhoEta->SetMarkerColor(2);
2437  hClusterPhoEta->SetMarkerStyle(20);
2438  hClusterPhoEta->Draw("same");
2439 
2440  hClusterPi0->Sumw2();
2441  hClusterPi0->SetMarkerColor(4);
2442  hClusterPi0->SetMarkerStyle(21);
2443  hClusterPi0->Draw("same");
2444 
2445  hClusterEta->Sumw2();
2446  hClusterEta->SetMarkerColor(2);
2447  hClusterEta->SetMarkerStyle(22);
2448  hClusterEta->Draw("same");
2449 
2450  hPrimPho->Sumw2();
2451  hPrimPho->SetMarkerColor(1);
2452  hPrimPho->SetMarkerStyle(24);
2453  hPrimPho->Draw("same");
2454 
2455  hPrimPhoPi0->Sumw2();
2456  hPrimPhoPi0->SetMarkerColor(4);
2457  hPrimPhoPi0->SetMarkerStyle(24);
2458  hPrimPhoPi0->Draw("same");
2459 
2460  hPrimPhoEta->Sumw2();
2461  hPrimPhoEta->SetMarkerColor(2);
2462  hPrimPhoEta->SetMarkerStyle(24);
2463  hPrimPhoEta->Draw("same");
2464 
2465  hPrimPi0->Sumw2();
2466  hPrimPi0->SetMarkerColor(4);
2467  hPrimPi0->SetMarkerStyle(25);
2468  hPrimPi0->Draw("same");
2469 
2470  hPrimEta->Sumw2();
2471  hPrimEta->SetMarkerColor(2);
2472  hPrimEta->SetMarkerStyle(26);
2473  hPrimEta->Draw("same");
2474 
2475  TLegend lR(0.5,0.5,0.7,0.89);
2476  lR.SetHeader("reco");
2477  lR.SetTextSize(0.04);
2478  lR.AddEntry(hClusterPho,"#gamma","P");
2479  lR.AddEntry(hClusterPhoPi0,"#gamma_{#pi^{0}}","P");
2480  lR.AddEntry(hClusterPhoEta,"#gamma_{#eta}","P");
2481  lR.AddEntry(hClusterPi0,"#pi^{0}","P");
2482  lR.AddEntry(hClusterEta,"#eta","P");
2483  lR.SetBorderSize(0);
2484  lR.SetFillColor(0);
2485  lR.Draw();
2486 
2487  TLegend lG(0.7,0.5,0.83,0.89);
2488  lG.SetHeader("gener");
2489  lG.SetTextSize(0.04);
2490  lG.AddEntry(hPrimPho,"#gamma","P");
2491  lG.AddEntry(hPrimPhoPi0,"#gamma_{#pi^{0}}","P");
2492  lG.AddEntry(hPrimPhoEta,"#gamma_{#eta}","P");
2493  lG.AddEntry(hPrimPi0,"#pi^{0}","P");
2494  lG.AddEntry(hPrimEta,"#eta","P");
2495  lG.SetBorderSize(0);
2496  lG.SetFillColor(0);
2497  lG.Draw();
2498 
2499  if(addOkFlag) text[3]->Draw();
2500 
2501  cmc->cd(2);
2502  gPad->SetLogy();
2503  TH1F* hRatPho = (TH1F*) hClusterPho ->Clone(Form("%s_hGenRecoRatPho" ,histoTag.Data()));
2504  TH1F* hRatPi0 = (TH1F*) hClusterPi0 ->Clone(Form("%s_hGenRecoRatPi0" ,histoTag.Data()));
2505  TH1F* hRatEta = (TH1F*) hClusterEta ->Clone(Form("%s_hGenRecoRatEta" ,histoTag.Data()));
2506  TH1F* hRatPhoPi0 = (TH1F*) hClusterPhoPi0->Clone(Form("%s_hGenRecoRatPhoPi0",histoTag.Data()));
2507  TH1F* hRatPhoEta = (TH1F*) hClusterPhoEta->Clone(Form("%s_hGenRecoRatPhoEta",histoTag.Data()));
2508 
2509  hRatPho ->Divide(hPrimPho);
2510  hRatPhoPi0->Divide(hPrimPhoPi0);
2511  hRatPhoEta->Divide(hPrimPhoEta);
2512  hRatPi0 ->Divide(hPrimPi0);
2513  hRatEta ->Divide(hPrimEta);
2514 
2515  hRatPho->SetTitle("Reconstructed cluster / Generated particle in Calo acc.");
2516  hRatPho->SetYTitle("#it{Ratio reco / gener}");
2517  hRatPho->SetXTitle("#it{E} (GeV)");
2518  hRatPho->SetMinimum(1e-3);
2519  hRatPho->SetMaximum(20);
2520  hRatPho->Draw("");
2521  hRatPhoPi0->Draw("same");
2522  hRatPhoEta->Draw("same");
2523  hRatPi0->Draw("same");
2524  hRatEta->Draw("same");
2525 
2526  TLegend l2(0.15,0.62,0.3,0.89);
2527  l2.SetTextSize(0.04);
2528  l2.AddEntry(hRatPho,"#gamma","P");
2529  l2.AddEntry(hRatPhoPi0,"#gamma_{#pi^{0}}","P");
2530  l2.AddEntry(hRatPhoEta,"#gamma_{#eta}","P");
2531  l2.AddEntry(hRatPi0,"#pi^{0}","P");
2532  l2.AddEntry(hRatEta,"#eta","P");
2533  l2.SetBorderSize(0);
2534  l2.SetFillColor(0);
2535  l2.Draw();
2536 
2537  if(addOkFlag) text[3]->Draw();
2538 
2539  cmc->cd(3);
2540  //gPad->SetLogy();
2541 
2542  TH2F* h2PrimPhoPhi = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hPhiPrim_MCPhoton",icalo));
2543  TH2F* h2PrimPi0Phi = (TH2F*) GetHisto(Form("AnaPi0_Calo%d_hPrimPi0Phi" ,icalo));
2544  TH2F* h2PrimEtaPhi = (TH2F*) GetHisto(Form("AnaPi0_Calo%d_hPrimEtaPhi" ,icalo));
2545 
2546  Int_t binMin = hPrimPho->FindBin(3);
2547 
2548  TH1F* hPrimPhoPhi = (TH1F*) h2PrimPhoPhi->ProjectionY(Form("%s_hPrimPhoPhi",histoTag.Data()),binMin,1000);
2549  TH1F* hPrimPi0Phi = (TH1F*) h2PrimPi0Phi->ProjectionY(Form("%s_hPrimPi0Phi",histoTag.Data()),binMin,1000);
2550  TH1F* hPrimEtaPhi = (TH1F*) h2PrimEtaPhi->ProjectionY(Form("%s_hPrimEtaPhi",histoTag.Data()),binMin,1000);
2551 
2552  hPrimPhoPhi->Sumw2();
2553  hPrimPi0Phi->Sumw2();
2554  hPrimEtaPhi->Sumw2();
2555 
2556  hPrimPhoPhi->Scale(1./hPrimPhoPhi->Integral(0,1000));
2557  hPrimPi0Phi->Scale(1./hPrimPi0Phi->Integral(0,1000));
2558  hPrimEtaPhi->Scale(1./hPrimEtaPhi->Integral(0,1000));
2559 
2560  Float_t maxPhi = hPrimPhoPhi->GetMaximum();
2561  if(maxPhi < hPrimPi0Phi->GetMaximum()) maxPhi = hPrimPi0Phi->GetMaximum();
2562  if(maxPhi < hPrimEtaPhi->GetMaximum()) maxPhi = hPrimEtaPhi->GetMaximum();
2563 
2564  Float_t minPhi = hPrimPhoPhi->GetMinimum();
2565  if(minPhi > hPrimPi0Phi->GetMinimum()) minPhi = hPrimPi0Phi->GetMinimum();
2566  if(minPhi > hPrimEtaPhi->GetMinimum()) minPhi = hPrimEtaPhi->GetMinimum();
2567 
2568  hPrimPi0Phi->SetMaximum(maxPhi*1.1);
2569  hPrimPi0Phi->SetMinimum(minPhi);
2570  TGaxis::SetMaxDigits(3);
2571 
2572  hPrimPi0Phi->SetYTitle("1/total entries d#it{N}/d#varphi");
2573  hPrimPi0Phi->SetTitle("Generated particles #varphi for #it{E} > 3 GeV");
2574  hPrimPi0Phi->SetTitleOffset(1.5,"Y");
2575  hPrimPi0Phi->SetMarkerColor(4);
2576  hPrimPi0Phi->SetMarkerStyle(21);
2577  hPrimPi0Phi->Draw("");
2578 
2579  hPrimPhoPhi->SetMarkerColor(1);
2580  hPrimPhoPhi->SetMarkerStyle(20);
2581  Float_t scale = TMath::RadToDeg();
2582  ScaleXaxis(hPrimPhoPhi, TMath::RadToDeg());
2583  hPrimPhoPhi->Draw("same");
2584 
2585  hPrimEtaPhi->SetMarkerColor(2);
2586  hPrimEtaPhi->SetMarkerStyle(22);
2587  hPrimEtaPhi->Draw("same");
2588 
2589  if(addOkFlag) text[3]->Draw();
2590 
2591  cmc->cd(4);
2592  //gPad->SetLogy();
2593 
2594  TH2F* h2PrimPhoEtaP = (TH2F*) GetHisto(Form("AnaPhoton_Calo%d_hYPrim_MCPhoton",icalo));
2595  TH2F* h2PrimPi0EtaP = (TH2F*) GetHisto(Form("AnaPi0_Calo%d_hPrimPi0Rapidity" ,icalo));
2596  TH2F* h2PrimEtaEtaP = (TH2F*) GetHisto(Form("AnaPi0_Calo%d_hPrimEtaRapidity" ,icalo));
2597 
2598  h2PrimPhoEtaP->Sumw2();
2599  h2PrimEtaEtaP->Sumw2();
2600  h2PrimPi0EtaP->Sumw2();
2601 
2602  binMin = hPrimPho->FindBin(3);
2603 
2604  TH1F* hPrimPhoEtaP = (TH1F*) h2PrimPhoEtaP->ProjectionY(Form("%s_hPrimPhoEtaP",histoTag.Data()),binMin,1000);
2605  TH1F* hPrimPi0EtaP = (TH1F*) h2PrimPi0EtaP->ProjectionY(Form("%s_hPrimPi0EtaP",histoTag.Data()),binMin,1000);
2606  TH1F* hPrimEtaEtaP = (TH1F*) h2PrimEtaEtaP->ProjectionY(Form("%s_hPrimEtaEtaP",histoTag.Data()),binMin,1000);
2607 
2608  hPrimPhoEtaP->Scale(1./hPrimPhoEtaP->Integral(0,1000));
2609  hPrimPi0EtaP->Scale(1./hPrimPi0EtaP->Integral(0,1000));
2610  hPrimEtaEtaP->Scale(1./hPrimEtaEtaP->Integral(0,1000));
2611 
2612  Float_t maxEta = hPrimPhoEtaP->GetMaximum();
2613  if(maxEta < hPrimPi0EtaP->GetMaximum()) maxEta = hPrimPi0EtaP->GetMaximum();
2614  if(maxEta < hPrimEtaEtaP->GetMaximum()) maxEta = hPrimEtaEtaP->GetMaximum();
2615 
2616  Float_t minEta = hPrimPhoEtaP->GetMinimum();
2617  if(minEta > hPrimPi0EtaP->GetMinimum()) minEta = hPrimPi0EtaP->GetMinimum();
2618  if(minEta > hPrimEtaEtaP->GetMinimum()) minEta = hPrimEtaEtaP->GetMinimum();
2619 
2620  hPrimPi0EtaP->SetMaximum(maxEta*1.1);
2621  hPrimPi0EtaP->SetMinimum(minEta);
2622  TGaxis::SetMaxDigits(3);
2623 
2624  hPrimPi0EtaP->SetYTitle("1/total entries d#it{N}/d#eta");
2625  hPrimPi0EtaP->SetTitle("Generated particles #eta for #it{E} > 3 GeV");
2626  hPrimPi0EtaP->SetTitleOffset(1.5,"Y");
2627  hPrimPi0EtaP->SetMarkerColor(4);
2628  hPrimPi0EtaP->SetMarkerStyle(21);
2629  hPrimPi0EtaP->Draw("");
2630 
2631  hPrimPhoEtaP->SetMarkerColor(1);
2632  hPrimPhoEtaP->SetMarkerStyle(20);
2633  scale = TMath::RadToDeg();
2634  hPrimPhoEtaP->Draw("same");
2635 
2636  hPrimEtaEtaP->SetMarkerColor(2);
2637  hPrimEtaEtaP->SetMarkerStyle(22);
2638  hPrimEtaEtaP->Draw("same");
2639 
2640  if(addOkFlag) text[3]->Draw();
2641 
2642  cmc->Print(Form("%s_MCHisto.%s",histoTag.Data(),format.Data()));
2643 
2644  // cleanup or save
2645  //
2646  if(exportToFile!=1)
2647  {
2648  delete hPrimPhoPhi ;
2649  delete hPrimPi0Phi ;
2650  delete hPrimEtaPhi ;
2651  delete hPrimPhoEtaP ;
2652  delete hPrimPi0EtaP ;
2653  delete hPrimEtaEtaP ;
2654 
2655  delete cmc ;
2656  }
2657  else
2658  {
2659  SaveHisto(hPrimPhoPhi ,kFALSE);
2660  SaveHisto(hPrimPi0Phi ,kFALSE);
2661  SaveHisto(hPrimEtaPhi ,kFALSE);
2662  SaveHisto(hPrimPhoEtaP,kFALSE);
2663  SaveHisto(hPrimPi0EtaP,kFALSE);
2664  SaveHisto(hPrimEtaEtaP,kFALSE);
2665 
2666  SaveCanvas(cmc) ;
2667  }
2668 }
2669 
2675 //____________________________________________________________________
2677 {
2678  if(histArr) delete histArr;
2679 
2680  histArr = (TList*) dir->Get(trigName);
2681 
2682  if ( !histArr )
2683  {
2684  printf("List not found, do nothing\n");
2685  return kFALSE;
2686  }
2687 
2688  if ( histArr->GetEntries() <= 0 )
2689  {
2690  printf("No histograms found <%d>, do nothing\n",histArr->GetEntries());
2691  return kFALSE;
2692  }
2693 
2694  if ( exportToFile == 2 )
2695  {
2696  fout = new TFile(Form("AnalysisResults%s.root",histoTag.Data()),"RECREATE");
2697  histArr->Write();
2698  fout->Close();
2699  }
2700 
2701  return kTRUE ;
2702 }
2703 
2711 //___________________________________
2713 {
2714  TObject *histo = 0x0;
2715 
2716  if ( histArr )
2717  histo = histArr->FindObject(histoName);
2718  else
2719  histo = file->Get (histoName);
2720 
2721  SaveHisto(histo);
2722 
2723  return histo;
2724 }
2725 
2732 //_________________________________________
2733 void SaveHisto(TObject* histo, Bool_t tag)
2734 {
2735  if(histo)
2736  {
2737  if(tag) histo->Write(Form("fig_ga_%s_%s",histoTag.Data(), histo->GetName()));
2738  else histo->Write(Form("fig_ga_%s" ,histo->GetName()));
2739  }
2740 // else
2741 // printf("Object not Available");
2742 }
2743 
2747 //_______________________________
2748 void SaveCanvas(TCanvas* canvas)
2749 {
2750  if(canvas) canvas->Write(Form("canvas_ga_%s",canvas->GetName()));
2751 }
2752 
2756 //___________________________________________________
2757 void ScaleAxis(TAxis *a, Double_t scale)
2758 {
2759  if (!a) return; // just a precaution
2760  if (a->GetXbins()->GetSize())
2761  {
2762  // an axis with variable bins
2763  // note: bins must remain in increasing order, hence the "Scale"
2764  // function must be strictly (monotonically) increasing
2765  TArrayD X(*(a->GetXbins()));
2766  for(Int_t i = 0; i < X.GetSize(); i++) X[i] = scale*X[i];
2767  a->Set((X.GetSize() - 1), X.GetArray()); // new Xbins
2768  }
2769  else
2770  {
2771  // an axis with fix bins
2772  // note: we modify Xmin and Xmax only, hence the "Scale" function
2773  // must be linear (and Xmax must remain greater than Xmin)
2774  a->Set(a->GetNbins(),
2775  - scale*a->GetXmin(), // new Xmin
2776  - scale*a->GetXmax()); // new Xmax
2777  }
2778  return;
2779 }
2780 
2784 //___________________________________________________
2785 void ScaleXaxis(TH1 *h, Double_t scale)
2786 {
2787  if (!h) return; // just a precaution
2788  ScaleAxis(h->GetXaxis(), scale);
2789  return;
2790 }
2791 
2792 
2795 //___________________________________________________
2796 void MyLatexMakeUp(TLatex *currentLatex, Int_t textfont, Double_t textsize, Int_t textcolor)
2797 {
2798  currentLatex->SetNDC();
2799  currentLatex->SetTextFont(textfont);
2800  currentLatex->SetTextSize(textsize);
2801  currentLatex->SetTextColor(textcolor);
2802  return;
2803 }
2804 
Int_t color[]
print message on plot with ok/not ok
void DrawAnaCaloTrackQA(TString listName="Pi0IM_GammaTrackCorr_EMCAL", TString fileName="AnalysisResults.root", Int_t exportTo=1, TString fileFormat="eps", TString outFileName="CaloTrackCorrQA_output", Bool_t okFlag=kTRUE)
double Double_t
Definition: External.C:58
Definition: External.C:236
TCanvas * canvas
Definition: DrawAnaELoss.C:28
void CaloQA(Int_t icalo)
TString fileName
TString format
file names tag, basically the trigger and calorimeter combination
void CorrelQA(Int_t icalo)
TObject * GetHisto(TString histoName)
TString histoTag
output file with plots or extracted histograms
void SaveCanvas(TCanvas *canvas)
void ProcessTrigger(TString trigName="default", Bool_t checkList=kTRUE)
TList * histArr
TDirectory file where lists per trigger are stored in train ouput.
void ScaleAxis(TAxis *a, Double_t scale)
TLatex * text[5]
option to what and if export to output file
int Int_t
Definition: External.C:63
void MCQA(Int_t icalo)
float Float_t
Definition: External.C:68
void IsolQA(Int_t icalo)
void SaveHisto(TObject *histo, Bool_t tag=kTRUE)
Definition: External.C:212
Double_t nEvents
plot quality messages
void MyLatexMakeUp(TLatex *currentLatex, Int_t textfont, Double_t textsize, Int_t textcolor)
Latex settings.
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)
Int_t exportToFile
plots format: eps, pdf, etc.
TFile * file
TList with histograms for a given trigger.
Bool_t GetList(TString trigName)
void Pi0QA(Int_t icalo)
bool Bool_t
Definition: External.C:53
TFile * fout
input train file
Bool_t addOkFlag
number of events analyzed
void ScaleXaxis(TH1 *h, Double_t scale)
void TrackQA()
Definition: External.C:196
TDirectoryFile * dir