AliPhysics  c7b8e89 (c7b8e89)
CompareInvMassGraphsDataProd.C
Go to the documentation of this file.
1 
13 #if !defined(__CINT__) || defined(__MAKECINT__)
14 
15 #include <TFile.h>
16 #include <TDirectoryFile.h>
17 #include <TList.h>
18 #include <TString.h>
19 #include <TROOT.h>
20 #include <TStyle.h>
21 #include <TH1D.h>
22 #include <TH2F.h>
23 #include <TH3F.h>
24 #include <TCanvas.h>
25 #include <TPad.h>
26 #include <TLegend.h>
27 #include <TObject.h>
28 #include <TAxis.h>
29 #include <TGaxis.h>
30 #include <TLine.h>
31 #include <TF1.h>
32 #include <TMath.h>
33 #include <TGraphErrors.h>
34 #include "PlotUtils.C"
35 
36 #endif
37 
38 TString fileFormat = ".eps";
39 
40 //------------------------------------------------------------------------------
46 //------------------------------------------------------------------------------
48 (
49  TString particle = "Eta",
50  TString calorimeter = "DCAL",
51  TString trigger = "default"
52 )
53 {
54  // SM
55  const Int_t nSM = 20;
56  Int_t ncol = 4;
57  Int_t nrow = 4;
58 
59  Int_t firstSM = 0;
60  Int_t lastSM = 11;
61 
62  if(calorimeter=="DCAL")
63  {
64  ncol = 3;
65  nrow = 3;
66 
67  firstSM = 12;
68  lastSM = 19;
69  }
70  //------------------------------------------------
71  // Put the path to the files in prod string array
72  // also the explanation title on prodLeg string array
73  //------------------------------------------------
74 
75  const Int_t nProd = 5;
76 
77  TString prod[] =
78  {
79  "LHC17_5TeV/CENTwoSDD"
80  ,"LHC17_5TeV/FAST"
81  ,"LHC17_13TeV"
82  ,"LHC16_13TeV"
83  ,"LHC15n"
84  };
85 
86  TString prodLeg[] =
87  {
88  "LHC17pq, #sqrt{s} = 5 TeV, CENT"
89  ,"LHC17pq, #sqrt{s} = 5 TeV, FAST"
90  ,"LHC17*, #sqrt{s} = 13 TeV"
91  ,"LHC16*, #sqrt{s} = 13 TeV"
92  ,"LHC15n, #sqrt{s} = 5 TeV"
93  };
94 
95  TString directory = Form("ComparisonGraphs/%s/%s/%s",calorimeter.Data(),trigger.Data(),particle.Data());
96 
97 // const Int_t nProd = 8;
98 //
99 // TString prod[] =
100 // {
101 // "LHC17_5TeV/LHC17p/CENTwoSDD"
102 // ,"LHC17_13TeV/LHC17i"
103 // ,"LHC17_13TeV/LHC17l"
104 // ,"LHC17_13TeV/LHC17o"
105 // ,"LHC16_13TeV/LHC16i"
106 // ,"LHC16_13TeV/LHC16k"
107 // ,"LHC16_13TeV/LHC16o"
108 // ,"LHC15n"
109 // };
110 //
111 // TString prodLeg[] =
112 // {
113 // "LHC17p, #sqrt{s} = 5 TeV, CENT"
114 // ,"LHC17i, #sqrt{s} = 13 TeV"
115 // ,"LHC17l, #sqrt{s} = 13 TeV"
116 // ,"LHC17o, #sqrt{s} = 13 TeV"
117 // ,"LHC16i, #sqrt{s} = 13 TeV"
118 // ,"LHC16l, #sqrt{s} = 13 TeV"
119 // ,"LHC16o, #sqrt{s} = 13 TeV"
120 // ,"LHC15n, #sqrt{s} = 5 TeV"
121 // };
122 //
123 // TString directory = Form("ComparisonGraphsPeriods/%s/%s/%s",calorimeter.Data(),trigger.Data(),particle.Data());
124 
125  // Create output directories
126  TString processline = Form(".! mkdir -p %s",directory.Data()) ;
127  gROOT->ProcessLine(processline.Data());
128 
129  // Init the graph arrays.
130  //
131  TGraphErrors* gMass[nSM][nProd];
132  TGraphErrors* gMassAllSM[nProd];
133  TGraphErrors* gMassAll [nProd];
134  TH1F* hAxis[nSM];
135 
136  TFile* file[nProd];
137 
138  TGraphErrors* gMassRat[nSM][nProd];
139  TGraphErrors* gMassRatAllSM[nProd];
140  TGraphErrors* gMassRatAll [nProd];
141 
142  // Init to 0 all arrays
143  for(Int_t iprod = 0; iprod < nProd; iprod++)
144  {
145  gMassAllSM [iprod] = 0;
146  gMassAll [iprod] = 0;
147  gMassRatAllSM[iprod] = 0;
148  gMassRatAll [iprod] = 0;
149  file [iprod] = 0;
150 
151  for(Int_t ism = 0; ism < nSM; ism++)
152  {
153  gMass [ism][iprod] = 0;
154  gMassRat[ism][iprod] = 0;
155  if ( iprod == 0 ) hAxis[ism] = 0;
156  }
157  }
158 
159  //===============================
160  // Open the input files and recover the graphs
161  // apply some style settings
162  //===============================
163 
164  Float_t colorProd[] = { 1,kBlue,kRed,kViolet,kRed,kOrange-2,kCyan, kBlue};
165  Float_t styleProd[] = {20, 24, 24, 24, 25, 25, 25, 20};
166  Float_t markerSize = 8;
167 
168  for(Int_t iprod = 0; iprod < nProd; iprod++)
169  {
170  //printf("iprod %d\n",iprod);
171  file[iprod] = new TFile(Form("%s/%s/%s/%s/MassWidthPtHistograms.root",
172  prod[iprod].Data(),calorimeter.Data(),trigger.Data(),particle.Data()),"read");
173 
174  printf("iprod %d, %s/%s/%s/%s/MassWidthPtHistograms.root %p\n",
175  iprod,prod[iprod].Data(),calorimeter.Data(),trigger.Data(),particle.Data(),file[iprod]);
176 
177  if(!file[iprod]) continue;
178 
179  gMassAllSM[iprod] = (TGraphErrors*) file[iprod]->Get("gMass_SameSM");
180  //printf("\t gMass_SameSM %p\n",gMassAllSM[iprod]);
181 
182  gMassAllSM[iprod]->SetMarkerColor(colorProd[iprod]);
183  gMassAllSM[iprod]->SetLineColor (colorProd[iprod]);
184  gMassAllSM[iprod]->SetMarkerStyle(styleProd[iprod]);
185  gMassAllSM[iprod]->SetMarkerSize (markerSize);
186 
187 
188  gMassAll[iprod] = (TGraphErrors*) file[iprod]->Get("gMass_AllSM");
189  //printf("\t gMass_SameSM %p\n",gMassAllSM[iprod]);
190 
191  gMassAll[iprod]->SetMarkerColor(colorProd[iprod]);
192  gMassAll[iprod]->SetLineColor (colorProd[iprod]);
193  gMassAll[iprod]->SetMarkerStyle(styleProd[iprod]);
194  gMassAll[iprod]->SetMarkerSize (markerSize);
195 
196  if ( iprod > 0 && gMassAllSM[0] )
197  {
198  gMassRatAllSM[iprod] = DivideGraphs(gMassAllSM[iprod],gMassAllSM[0]);
199  gMassRatAllSM[iprod]->SetMarkerColor(colorProd[iprod]);
200  gMassRatAllSM[iprod]->SetLineColor (colorProd[iprod]);
201  gMassRatAllSM[iprod]->SetMarkerStyle(styleProd[iprod]);
202  gMassRatAllSM[iprod]->SetMarkerSize (markerSize);
203 
204  gMassRatAll[iprod] = DivideGraphs(gMassAll[iprod],gMassAll[0]);
205  gMassRatAll[iprod]->SetMarkerColor(colorProd[iprod]);
206  gMassRatAll[iprod]->SetLineColor (colorProd[iprod]);
207  gMassRatAll[iprod]->SetMarkerStyle(styleProd[iprod]);
208  gMassRatAll[iprod]->SetMarkerSize (markerSize);
209  }
210 
211  for(Int_t ism = firstSM; ism <= lastSM; ism++)
212  {
213  gMass[ism][iprod] = (TGraphErrors*) file[iprod]->Get(Form("gMass_SM%d",ism));
214  //printf("\t gMass_SM%d %p\n",ism,gMass[ism][iprod]);
215 
216  if(!gMass[ism][iprod] ) continue;
217 
218  gMass[ism][iprod]->SetMarkerColor(colorProd[iprod]);
219  gMass[ism][iprod]->SetLineColor (colorProd[iprod]);
220  gMass[ism][iprod]->SetMarkerStyle(styleProd[iprod]);
221  gMass[ism][iprod]->SetMarkerSize (markerSize);
222 
223  if ( iprod > 0 && gMass[ism][0] )
224  {
225  gMassRat[ism][iprod] = DivideGraphs(gMass[ism][iprod],gMass[ism][0]);
226  gMassRat[ism][iprod]->SetMarkerColor(colorProd[iprod]);
227  gMassRat[ism][iprod]->SetLineColor (colorProd[iprod]);
228  gMassRat[ism][iprod]->SetMarkerStyle(styleProd[iprod]);
229  gMassRat[ism][iprod]->SetMarkerSize (markerSize);
230  }
231 
232  } // ism
233 
234  } // iprod
235 
236  //===============================
237  // PLOTS
238  //===============================
239  gStyle->SetPadRightMargin(0.01);
240  gStyle->SetPadLeftMargin(0.12);
241  gStyle->SetTitleFontSize(0.06);
242  // gStyle->SetStatFontSize(0.06);
243  gStyle->SetOptStat(0);
244  gStyle->SetOptTitle(1);
245  gStyle->SetPadTopMargin(0.07);
246 
247  //-----------------------------------------------------
248  // Plot per SM the Mass distributions and their ratios
249  //-----------------------------------------------------
251 
252  TLegend *lE = new TLegend(-0.04,0.1,1,1);
253  lE->SetFillColor(0);
254  lE->SetFillStyle(0);
255  lE->SetLineColor(0);
256  lE->SetBorderSize(0);
257  lE->SetTextSize(0.07);
258  if ( trigger == "default" ) lE->SetHeader(" kINT7");
259  else if ( trigger == "EMCAL_L0") lE->SetHeader(" kEMC7");
260  else if ( trigger == "DCAL_L0" ) lE->SetHeader(" kDMC7");
261  else if ( trigger == "EMCAL_L1") lE->SetHeader(" kEMCEGA: EG1");
262  else if ( trigger == "DCAL_L1" ) lE->SetHeader(" kEMCEGA: DG1");
263  else if ( trigger == "EMCAL_L2") lE->SetHeader(" kEMCEGA: EG2");
264  else if ( trigger == "DCAL_L2" ) lE->SetHeader(" kEMCEGA: DG2");
265 
266  // Axis E range
267  Float_t minE = 1.5;
268  Float_t maxE = 12;
269  if(particle=="Eta")
270  {
271  minE = 2;
272  maxE = 20;
273  }
274 
275  TCanvas * cGraphSM = new TCanvas("cGraphSM","Per SM",
276  ncol*2000,nrow*2000);
277 
278  cGraphSM->Divide(ncol,nrow);
279 
280  for(Int_t ism = firstSM; ism <= lastSM; ism++)
281  {
282  cGraphSM->cd(ism+1-firstSM);
283 
284  //gPad->SetGridx();
285  gPad->SetGridy();
286 
287  hAxis[ism] = new TH1F(Form("hAxisBisRat_SM%d",ism),Form("SM %d",ism),1000,minE,maxE);
288 
289  hAxis[ism]->GetYaxis()->SetTitle("Mass (MeV/#it{c}^{2})");
290  hAxis[ism]->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
291  hAxis[ism]->SetMaximum(145);
292  hAxis[ism]->SetMinimum(125);
293 // hAxis[ism]->SetMaximum(137);
294 // hAxis[ism]->SetMinimum(123);
295 
296 // hAxis[ism]->SetMaximum(140);
297 // hAxis[ism]->SetMinimum(128);
298 
299  if(particle=="Eta")
300  {
301  hAxis[ism]->SetMaximum(600);
302  hAxis[ism]->SetMinimum(450);
303  }
304 
305  hAxis[ism]->Draw("");
306 
307  for(Int_t iprod = 0; iprod < nProd; iprod++)
308  {
309  if ( !file[iprod] || !gMass[ism][iprod] ) continue;
310 
311  if ( ism == firstSM )
312  lE->AddEntry(gMass[ism][iprod],prodLeg[iprod],"PL");
313 
314  gMass[ism][iprod]->Draw("P");
315  }
316 
317  }
318  cGraphSM->cd(lastSM-firstSM+2);
319  lE->Draw();
320 
321  fileName = Form("%s/Mass_PerSM",
322  directory.Data());
323 
324  fileName+=fileFormat;
325  cGraphSM->Print(fileName);
326 
327  // Ratio
328  //
329  TCanvas * cGraphSMRat = new TCanvas("cGraphSMRat","Ratio Per SM",
330  ncol*2000,nrow*2000);
331 
332  cGraphSMRat->Divide(ncol,nrow);
333 
334  for(Int_t ism = firstSM; ism <= lastSM; ism++)
335  {
336  cGraphSMRat->cd(ism+1-firstSM);
337 
338  //gPad->SetGridx();
339  gPad->SetGridy();
340 
341  hAxis[ism] = new TH1F(Form("hAxisBis_SM%d",ism),Form("SM %d",ism),1000,minE,maxE);
342 
343  hAxis[ism]->GetYaxis()->SetTitle("Ratio to LHC17pq");
344  hAxis[ism]->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
345 // hAxis[ism]->SetMaximum(1.002);
346 // hAxis[ism]->SetMinimum(0.932);
347 // hAxis[ism]->SetMaximum(1.05);
348 // hAxis[ism]->SetMinimum(0.90);
349  hAxis[ism]->SetMaximum(1.1);
350  hAxis[ism]->SetMinimum(0.9);
351  hAxis[ism]->Draw("");
352 
353  for(Int_t iprod = 1; iprod < nProd; iprod++)
354  {
355  if ( !gMassRat[ism][iprod] ) continue;
356 
357  //printf("%d %p\n",iprod,gMassRat[ism][iprod]);
358  gMassRat[ism][iprod]->Draw("P");
359  }
360 
361  }
362 
363  cGraphSMRat->cd(lastSM-firstSM+2);
364 
365  lE->Draw();
366 
367  fileName = Form("%s/MassRatio_PerSM",
368  directory.Data());
369 
370  fileName+=fileFormat;
371  cGraphSMRat->Print(fileName);
372 
373 
374  //-----------------------------------------------------
375  // Plot sum of all SM or group of SM
376  // the Mass distributions and their ratios
377  //-----------------------------------------------------
378 
379  ncol = 2;
380  nrow = 2;
381 
382  TH1F* hAxisAll = new TH1F("hAxisAll","All pairs in same SM",1000,minE,maxE);
383 
384  hAxisAll->GetYaxis()->SetTitle("Mass (MeV/#it{c}^{2})");
385  hAxisAll->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
386  hAxisAll->SetMaximum(140);
387  hAxisAll->SetMinimum(125);
388 
389  if(particle=="Eta")
390  {
391  hAxisAll->SetMaximum(600);
392  hAxisAll->SetMinimum(450);
393  }
394 
395 // hAxisAll->SetMaximum(140);
396 // hAxisAll->SetMinimum(115);
397 
398  TH1F* hAxisAllRat = new TH1F("hAxisAllRat","Ratio to LHC17pq CENT, all SM, pairs in same SM",1000,minE,maxE);
399  hAxisAllRat->GetYaxis()->SetTitle("Ratio to LHC17pq CENT");
400  hAxisAllRat->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
401 // hAxisAllRat->SetMaximum(1.002);
402 // hAxisAllRat->SetMinimum(0.932);
403 
404 // hAxisAllRat->SetMaximum(1.02);
405 // hAxisAllRat->SetMinimum(0.92);
406  hAxisAllRat->SetMaximum(1.05);
407  hAxisAllRat->SetMinimum(0.95);
408 
409 
410  TCanvas * cGraph = new TCanvas("cGraphAllSameSM","All same SM",
411  ncol*2000,nrow*2000);
412 
413  cGraph->Divide(ncol,nrow);
414 
415  cGraph->cd(1);
416 
417  //gPad->SetGridx();
418  gPad->SetGridy();
419 
420  hAxisAll->Draw("");
421 
422  for(Int_t iprod = 0; iprod < nProd; iprod++)
423  {
424  if(!file[iprod]) continue;
425  gMassAllSM[iprod]->Draw("P");
426  }
427 
428  gMassAllSM[1]->Draw("P");
429 
430  cGraph->cd(2);
431 
432  //gPad->SetGridx();
433  gPad->SetGridy();
434 
435  TLegend *lA = new TLegend(0.2,0.8,0.9,0.9);
436  lA->SetFillColor(0);
437  lA->SetFillStyle(0);
438  lA->SetLineColor(0);
439  lA->SetBorderSize(0);
440  lA->SetTextSize(0.05);
441 
442  hAxisAllRat->Draw("");
443 
444  for(Int_t iprod = 1; iprod < nProd; iprod++)
445  {
446  if(!file[iprod]) continue;
447 
448  gMassRatAllSM[iprod]->Draw("P");
449  //gMassRatAllSM[iprod]->Fit("pol0","R","",2,8);
450  //TF1* func = gMassRatAllSM[iprod]->GetFunction("pol0");
451  //lA->AddEntry(func,Form("Fit %2.3f #pm %2.3f",func->GetParameter(0),func->GetParError(0)),"L");
452  }
453  gMassRatAllSM[1]->Draw("P");
454 
455  lA->Draw();
456 
457  cGraph->cd(3);
458 
459  lE->Draw();
460 
461  //cGraph->cd(ncol*nrow);
462  //lE->Draw();
463 
464  fileName = Form("%s/Mass_AllSameSM",
465  directory.Data());
466 
467  fileName+=fileFormat;
468  cGraph->Print(fileName);
469 
470 
471  TCanvas * cGraphA = new TCanvas("cGraphAllSM","All SM",
472  ncol*2000,nrow*2000);
473 
474  cGraphA->Divide(ncol,nrow);
475 
476  cGraphA->cd(1);
477 
478  //gPad->SetGridx();
479  gPad->SetGridy();
480 
481  hAxisAll->SetTitle("all SM combinations");
482  hAxisAll->Draw("");
483 
484  for(Int_t iprod = 0; iprod < nProd; iprod++)
485  {
486  if(!file[iprod]) continue;
487  gMassAll[iprod]->Draw("P");
488  }
489 
490  gMassAll[1]->Draw("P");
491 
492  cGraphA->cd(2);
493 
494  //gPad->SetGridx();
495  gPad->SetGridy();
496 
497  hAxisAllRat->SetTitle("Ratio to LHC17pq CENT, all SM combinations");
498  hAxisAllRat->Draw("");
499 
500  for(Int_t iprod = 1; iprod < nProd; iprod++)
501  {
502  if(!file[iprod]) continue;
503 
504  gMassRatAllSM[iprod]->Draw("P");
505  //gMassRatAllSM[iprod]->Fit("pol0","R","",2,8);
506  //TF1* func = gMassRatAllSM[iprod]->GetFunction("pol0");
507  //lA->AddEntry(func,Form("Fit %2.3f #pm %2.3f",func->GetParameter(0),func->GetParError(0)),"L");
508  }
509 
510  gMassRatAllSM[1]->Draw("P");
511 
512  lA->Draw();
513 
514  cGraphA->cd(3);
515 
516  lE->Draw();
517 
518  //cGraph->cd(ncol*nrow);
519  //lE->Draw();
520 
521  fileName = Form("%s/Mass_AllSM",
522  directory.Data());
523 
524  fileName+=fileFormat;
525  cGraphA->Print(fileName);
526 }
527 
528 
TString fileName
int lastSM
Plotting utilities.
const TString calorimeter
Definition: anaM.C:36
static TGraphErrors * DivideGraphs(TGraphErrors *gNum, TGraphErrors *gDen)
Definition: PlotUtils.C:151
int Int_t
Definition: External.C:63
void CompareInvMassGraphsDataProd(TString particle="Eta", TString calorimeter="DCAL", TString trigger="default")
File format: ".eps" ".pdf" ...
const Int_t nProd
float Float_t
Definition: External.C:68
TString prodLeg[]
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)
TFile * file
TList with histograms for a given trigger.
TString fileFormat
TString prod[]
productions to be compared, directory name