Physics selection QA
[u/mrichter/AliRoot.git] / PWG1 / macros / PhysSel / PlotQAPhysicsSelection.C
1 //using namespace std;
2
3 #include <iostream>
4
5 #include "QALHC11c_QA69.h"
6 #include "DefinePlots.h"
7 #include "TGraphErrors.h"
8
9
10 const char * GetLocalFileName(Int_t run, const char * suffix, const char * path);
11 TH1F * GetEmpty(const char * name, Int_t nfile);
12 TGraphErrors * GetGraphRej(TGraphErrors * gr, TList * rejRunList, const char * reason, Float_t &mean, Bool_t doDraw) ;
13 double GetMedian(double* arr, int n);
14 double meanMed(double* vec, int np, double nsigmaCut, int &nrej, int *rejList);
15 TGraphErrors * graph[kNGraphs];
16
17 //---------------------------------------------------------
18 TH1F *fAccOverAll[3];
19 TH1F *fBGOverAll[3];
20 TH1F *fV0BGOverAll[3];
21 TH1F *fV0ABGOverAll[3];
22 TH1F *fV0CBGOverAll[3];
23 TH1F *fAfterOverBefore[3];
24 TH1F *fV0AOverV0C[3];
25 TH1F *fF02OverAll[3];
26 TH1F *fF01OverAll[3];
27
28 void PlotQAPhysicsSelection() {
29
30
31   for(Int_t i=0;i<3;++i){
32     fAccOverAll[i]=new TH1F(Form("fAccOverAll_%d",i),Form("fAccOverAll_%d",i),100,0,1);
33     fBGOverAll[i]=new TH1F(Form("fBGOverAll_%d",i),Form("fBGOverAll_%d",i),50,0,0.5);
34     fV0BGOverAll[i]=new TH1F(Form("fV0BGOverAll_%d",i),Form("fV0BGOverAll_%d",i),50,0,0.1);
35     fV0ABGOverAll[i]=new TH1F(Form("fV0ABGOverAll_%d",i),Form("fV0ABGOverAll_%d",i),50,0,0.1);
36     fV0CBGOverAll[i]=new TH1F(Form("fV0CBGOverAll_%d",i),Form("fV0CBGOverAll_%d",i),50,0,0.1);
37     fAfterOverBefore[i]=new TH1F(Form("fAfterOverBefore_%d",i),Form("fAfterOverBefore_%d",i),100,0.5,1);
38     fV0AOverV0C[i]=new TH1F(Form("fV0AOverV0C_%d",i),Form("fV0AOverV0C_%d",i),100,0.7,1.2);
39     fF01OverAll[i]=new TH1F(Form("fF01OverAll_%d",i),Form("fF01OverAll_%d",i),200,0.6,1.1);
40     fF02OverAll[i]=new TH1F(Form("fF02OverAll_%d",i),Form("fF02OverAll_%d",i),200,0.6,1.1);
41   }
42
43
44   gStyle->SetOptStat(0);
45   gStyle->SetPalette(1);
46   gStyle->SetCanvasColor(10);
47   gStyle->SetFrameFillColor(10);
48   gStyle->SetOptStat(0);
49   gStyle->SetOptTitle(0);
50
51   // Book graphs
52   TGraphErrors * graph[kNGraphs] = {0};
53   for(Int_t igraph = 0; igraph < kNGraphs; igraph++){
54     graph[igraph] = new TGraphErrors;
55     graph[igraph]->SetName(gnames[igraph]);
56     graph[igraph]->GetXaxis()->SetTitle("run");
57     graph[igraph]->GetYaxis()->SetTitle(ylabels[igraph]);
58     graph[igraph]->SetMarkerStyle(20);
59   }
60
61   //loading libraries
62   loadlibs();
63
64   TList * listEmptyRuns    = new TList(); // keep track of empty runs
65   TList * listRejectedRuns = new TList(); // keep track of rejected runes
66
67   // Count ncycles
68   Int_t cycle=0;
69   while(QAcycle[++cycle]>0) {
70     
71     cout << "." ;
72   }  
73   cout << "  Ncycles " << cycle <<endl;
74   cout << QAcycle[0] ;
75
76   // connect to grid
77   if(!localMode) TGrid::Connect("alien://");  
78
79   // do not use scientific notation for run number
80   TGaxis::SetMaxDigits(7)  ;
81
82
83
84
85   
86   // loop over all files
87   Int_t ifile =-1;
88   Int_t ifileGood = 0;
89   Int_t ifileNotEmpty  = 0;
90   while (runs[++ifile] > 0) {
91     
92     Long_t *id,*size,*flags,*mt;
93     
94
95
96
97     // check if QA file is available and open it
98     // try different QA train outputs
99     TString file ;
100     TFile *fr=0;
101     TFile *fc=0; // centrality, only in local mode for the time being
102     for(Int_t c=0;c<cycle;c++){ 
103       if(fr)break; 
104       if(output.Contains("AOD"))
105         file.Form("alien://%s/%09d/%s%3.3d/event_stat.root",location.Data(),runs[ifile],
106                   output.Data(),QAcycle[c] );
107       else
108         file.Form("alien://%s/%09d/%s%d/event_stat.root",location.Data(),runs[ifile],
109                   output.Data(),QAcycle[c] );
110       if(localMode){ 
111         if(!GetLocalFileName(runs[ifile], localSuffix, localPath))continue;
112         file = GetLocalFileName(runs[ifile], localSuffix, localPath);
113         if(gSystem->GetPathInfo(file,id,size,flags,mt)){Printf("not found");continue;}
114       }
115   
116       Printf("\nBegin of reading: %s", file.Data());
117       fr=TFile::Open(file);     
118
119     }
120  
121
122
123     // set value to -1 by default:
124     for(Int_t igraph = 0; igraph < kNGraphs; igraph++){
125       graph[igraph]->SetPoint(ifile, runs[ifile], -1);
126     }
127
128     // If the file is not available, continue
129     if(!fr){
130       Printf("File %d is not available.\n",runs[ifile]);
131       listEmptyRuns->Add(new TObjString(Form("File not Existing [%d]",runs[ifile])));
132       continue;
133     }
134
135
136     // get stats and fill graphs
137     TH2F * hStats = (TH2F*) fr->Get("fHistStatistics"); 
138     if(!localMode) {
139       gSystem->Exec(Form("alien_cp %s %s",file.Data(), GetLocalFileName(runs[ifile], localSuffix, localPath)));
140       cout << Form("alien_cp %s %s",file.Data(), GetLocalFileName(runs[ifile], localSuffix, localPath)) <<endl;
141     }
142
143  
144
145     Int_t rowC0SMH   = -1;
146     Int_t rowCMBAC   = -1;
147     Int_t rowCMBACS2 = -1;
148     Int_t nbiny = hStats->GetNbinsY();
149     for(Int_t ibiny = 1; ibiny <= nbiny; ibiny++){
150       TString label = hStats->GetYaxis()->GetBinLabel(ibiny);
151       if(label.Contains(trigger3))   rowC0SMH   = ibiny; 
152       if(label.Contains(trigger1))   rowCMBAC   = ibiny; 
153       if(label.Contains(trigger2)) rowCMBACS2 = ibiny; 
154     }
155
156     //Number of events in the selected trigger class
157     Float_t C0SMH   = hStats->GetBinContent(AliPhysicsSelection::kStatTriggerClass,rowC0SMH);    
158     Float_t CMBAC   = hStats->GetBinContent(AliPhysicsSelection::kStatTriggerClass,rowCMBAC);
159     Float_t CMBACS2 = hStats->GetBinContent(AliPhysicsSelection::kStatTriggerClass,rowCMBACS2);
160
161
162     //Number of events after physics selection
163     Float_t C0SMH_APS   = hStats->GetBinContent(AliPhysicsSelection::kStatOffline,rowC0SMH);    
164     Float_t CMBAC_APS   = hStats->GetBinContent(AliPhysicsSelection::kStatOffline,rowCMBAC);
165     Float_t CMBACS2_APS = hStats->GetBinContent(AliPhysicsSelection::kStatOffline,rowCMBACS2);
166
167
168     //Fraction of Events rejected as background because out of the correlation tracklets vs clusters
169     Float_t C0SMHBG   = hStats->GetBinContent(AliPhysicsSelection::kStatBG,rowC0SMH);    
170     Float_t CMBACBG   = hStats->GetBinContent(AliPhysicsSelection::kStatBG,rowCMBAC);
171     Float_t CMBACS2BG = hStats->GetBinContent(AliPhysicsSelection::kStatBG,rowCMBACS2);
172
173  
174
175
176     //Events with a signal from V0A in the collision time window
177     Float_t C0SMHV0A   = hStats->GetBinContent(AliPhysicsSelection::kStatV0A,rowC0SMH);    
178     Float_t CMBACV0A   = hStats->GetBinContent(AliPhysicsSelection::kStatV0A,rowCMBAC);
179     Float_t CMBACS2V0A = hStats->GetBinContent(AliPhysicsSelection::kStatV0A,rowCMBACS2);
180     
181     //Events with a signal from V0C in the collision time window
182     Float_t C0SMHV0C   = hStats->GetBinContent(AliPhysicsSelection::kStatV0C,rowC0SMH);    
183     Float_t CMBACV0C   = hStats->GetBinContent(AliPhysicsSelection::kStatV0C,rowCMBAC);
184     Float_t CMBACS2V0C = hStats->GetBinContent(AliPhysicsSelection::kStatV0C,rowCMBACS2);
185
186     if(C0SMHV0C>0)graph[kGraphV0AOverV0CSMH]->SetPoint(ifile,runs[ifile], C0SMHV0A/C0SMHV0C);
187     if(CMBACV0C>0)graph[kGraphV0AOverV0CAC]->SetPoint(ifile,runs[ifile], CMBACV0A/CMBACV0C);
188     if(CMBACS2V0C>0)graph[kGraphV0AOverV0CAC2]->SetPoint(ifile,runs[ifile], CMBACS2V0A/CMBACS2V0C);
189
190  
191     //Events flagged as BG from V0 side A
192     Float_t C0SMHV0ABG   = hStats->GetBinContent(AliPhysicsSelection::kStatV0ABG,rowC0SMH);    
193     Float_t CMBACV0ABG   = hStats->GetBinContent(AliPhysicsSelection::kStatV0ABG,rowCMBAC);
194     Float_t CMBACS2V0ABG = hStats->GetBinContent(AliPhysicsSelection::kStatV0ABG,rowCMBACS2);
195     
196     //Events flagged as BG from V0 side C
197     Float_t C0SMHV0CBG   = hStats->GetBinContent(AliPhysicsSelection::kStatV0CBG,rowC0SMH);    
198     Float_t CMBACV0CBG   = hStats->GetBinContent(AliPhysicsSelection::kStatV0CBG,rowCMBAC);
199     Float_t CMBACS2V0CBG = hStats->GetBinContent(AliPhysicsSelection::kStatV0CBG,rowCMBACS2);
200
201
202     //Number of events with more than 1 chip hit in the pixels, computed offline
203     Float_t C0SMHF01   = hStats->GetBinContent(AliPhysicsSelection::kStatFO1,rowC0SMH);    
204     Float_t CMBACF01   = hStats->GetBinContent(AliPhysicsSelection::kStatFO1,rowCMBAC);
205     Float_t CMBACS2F01 = hStats->GetBinContent(AliPhysicsSelection::kStatFO1,rowCMBACS2);
206
207
208     //Number of events with more than 2 chip hit in the pixels, computed offline
209     Float_t C0SMHF02   = hStats->GetBinContent(AliPhysicsSelection::kStatFO2,rowC0SMH);    
210     Float_t CMBACF02   = hStats->GetBinContent(AliPhysicsSelection::kStatFO2,rowCMBAC);
211     Float_t CMBACS2F02 = hStats->GetBinContent(AliPhysicsSelection::kStatFO2,rowCMBACS2);
212
213
214     //Accepted by V0
215     Float_t C0SMHV0   = hStats->GetBinContent(AliPhysicsSelection::kStatV0,rowC0SMH);    
216     Float_t CMBACV0   = hStats->GetBinContent(AliPhysicsSelection::kStatV0,rowCMBAC);
217     Float_t CMBACS2V0 = hStats->GetBinContent(AliPhysicsSelection::kStatV0,rowCMBACS2);
218
219
220
221     //Events flagged as BG from V0 and the A or C side
222     Float_t C0SMHV0BG   = hStats->GetBinContent(AliPhysicsSelection::kStatV0ABG,rowC0SMH)  +hStats->GetBinContent(AliPhysicsSelection::kStatV0CBG,rowC0SMH)  ;
223     Float_t CMBACV0BG   = hStats->GetBinContent(AliPhysicsSelection::kStatV0ABG,rowCMBAC)  +hStats->GetBinContent(AliPhysicsSelection::kStatV0CBG,rowCMBAC)  ;
224     Float_t CMBACS2V0BG = hStats->GetBinContent(AliPhysicsSelection::kStatV0ABG,rowCMBACS2)+hStats->GetBinContent(AliPhysicsSelection::kStatV0CBG,rowCMBACS2);
225
226
227    //Events passing the ZDC time cut on the correlation between the sum and the difference of the timing (rejects the "debunched" events)
228     
229     Float_t C0SMHZDC   = hStats->GetBinContent(AliPhysicsSelection::kStatZDCTime,rowC0SMH);
230     Float_t CMBACZDC   = hStats->GetBinContent(AliPhysicsSelection::kStatZDCTime,rowCMBAC);
231     Float_t CMBACS2ZDC = hStats->GetBinContent(AliPhysicsSelection::kStatZDCTime,rowCMBACS2);
232                                                                     
233
234     //Accepted events                                               kStatAccepted
235     Float_t C0SMHACC   = hStats->GetBinContent(AliPhysicsSelection::kStatAccepted,rowC0SMH);
236     Float_t CMBACACC   = hStats->GetBinContent(AliPhysicsSelection::kStatAccepted,rowCMBAC);
237     Float_t CMBACS2ACC = hStats->GetBinContent(AliPhysicsSelection::kStatAccepted,rowCMBACS2);
238
239
240     cout<<"C0SMHBG="<<C0SMHBG<<"   C0SMHACC="<<C0SMHACC<<endl;
241
242
243
244
245
246     cout << runs[ifile] << " " << CMBACS2 << " " << CMBACS2V0BG << " " << " " << CMBACS2ACC << endl;
247     
248   
249    
250
251    
252     
253
254     if(CMBAC>0) {
255       graph[kGraphBGOverAllAC]->SetPoint(ifile,runs[ifile], CMBACBG / CMBAC);
256       graph[kGraphV0BGOverAllAC]->SetPoint(ifile,runs[ifile], CMBACV0BG / CMBAC);
257       graph[kGraphV0ABGOverAllAC]->SetPoint(ifile,runs[ifile], CMBACV0ABG / CMBAC);//
258       graph[kGraphV0CBGOverAllAC]->SetPoint(ifile,runs[ifile], CMBACV0CBG / CMBAC);//
259       graph[kGraphACCOverAllAC]    ->SetPoint(ifile,runs[ifile],CMBACACC/CMBAC);
260       graph[kGraphNevACratioAfter]    ->SetPoint(ifile,runs[ifile],CMBAC_APS/CMBAC);
261       //F0
262       graph[kGraphFO1OverAllAC]    ->SetPoint(ifile,runs[ifile],CMBACF01/CMBAC);
263       graph[kGraphFO2OverAllAC]    ->SetPoint(ifile,runs[ifile],CMBACF02/CMBAC);
264
265
266
267     }
268     if(CMBACS2>0){
269       graph[kGraphBGOverAllAC2]->SetPoint(ifile,runs[ifile], CMBACS2BG / CMBACS2);
270       graph[kGraphV0BGOverAllAC2]->SetPoint(ifile,runs[ifile], CMBACS2V0BG / CMBACS2);
271       graph[kGraphV0ABGOverAllAC2]->SetPoint(ifile,runs[ifile], CMBACS2V0ABG / CMBACS2);//
272       graph[kGraphV0CBGOverAllAC2]->SetPoint(ifile,runs[ifile], CMBACS2V0CBG / CMBACS2);//
273       graph[kGraphACCOverAllACS2]->SetPoint(ifile,runs[ifile],CMBACS2ACC/CMBACS2);
274       graph[kGraphNevAC2ratioAfter]->SetPoint(ifile,runs[ifile],CMBACS2_APS/CMBACS2);
275
276       //FO
277       graph[kGraphFO1OverAllAC2]    ->SetPoint(ifile,runs[ifile],CMBACS2F01/CMBACS2);
278       graph[kGraphFO2OverAllAC2]    ->SetPoint(ifile,runs[ifile],CMBACS2F02/CMBACS2);
279
280
281
282     }
283     if(C0SMH>0){
284       graph[kGraphBGOverAllSMH]->SetPoint(ifile,runs[ifile], C0SMHBG / C0SMH);
285       graph[kGraphV0BGOverAllSMH]->SetPoint(ifile,runs[ifile], C0SMHV0BG / C0SMH);
286       graph[kGraphV0ABGOverAllSMH]->SetPoint(ifile,runs[ifile],  C0SMHV0ABG/ C0SMH);
287       graph[kGraphV0CBGOverAllSMH]->SetPoint(ifile,runs[ifile],  C0SMHV0CBG/ C0SMH);
288       graph[kGraphACCOverAllSMH]   ->SetPoint(ifile,runs[ifile],C0SMHACC/C0SMH);
289       graph[kGraphNevSMHratioAfter]   ->SetPoint(ifile,runs[ifile],C0SMH_APS/C0SMH);
290
291       graph[kGraphFO1OverAllSMH]    ->SetPoint(ifile,runs[ifile], C0SMHF01/C0SMH);
292       graph[kGraphFO2OverAllSMH]    ->SetPoint(ifile,runs[ifile], C0SMHF02/C0SMH);
293
294
295     } 
296     
297     
298     
299     
300
301     graph[kGraphNevSMH]->SetPoint(ifile,runs[ifile],C0SMHACC);
302     graph[kGraphNevAC ]->SetPoint(ifile,runs[ifile],CMBACACC);
303     graph[kGraphNevAC2]->SetPoint(ifile,runs[ifile],CMBACS2ACC);
304     graph[kGraphNevSMHBefore]->SetPoint(ifile,runs[ifile],C0SMH);
305     graph[kGraphNevACBefore ]->SetPoint(ifile,runs[ifile],CMBAC);
306     graph[kGraphNevAC2Before]->SetPoint(ifile,runs[ifile],CMBACS2);
307     graph[kGraphNevSMHAfter]->SetPoint(ifile,runs[ifile],C0SMH_APS);
308     graph[kGraphNevACAfter ]->SetPoint(ifile,runs[ifile],CMBAC_APS);
309     graph[kGraphNevAC2After]->SetPoint(ifile,runs[ifile],CMBACS2_APS);
310
311  
312
313     ifileGood++;
314     if(C0SMH>0 || CMBACS2>0 || CMBAC>0) ifileNotEmpty++;
315     else  listEmptyRuns->Add(new TObjString(Form("No events [%d]",runs[ifile])));
316  
317
318   }
319
320
321  
322
323
324   // Set bin labels with run number and save graphs
325   // also prepare a table with number of events
326   AliLatexTable tableEvts (7,"c|ccc|ccc");
327   //tableEvts.InsertCustomRow("& \\multicolumn{3}{c}{Before Phys. Sel.} & \\multicolumn{3}{c}{After Phys. Sel.}\\\\");
328   tableEvts.InsertCustomRow("Run & CINT1B & CEMC1B & CSH1B &  CINT1B & CEMC1B & CSH1B \\\\"); 
329  
330
331
332
333   fou = new TFile("qaLHC10apass2.root","recreate");
334   fou->cd();
335
336   for(Int_t igraph = 0; igraph < kNGraphs; igraph++){
337     //    TAxis * ax = graph[igraph]->GetHistogram()->GetXaxis();
338     graph[igraph]->GetXaxis()->SetTitle("run");
339     graph[igraph]->GetYaxis()->SetTitle(ylabels[igraph]);
340
341     graph[igraph]->GetXaxis()->SetLabelSize(0.06);
342     graph[igraph]->GetXaxis()->SetTitleSize(0.07);
343     graph[igraph]->GetXaxis()->SetTitleOffset(0.5);
344     graph[igraph]->GetYaxis()->SetLabelSize(0.08);
345     graph[igraph]->GetYaxis()->SetTitleSize(0.06);
346     graph[igraph]->GetYaxis()->SetTitleOffset(0.7);
347     
348     
349     graph[igraph]->GetYaxis()->CenterTitle(0);
350     graph[igraph]->GetYaxis()->SetLabelFont(42);
351     graph[igraph]->GetYaxis()->SetTitleFont(42);
352     graph[igraph]->GetYaxis()->SetNoExponent(kTRUE);
353     graph[igraph]->GetXaxis()->SetLabelFont(42);
354     graph[igraph]->GetXaxis()->SetTitleFont(52);
355
356     graph[igraph]->SetMarkerStyle(30);
357     graph[igraph]->SetMarkerColor(4);
358
359
360     graph[igraph]->SetMinimum(0);
361     graph[igraph]->SetName(gnames[igraph]);
362     graph[igraph]->SetTitle(ylabels[igraph]);
363     //    graph[igraph]->SetMaximum(1);
364     graph[igraph]->Write();
365   }
366
367  
368  
369   fou->Close();
370
371
372   
373   Int_t npoint = graph[kGraphNevAC2] -> GetN();
374   for(Int_t ipoint = 0; ipoint < npoint; ipoint++){      
375     tableEvts.SetNextCol(TMath::Nint(graph[kGraphNevAC2]->GetX()[ipoint]));
376     tableEvts.SetNextCol(TMath::Nint(graph[kGraphNevACBefore]->GetY()[ipoint]));
377     tableEvts.SetNextCol(TMath::Nint(graph[kGraphNevAC2Before]->GetY()[ipoint]));
378     tableEvts.SetNextCol(TMath::Nint(graph[kGraphNevSMHBefore]->GetY()[ipoint]));
379     tableEvts.SetNextCol(TMath::Nint(graph[kGraphNevAC]->GetY()[ipoint]));
380     tableEvts.SetNextCol(TMath::Nint(graph[kGraphNevAC2]->GetY()[ipoint]));
381     tableEvts.SetNextCol(TMath::Nint(graph[kGraphNevSMH]->GetY()[ipoint]));
382     //tableEvts.SetNextCol(TMath::Nint(graph[kGraphNev90]->GetY()[ipoint]));
383     tableEvts.InsertRow();
384     //   ax->SetBinLabel(ax->FindBin(ipoint), Form("%d", runs[ipoint-1]));
385   }    
386
387
388
389
390   Float_t meanDummy;
391   
392
393   /*
394     THREE TRIGGERS IN THE SAME CANVAS
395     
396   */
397   c4a = new TCanvas("GraphACCOverAll","GraphACCOverAll",1000,700);
398   
399
400   TPad * pad1=new TPad("pad1","pad1",0.01,0.635,0.7,0.94,0);
401   TPad * pad2=new TPad("pad2","pad2",0.01,0.33,0.7,0.635,0);
402   TPad * pad3=new TPad("pad3","pad3",0.01,0.01,0.7,0.33,0);
403
404   TPad * pad1a=new TPad("pad1a","pad1a",0.7,0.635,0.99,0.94,0);
405   TPad * pad2a=new TPad("pad2a","pad2a",0.7,0.33,0.99,0.635,0);
406   TPad * pad3a=new TPad("pad3a","pad3a",0.7,0.01,0.99,0.33,0);
407
408   pad1->Draw();
409   pad2->Draw();
410   pad3->Draw();
411   pad1a->Draw();
412   pad2a->Draw();
413   pad3a->Draw();
414
415
416   pad1->SetBottomMargin(0);
417   pad1->SetBorderSize(0);
418   pad1->SetRightMargin(0.01);
419
420   pad2->SetBottomMargin(0.0);
421   pad2->SetTopMargin(0);
422   pad2->SetRightMargin(0.01);
423   pad2->SetBorderSize(0);
424
425   pad3->SetBottomMargin(0.2);
426   pad3->SetTopMargin(0);
427   pad3->SetRightMargin(0.01);
428   pad3->SetBorderSize(0);
429
430   pad1a->SetBottomMargin(0);
431   pad1a->SetBorderSize(0);
432   pad1a->SetRightMargin(0.01);
433
434   pad2a->SetBottomMargin(0.0);
435   pad2a->SetTopMargin(0);
436   pad2a->SetRightMargin(0.01);
437   pad2a->SetBorderSize(0);
438
439   pad3a->SetBottomMargin(0.2);
440   pad3a->SetTopMargin(0);
441   pad3a->SetRightMargin(0.01);
442   pad3a->SetBorderSize(0);
443
444
445   // Draw a global picture title
446   TPaveLabel *title = new TPaveLabel(0.01,0.95,0.99,0.99,
447                                      Form("%s",period));
448
449   title->SetFillColor(0);
450   title->SetTextFont(52);
451   title->SetTextColor(4);
452   title->SetTextSize(0.7);
453   title->Draw();
454
455
456   c4a->cd();
457   pad1->Draw();
458
459
460   pad1->cd();
461
462
463   graph[kGraphACCOverAllAC]->Draw("PA");
464   graph[kGraphACCOverAllACRej] = GetGraphRej(graph[kGraphACCOverAllAC] , listRejectedRuns, "Acc/All AC" ,  meanDummy, 1);
465
466   pad1a->cd();
467   pad1a->cd()->SetGridx(1);
468   fAccOverAll[0]->Draw();
469
470
471   pad2->cd();
472
473
474   graph[kGraphACCOverAllACS2]->Draw("PA");
475   graph[kGraphACCOverAllACS2Rej] = GetGraphRej(graph[kGraphACCOverAllACS2] , listRejectedRuns, "Acc/All [ACS2]" ,  meanDummy, 1);
476
477   pad2a->cd();
478   pad2a->cd()->SetGridx(1);
479   //  GetHisto(graph[kGraphACCOverAllACS2],fAccOverAll[1]);
480   fAccOverAll[1]->Draw();
481
482
483
484   pad3->cd();
485
486
487   graph[kGraphACCOverAllSMH]->Draw("PA");
488   graph[kGraphACCOverAllSMHRej] = GetGraphRej(graph[kGraphACCOverAllSMH] , listRejectedRuns, "Acc/All SMH" ,  meanDummy, 1);
489
490   //fAccOverAll[i]
491
492   pad3a->cd();
493   pad3a->cd()->SetGridx(1);
494   //  GetHisto(graph[kGraphACCOverAllSMH],fAccOverAll[2]);
495   fAccOverAll[2]->Draw();
496   fAccOverAll[2]->GetXaxis()->SetTitle("Accepted / All");
497
498   /*
499   for(Int_t ipoint = 0; ipoint < graph[kGraphACCOverAllAC]->GetN(); ipoint++){
500     fhistTest->Fill(graph[kGraphACCOverAllAC]->GetY()[ipoint]);
501   }
502
503   pad1->cd(4);
504   fhistTest->Sumw2();
505   fhistTest->SetMarkerStyle(25);
506   fhistTest->GetXaxis()->SetTitle(ylabels[kGraphACCOverAllAC]);
507   fhistTest->GetYaxis()->SetTitle("Entries");
508   fhistTest->Draw();
509   */
510   //pad1->cd(2);
511
512
513
514
515
516   c4a->Update();
517   gSystem->ProcessEvents();
518   c4a->SaveAs(Form("picturesLHC11hAOD50/c4a_%s.png",c4a->GetName()));
519
520
521
522   c5a = new TCanvas("GraphNev","GraphNev",1000,700);
523   
524   TPad *pad1 = new TPad("pad1",
525                         "The pad with the function",0.01,0.01,0.99,0.94,0);
526   pad1->Draw();
527
528
529   // Draw a global picture title
530   TPaveLabel *title = new TPaveLabel(0.01,0.95,0.99,0.99,
531                                      Form("%s",period));
532
533   title->SetFillColor(0);
534   title->SetTextFont(52);
535   title->SetTextColor(4);
536   title->SetTextSize(0.7);
537   title->Draw();
538
539
540
541
542   pad1->Divide(1,3);
543   pad1->cd(1);
544   graph[kGraphNevAC]->Draw("PA");
545
546   pad1->cd(2);
547   graph[kGraphNevAC2]->Draw("PA");
548
549
550   pad1->cd(3);
551   graph[kGraphNevSMH]->Draw("PA");
552  
553
554
555   c5a->Update();
556   gSystem->ProcessEvents();
557   c5a->SaveAs(Form("picturesLHC11hAOD50/c5a_%s.png",c5a->GetName()));
558
559
560
561
562
563   c7a = new TCanvas("GraphBGOverAll","GraphACCOverAll",1000,700);
564   
565
566   TPad * pad1=new TPad("pad1","pad1",0.01,0.635,0.7,0.94,0);
567   TPad * pad2=new TPad("pad2","pad2",0.01,0.33,0.7,0.635,0);
568   TPad * pad3=new TPad("pad3","pad3",0.01,0.01,0.7,0.33,0);
569
570   TPad * pad1a=new TPad("pad1a","pad1a",0.7,0.635,0.99,0.94,0);
571   TPad * pad2a=new TPad("pad2a","pad2a",0.7,0.33,0.99,0.635,0);
572   TPad * pad3a=new TPad("pad3a","pad3a",0.7,0.01,0.99,0.33,0);
573
574   pad1->Draw();
575   pad2->Draw();
576   pad3->Draw();
577   pad1a->Draw();
578   pad2a->Draw();
579   pad3a->Draw();
580
581
582   pad1->SetBottomMargin(0);
583   pad1->SetBorderSize(0);
584   pad1->SetRightMargin(0.01);
585
586   pad2->SetBottomMargin(0.0);
587   pad2->SetTopMargin(0);
588   pad2->SetRightMargin(0.01);
589   pad2->SetBorderSize(0);
590
591   pad3->SetBottomMargin(0.2);
592   pad3->SetTopMargin(0);
593   pad3->SetRightMargin(0.01);
594   pad3->SetBorderSize(0);
595
596   pad1a->SetBottomMargin(0);
597   pad1a->SetBorderSize(0);
598   pad1a->SetRightMargin(0.01);
599
600   pad2a->SetBottomMargin(0.0);
601   pad2a->SetTopMargin(0);
602   pad2a->SetRightMargin(0.01);
603   pad2a->SetBorderSize(0);
604
605   pad3a->SetBottomMargin(0.2);
606   pad3a->SetTopMargin(0);
607   pad3a->SetRightMargin(0.01);
608   pad3a->SetBorderSize(0);
609
610
611   // Draw a global picture title
612   TPaveLabel *title = new TPaveLabel(0.01,0.95,0.99,0.99,
613                                      Form("%s",period));
614
615   title->SetFillColor(0);
616   title->SetTextFont(52);
617   title->SetTextColor(4);
618   title->SetTextSize(0.7);
619   title->Draw();
620
621
622   c7a->cd();
623   pad1->Draw();
624
625
626   pad1->cd();
627
628   graph[kGraphBGOverAllAC]->Draw("PA");
629   graph[kGraphBGOverAllAC] = GetGraphRej(graph[kGraphBGOverAllAC] , listRejectedRuns, "BG/All AC" ,  meanDummy, 1);
630
631   pad1a->cd();
632   pad1a->cd()->SetGridx(1);
633   fBGOverAll[0]->Draw();
634  
635
636
637
638
639
640
641   pad2->cd();
642
643   graph[kGraphBGOverAllAC2]->Draw("PA");
644   graph[kGraphBGOverAllAC2] = GetGraphRej(graph[kGraphBGOverAllAC2] , listRejectedRuns, "BG/All AC2" ,  meanDummy, 1);
645
646
647   pad2a->cd();
648   pad2a->cd()->SetGridx(1);
649   fBGOverAll[1]->Draw();
650
651
652
653   pad3->cd();
654
655  
656   graph[kGraphBGOverAllSMH]->Draw("PA");
657   graph[kGraphBGOverAllSMHRej] = GetGraphRej(graph[kGraphBGOverAllSMH] , listRejectedRuns, "BG/All SMH" ,  meanDummy, 1);
658
659
660   pad3a->cd();
661   pad3a->cd()->SetGridx(1);
662
663   fBGOverAll[2]->GetXaxis()->SetTitle("BG / All");
664   fBGOverAll[2]->Draw();
665
666
667
668
669   c7a->Update();
670   gSystem->ProcessEvents();
671   c7a->SaveAs(Form("picturesLHC11hAOD50/c7a_%s.png",c7a->GetName()));
672
673
674
675   c7b = new TCanvas("GraphV0BGOverAll","GraphV0BGOverAll",1000,700);
676     TPad * pad1=new TPad("pad1","pad1",0.01,0.635,0.7,0.94,0);
677   TPad * pad2=new TPad("pad2","pad2",0.01,0.33,0.7,0.635,0);
678   TPad * pad3=new TPad("pad3","pad3",0.01,0.01,0.7,0.33,0);
679
680   TPad * pad1a=new TPad("pad1a","pad1a",0.7,0.635,0.99,0.94,0);
681   TPad * pad2a=new TPad("pad2a","pad2a",0.7,0.33,0.99,0.635,0);
682   TPad * pad3a=new TPad("pad3a","pad3a",0.7,0.01,0.99,0.33,0);
683
684   pad1->Draw();
685   pad2->Draw();
686   pad3->Draw();
687   pad1a->Draw();
688   pad2a->Draw();
689   pad3a->Draw();
690
691
692   pad1->SetBottomMargin(0);
693   pad1->SetBorderSize(0);
694   pad1->SetRightMargin(0.01);
695
696   pad2->SetBottomMargin(0.0);
697   pad2->SetTopMargin(0);
698   pad2->SetRightMargin(0.01);
699   pad2->SetBorderSize(0);
700
701   pad3->SetBottomMargin(0.2);
702   pad3->SetTopMargin(0);
703   pad3->SetRightMargin(0.01);
704   pad3->SetBorderSize(0);
705
706   pad1a->SetBottomMargin(0);
707   pad1a->SetBorderSize(0);
708   pad1a->SetRightMargin(0.01);
709
710   pad2a->SetBottomMargin(0.0);
711   pad2a->SetTopMargin(0);
712   pad2a->SetRightMargin(0.01);
713   pad2a->SetBorderSize(0);
714
715   pad3a->SetBottomMargin(0.2);
716   pad3a->SetTopMargin(0);
717   pad3a->SetRightMargin(0.01);
718   pad3a->SetBorderSize(0);
719
720
721   // Draw a global picture title
722   TPaveLabel *title = new TPaveLabel(0.01,0.95,0.99,0.99,
723                                      Form("%s",period));
724
725   title->SetFillColor(0);
726   title->SetTextFont(52);
727   title->SetTextColor(4);
728   title->SetTextSize(0.7);
729   title->Draw();
730
731
732   c7b->cd();
733
734
735
736   pad1->cd();
737
738   graph[kGraphV0BGOverAllAC]->Draw("PA");
739   graph[kGraphV0BGOverAllAC] = GetGraphRej(graph[kGraphV0BGOverAllAC] , listRejectedRuns, "V0BG/All AC" ,  meanDummy, 1);
740
741   pad1a->cd();
742   pad1a->cd()->SetGridx(1);
743   fV0BGOverAll[0]->Draw();
744
745
746   pad2->cd();
747   graph[kGraphV0BGOverAllAC2]->Draw("PA");
748   graph[kGraphV0BGOverAllAC2] = GetGraphRej(graph[kGraphV0BGOverAllAC2] , listRejectedRuns, "V0BG/All AC2" ,  meanDummy, 1);
749
750   pad2a->cd();
751   pad2a->cd()->SetGridx(1);
752   fV0BGOverAll[1]->Draw();
753
754   pad3->cd();
755   graph[kGraphV0BGOverAllSMH]->Draw("PA");
756   graph[kGraphV0BGOverAllSMHRej] = GetGraphRej(graph[kGraphV0BGOverAllSMH] , listRejectedRuns, "V0BG/All SMH" ,  meanDummy, 1);
757
758   pad3a->cd();
759   pad3a->cd()->SetGridx(1);
760   fV0BGOverAll[2]->GetXaxis()->SetTitle("V0BG / All");
761   fV0BGOverAll[2]->Draw();
762
763
764
765   c7b->Update();
766   gSystem->ProcessEvents();
767   c7b->SaveAs(Form("picturesLHC11hAOD50/c7b_%s.png",c7b->GetName()));
768
769
770
771
772
773
774
775
776
777
778   c7c = new TCanvas("GraphV0A_BGOverAll","GraphV0A_BGOverAll",1000,700);
779    TPad * pad1=new TPad("pad1","pad1",0.01,0.635,0.7,0.94,0);
780   TPad * pad2=new TPad("pad2","pad2",0.01,0.33,0.7,0.635,0);
781   TPad * pad3=new TPad("pad3","pad3",0.01,0.01,0.7,0.33,0);
782
783   TPad * pad1a=new TPad("pad1a","pad1a",0.7,0.635,0.99,0.94,0);
784   TPad * pad2a=new TPad("pad2a","pad2a",0.7,0.33,0.99,0.635,0);
785   TPad * pad3a=new TPad("pad3a","pad3a",0.7,0.01,0.99,0.33,0);
786
787   pad1->Draw();
788   pad2->Draw();
789   pad3->Draw();
790   pad1a->Draw();
791   pad2a->Draw();
792   pad3a->Draw();
793
794
795   pad1->SetBottomMargin(0);
796   pad1->SetBorderSize(0);
797   pad1->SetRightMargin(0.01);
798
799   pad2->SetBottomMargin(0.0);
800   pad2->SetTopMargin(0);
801   pad2->SetRightMargin(0.01);
802   pad2->SetBorderSize(0);
803
804   pad3->SetBottomMargin(0.2);
805   pad3->SetTopMargin(0);
806   pad3->SetRightMargin(0.01);
807   pad3->SetBorderSize(0);
808
809   pad1a->SetBottomMargin(0);
810   pad1a->SetBorderSize(0);
811   pad1a->SetRightMargin(0.01);
812
813   pad2a->SetBottomMargin(0.0);
814   pad2a->SetTopMargin(0);
815   pad2a->SetRightMargin(0.01);
816   pad2a->SetBorderSize(0);
817
818   pad3a->SetBottomMargin(0.2);
819   pad3a->SetTopMargin(0);
820   pad3a->SetRightMargin(0.01);
821   pad3a->SetBorderSize(0);
822
823
824   // Draw a global picture title
825   TPaveLabel *title = new TPaveLabel(0.01,0.95,0.99,0.99,
826                                      Form("%s",period));
827
828   title->SetFillColor(0);
829   title->SetTextFont(52);
830   title->SetTextColor(4);
831   title->SetTextSize(0.7);
832   title->Draw();
833
834
835   c7c->cd();
836
837   pad1->cd(); 
838
839
840   graph[kGraphV0ABGOverAllAC]->Draw("PA");
841   graph[kGraphV0ABGOverAllAC] = GetGraphRej(graph[kGraphV0ABGOverAllAC] , listRejectedRuns, "V0A_BG/All AC" ,  meanDummy, 1);
842
843   pad1a->cd();
844   pad1a->cd()->SetGridx(1);
845   fV0ABGOverAll[0]->Draw();
846
847
848
849   pad2->cd(); 
850   graph[kGraphV0ABGOverAllAC2]->Draw("PA");
851   graph[kGraphV0ABGOverAllAC2] = GetGraphRej(graph[kGraphV0ABGOverAllAC2] , listRejectedRuns, "V0A_BG/All AC2" ,  meanDummy, 1);
852
853   pad2a->cd();
854   pad2a->cd()->SetGridx(1);
855   fV0ABGOverAll[1]->Draw();
856
857   pad3->cd(); 
858   graph[kGraphV0ABGOverAllSMH]->Draw("PA");
859   graph[kGraphV0ABGOverAllSMHRej] = GetGraphRej(graph[kGraphV0ABGOverAllSMH] , listRejectedRuns, "V0A_BG/All SMH" ,  meanDummy, 1);
860
861
862   pad3a->cd();
863   pad3a->cd()->SetGridx(1);
864   fV0ABGOverAll[2]->Draw();
865   fV0ABGOverAll[2]->GetXaxis()->SetTitle("V0ABG / All");
866
867
868   c7c->Update();
869   gSystem->ProcessEvents();
870   c7c->SaveAs(Form("picturesLHC11hAOD50/c7c_%s.png",c7c->GetName()));
871
872
873
874
875
876   c7d = new TCanvas("GraphV0C_BGOverAll","GraphV0C_BGOverAll",1000,700);
877   TPad * pad1=new TPad("pad1","pad1",0.01,0.635,0.7,0.94,0);
878   TPad * pad2=new TPad("pad2","pad2",0.01,0.33,0.7,0.635,0);
879   TPad * pad3=new TPad("pad3","pad3",0.01,0.01,0.7,0.33,0);
880
881   TPad * pad1a=new TPad("pad1a","pad1a",0.7,0.635,0.99,0.94,0);
882   TPad * pad2a=new TPad("pad2a","pad2a",0.7,0.33,0.99,0.635,0);
883   TPad * pad3a=new TPad("pad3a","pad3a",0.7,0.01,0.99,0.33,0);
884
885   pad1->Draw();
886   pad2->Draw();
887   pad3->Draw();
888   pad1a->Draw();
889   pad2a->Draw();
890   pad3a->Draw();
891
892
893   pad1->SetBottomMargin(0);
894   pad1->SetBorderSize(0);
895   pad1->SetRightMargin(0.01);
896
897   pad2->SetBottomMargin(0.0);
898   pad2->SetTopMargin(0);
899   pad2->SetRightMargin(0.01);
900   pad2->SetBorderSize(0);
901
902   pad3->SetBottomMargin(0.2);
903   pad3->SetTopMargin(0);
904   pad3->SetRightMargin(0.01);
905   pad3->SetBorderSize(0);
906
907   pad1a->SetBottomMargin(0);
908   pad1a->SetBorderSize(0);
909   pad1a->SetRightMargin(0.01);
910
911   pad2a->SetBottomMargin(0.0);
912   pad2a->SetTopMargin(0);
913   pad2a->SetRightMargin(0.01);
914   pad2a->SetBorderSize(0);
915
916   pad3a->SetBottomMargin(0.2);
917   pad3a->SetTopMargin(0);
918   pad3a->SetRightMargin(0.01);
919   pad3a->SetBorderSize(0);
920
921
922   // Draw a global picture title
923   TPaveLabel *title = new TPaveLabel(0.01,0.95,0.99,0.99,
924                                      Form("%s",period));
925
926   title->SetFillColor(0);
927   title->SetTextFont(52);
928   title->SetTextColor(4);
929   title->SetTextSize(0.7);
930   title->Draw();
931
932
933   c7d->cd();
934
935   pad1->cd();
936   graph[kGraphV0CBGOverAllAC]->Draw("PA");
937   graph[kGraphV0CBGOverAllAC] = GetGraphRej(graph[kGraphV0CBGOverAllAC] , listRejectedRuns, "V0C_BG/All AC" ,  meanDummy, 1);
938
939   pad1a->cd();
940   pad1a->cd()->SetGridx(1);
941   fV0CBGOverAll[0]->Draw();
942
943
944   pad2->cd();
945   graph[kGraphV0CBGOverAllAC2]->Draw("PA");
946   graph[kGraphV0CBGOverAllAC2] = GetGraphRej(graph[kGraphV0CBGOverAllAC2] , listRejectedRuns, "V0C_BG/All AC2" ,  meanDummy, 1);
947
948   pad2a->cd();
949   pad2a->cd()->SetGridx(1);
950   fV0CBGOverAll[1]->Draw();
951
952
953            pad3->cd();
954   graph[kGraphV0CBGOverAllSMH]->Draw("PA");
955   graph[kGraphV0CBGOverAllSMHRej] = GetGraphRej(graph[kGraphV0CBGOverAllSMH] , listRejectedRuns, "V0A_BG/All SMH" ,  meanDummy, 1);
956
957
958   pad3a->cd();
959   pad3a->cd()->SetGridx(1);
960   fV0CBGOverAll[2]->Draw();
961   fV0CBGOverAll[2]->GetXaxis()->SetTitle("V0CBG / All");
962
963
964
965   c7d->Update();
966   gSystem->ProcessEvents();
967   c7d->SaveAs(Form("picturesLHC11hAOD50/c7d_%s.png",c7d->GetName()));
968
969
970
971
972
973   c8a = new TCanvas("GraphAfterOverBefore","GraphAfterOverBefore",1000,700);
974   TPad * pad1=new TPad("pad1","pad1",0.01,0.635,0.7,0.94,0);
975   TPad * pad2=new TPad("pad2","pad2",0.01,0.33,0.7,0.635,0);
976   TPad * pad3=new TPad("pad3","pad3",0.01,0.01,0.7,0.33,0);
977
978   TPad * pad1a=new TPad("pad1a","pad1a",0.7,0.635,0.99,0.94,0);
979   TPad * pad2a=new TPad("pad2a","pad2a",0.7,0.33,0.99,0.635,0);
980   TPad * pad3a=new TPad("pad3a","pad3a",0.7,0.01,0.99,0.33,0);
981
982   pad1->Draw();
983   pad2->Draw();
984   pad3->Draw();
985   pad1a->Draw();
986   pad2a->Draw();
987   pad3a->Draw();
988
989
990   pad1->SetBottomMargin(0);
991   pad1->SetBorderSize(0);
992   pad1->SetRightMargin(0.01);
993
994   pad2->SetBottomMargin(0.0);
995   pad2->SetTopMargin(0);
996   pad2->SetRightMargin(0.01);
997   pad2->SetBorderSize(0);
998
999   pad3->SetBottomMargin(0.2);
1000   pad3->SetTopMargin(0);
1001   pad3->SetRightMargin(0.01);
1002   pad3->SetBorderSize(0);
1003
1004   pad1a->SetBottomMargin(0);
1005   pad1a->SetBorderSize(0);
1006   pad1a->SetRightMargin(0.01);
1007
1008   pad2a->SetBottomMargin(0.0);
1009   pad2a->SetTopMargin(0);
1010   pad2a->SetRightMargin(0.01);
1011   pad2a->SetBorderSize(0);
1012
1013   pad3a->SetBottomMargin(0.2);
1014   pad3a->SetTopMargin(0);
1015   pad3a->SetRightMargin(0.01);
1016   pad3a->SetBorderSize(0);
1017
1018
1019   // Draw a global picture title
1020   TPaveLabel *title = new TPaveLabel(0.01,0.95,0.99,0.99,
1021                                      Form("%s",period));
1022
1023   title->SetFillColor(0);
1024   title->SetTextFont(52);
1025   title->SetTextColor(4);
1026   title->SetTextSize(0.7);
1027   title->Draw();
1028   
1029   pad1->cd();
1030   
1031   graph[kGraphNevACratioAfter]->Draw("PA");
1032   graph[kGraphNevACratioAfterRej] = GetGraphRej(graph[kGraphNevACratioAfter] , listRejectedRuns, "V0BG/All AC" ,  meanDummy, 1);
1033   
1034   pad1a->cd();
1035   pad1a->cd()->SetGridx(1);
1036   fAfterOverBefore[0]->Draw();
1037
1038
1039   pad2->cd();
1040   graph[kGraphNevAC2ratioAfter]->Draw("PA");
1041   graph[kGraphNevAC2ratioAfterRej] = GetGraphRej(graph[kGraphNevAC2ratioAfter] , listRejectedRuns, "V0BG/All AC2" ,  meanDummy, 1);
1042
1043   pad2a->cd();
1044   pad2a->cd()->SetGridx(1);
1045   fAfterOverBefore[1]->Draw();
1046
1047   pad3->cd();
1048   graph[kGraphNevSMHratioAfter]->Draw("PA");
1049   graph[kGraphNevSMHratioAfterRej] = GetGraphRej(graph[kGraphNevSMHratioAfter] , listRejectedRuns, "V0BG/All SMH" ,  meanDummy, 1);
1050
1051   pad3a->cd();
1052   pad3a->cd()->SetGridx(1);
1053   fAfterOverBefore[2]->Draw();
1054   fAfterOverBefore[2]->GetXaxis()->SetTitle("After_PS / Before_PS");
1055
1056
1057   c8a->Update();
1058   gSystem->ProcessEvents();
1059   c8a->SaveAs(Form("picturesLHC11hAOD50/c8a_%s.png",c8a->GetName()));
1060
1061
1062
1063
1064
1065   c8b = new TCanvas("GraphNevBefore","GraphNevBefore",1000,700);
1066   
1067   TPad *pad1 = new TPad("pad1",
1068                         "The pad with the function",0.01,0.01,0.99,0.94,0);
1069   pad1->Draw();
1070
1071
1072   // Draw a global picture title
1073   TPaveLabel *title = new TPaveLabel(0.01,0.95,0.99,0.99,
1074                                      Form("%s",period));
1075
1076   title->SetFillColor(0);
1077   title->SetTextFont(52);
1078   title->SetTextColor(4);
1079   title->SetTextSize(0.7);
1080   title->Draw();
1081
1082
1083
1084
1085   pad1->Divide(1,3);
1086
1087   pad1->cd(1);
1088   graph[kGraphNevACBefore]->Draw("PA");
1089
1090   pad1->cd(2);
1091   graph[kGraphNevAC2Before]->Draw("PA");
1092
1093   pad1->cd(3);
1094   graph[kGraphNevSMHBefore]->Draw("PA");
1095
1096
1097
1098   c8b->Update();
1099   gSystem->ProcessEvents();
1100   c8b->SaveAs(Form("picturesLHC11hAOD50/c8b_%s.png",c8b->GetName()));
1101
1102
1103
1104
1105
1106   c8c = new TCanvas("GraphNevAfter","GraphNevAfter",1000,700);
1107   
1108   TPad *pad1 = new TPad("pad1",
1109                         "The pad with the function",0.01,0.01,0.99,0.94,0);
1110   pad1->Draw();
1111
1112
1113   // Draw a global picture title
1114   TPaveLabel *title = new TPaveLabel(0.01,0.95,0.99,0.99,
1115                                      Form("%s",period));
1116
1117   title->SetFillColor(0);
1118   title->SetTextFont(52);
1119   title->SetTextColor(4);
1120   title->SetTextSize(0.7);
1121   title->Draw();
1122
1123
1124
1125
1126   pad1->Divide(1,3);
1127
1128   pad1->cd(1);
1129   graph[kGraphNevACAfter]->Draw("PA");
1130
1131
1132   pad1->cd(2);
1133   graph[kGraphNevAC2After]->Draw("PA");
1134
1135
1136   pad1->cd(3);
1137
1138   graph[kGraphNevSMHAfter]->Draw("PA");
1139
1140
1141
1142   c8c->Update();
1143   gSystem->ProcessEvents();
1144   c8c->SaveAs(Form("picturesLHC11hAOD50/c8c_%s.png",c8c->GetName()));
1145
1146
1147
1148
1149
1150
1151
1152
1153   c9a = new TCanvas("GraphV0AOverV0C","GraphV0AOverV0C",1000,700);
1154   TPad * pad1=new TPad("pad1","pad1",0.01,0.635,0.7,0.94,0);
1155   TPad * pad2=new TPad("pad2","pad2",0.01,0.33,0.7,0.635,0);
1156   TPad * pad3=new TPad("pad3","pad3",0.01,0.01,0.7,0.33,0);
1157
1158   TPad * pad1a=new TPad("pad1a","pad1a",0.7,0.635,0.99,0.94,0);
1159   TPad * pad2a=new TPad("pad2a","pad2a",0.7,0.33,0.99,0.635,0);
1160   TPad * pad3a=new TPad("pad3a","pad3a",0.7,0.01,0.99,0.33,0);
1161
1162   pad1->Draw();
1163   pad2->Draw();
1164   pad3->Draw();
1165   pad1a->Draw();
1166   pad2a->Draw();
1167   pad3a->Draw();
1168
1169
1170   pad1->SetBottomMargin(0);
1171   pad1->SetBorderSize(0);
1172   pad1->SetRightMargin(0.01);
1173
1174   pad2->SetBottomMargin(0.0);
1175   pad2->SetTopMargin(0);
1176   pad2->SetRightMargin(0.01);
1177   pad2->SetBorderSize(0);
1178
1179   pad3->SetBottomMargin(0.2);
1180   pad3->SetTopMargin(0);
1181   pad3->SetRightMargin(0.01);
1182   pad3->SetBorderSize(0);
1183
1184   pad1a->SetBottomMargin(0);
1185   pad1a->SetBorderSize(0);
1186   pad1a->SetRightMargin(0.01);
1187
1188   pad2a->SetBottomMargin(0.0);
1189   pad2a->SetTopMargin(0);
1190   pad2a->SetRightMargin(0.01);
1191   pad2a->SetBorderSize(0);
1192
1193   pad3a->SetBottomMargin(0.2);
1194   pad3a->SetTopMargin(0);
1195   pad3a->SetRightMargin(0.01);
1196   pad3a->SetBorderSize(0);
1197
1198
1199   // Draw a global picture title
1200   TPaveLabel *title = new TPaveLabel(0.01,0.95,0.99,0.99,
1201                                      Form("%s",period));
1202
1203   title->SetFillColor(0);
1204   title->SetTextFont(52);
1205   title->SetTextColor(4);
1206   title->SetTextSize(0.7);
1207   title->Draw();
1208
1209   pad1->cd();
1210   graph[kGraphV0AOverV0CAC]->Draw("PA");
1211   graph[kGraphV0AOverV0CACRej] = GetGraphRej(graph[kGraphV0AOverV0CAC] , listRejectedRuns, "Nev_V0A/Nev_V0C AC" ,  meanDummy, 1);
1212
1213   pad1a->cd();
1214   pad1a->cd()->SetGridx(1);
1215   fV0AOverV0C[0]->Draw();
1216
1217
1218   pad2->cd();
1219   graph[kGraphV0AOverV0CAC2]->Draw("PA");
1220   graph[kGraphV0AOverV0CAC2Rej] = GetGraphRej(graph[kGraphV0AOverV0CAC2] , listRejectedRuns, "Nev_V0A/Nev_V0C AC2" ,  meanDummy, 1);
1221
1222   pad2a->cd();
1223   pad2a->cd()->SetGridx(1);
1224   fV0AOverV0C[1]->Draw();
1225
1226
1227   pad3->cd();
1228   graph[kGraphV0AOverV0CSMH]->Draw("PA");
1229   graph[kGraphV0AOverV0CSMHRej] = GetGraphRej(graph[kGraphV0AOverV0CSMH] , listRejectedRuns, "Nev_V0A/Nev_V0C SMH" ,  meanDummy, 1);
1230
1231   pad3a->cd();
1232   pad3a->cd()->SetGridx(1);
1233   fV0AOverV0C[2]->Draw();
1234   fV0AOverV0C[2]->GetXaxis()->SetTitle("V0A / VOC");
1235
1236
1237
1238   c9a->Update();
1239   gSystem->ProcessEvents();
1240   c9a->SaveAs(Form("picturesLHC11hAOD50/c9a_%s.png",c9a->GetName()));
1241
1242
1243
1244   c10a = new TCanvas("GraphFO1OverAll","GraphFO1OverAll",1000,700);
1245   TPad * pad1=new TPad("pad1","pad1",0.01,0.635,0.7,0.94,0);
1246   TPad * pad2=new TPad("pad2","pad2",0.01,0.33,0.7,0.635,0);
1247   TPad * pad3=new TPad("pad3","pad3",0.01,0.01,0.7,0.33,0);
1248
1249   TPad * pad1a=new TPad("pad1a","pad1a",0.7,0.635,0.99,0.94,0);
1250   TPad * pad2a=new TPad("pad2a","pad2a",0.7,0.33,0.99,0.635,0);
1251   TPad * pad3a=new TPad("pad3a","pad3a",0.7,0.01,0.99,0.33,0);
1252
1253   pad1->Draw();
1254   pad2->Draw();
1255   pad3->Draw();
1256   pad1a->Draw();
1257   pad2a->Draw();
1258   pad3a->Draw();
1259
1260
1261   pad1->SetBottomMargin(0);
1262   pad1->SetBorderSize(0);
1263   pad1->SetRightMargin(0.01);
1264
1265   pad2->SetBottomMargin(0.0);
1266   pad2->SetTopMargin(0);
1267   pad2->SetRightMargin(0.01);
1268   pad2->SetBorderSize(0);
1269
1270   pad3->SetBottomMargin(0.2);
1271   pad3->SetTopMargin(0);
1272   pad3->SetRightMargin(0.01);
1273   pad3->SetBorderSize(0);
1274
1275   pad1a->SetBottomMargin(0);
1276   pad1a->SetBorderSize(0);
1277   pad1a->SetRightMargin(0.01);
1278
1279   pad2a->SetBottomMargin(0.0);
1280   pad2a->SetTopMargin(0);
1281   pad2a->SetRightMargin(0.01);
1282   pad2a->SetBorderSize(0);
1283
1284   pad3a->SetBottomMargin(0.2);
1285   pad3a->SetTopMargin(0);
1286   pad3a->SetRightMargin(0.01);
1287   pad3a->SetBorderSize(0);
1288
1289
1290   // Draw a global picture title
1291   TPaveLabel *title = new TPaveLabel(0.01,0.95,0.99,0.99,
1292                                      Form("%s",period));
1293
1294   title->SetFillColor(0);
1295   title->SetTextFont(52);
1296   title->SetTextColor(4);
1297   title->SetTextSize(0.7);
1298   title->Draw();  
1299
1300
1301   pad1->cd();
1302   graph[kGraphFO1OverAllAC]->Draw("PA");
1303   graph[kGraphFO1OverAllACRej] = GetGraphRej(graph[kGraphFO1OverAllAC] , listRejectedRuns, "Nev_FO1/All AC" ,  meanDummy, 1);
1304
1305   pad1a->cd();
1306   pad1a->cd()->SetGridx(1);
1307   fF01OverAll[0]->Draw();
1308
1309
1310   pad2->cd();
1311   graph[kGraphFO1OverAllAC2]->Draw("PA");
1312   graph[kGraphFO1OverAllAC2Rej] = GetGraphRej(graph[kGraphFO1OverAllAC2] , listRejectedRuns, "Nev_FO1/All AC2" ,  meanDummy, 1);
1313
1314   pad2a->cd();
1315   pad2a->cd()->SetGridx(1);
1316   fF01OverAll[1]->Draw();
1317
1318
1319   pad3->cd();
1320   graph[kGraphFO1OverAllSMH]->Draw("PA");
1321   graph[kGraphFO1OverAllSMHRej] = GetGraphRej(graph[kGraphFO1OverAllSMH] , listRejectedRuns, "Nev_FO1/All SMH" ,  meanDummy, 1);
1322
1323   pad3a->cd();
1324   pad3a->cd()->SetGridx(1);
1325   fF01OverAll[2]->Draw();
1326   fF01OverAll[2]->GetXaxis()->SetTitle("F01 / All");
1327
1328
1329
1330
1331   c10a->Update();
1332   gSystem->ProcessEvents();
1333   c10a->SaveAs(Form("picturesLHC11hAOD50/c10a_%s.png",c10a->GetName()));
1334
1335
1336
1337
1338   c10b = new TCanvas("GraphFO2OverAll","GraphFO2OverAll",1000,700);
1339   TPad * pad1=new TPad("pad1","pad1",0.01,0.635,0.7,0.94,0);
1340   TPad * pad2=new TPad("pad2","pad2",0.01,0.33,0.7,0.635,0);
1341   TPad * pad3=new TPad("pad3","pad3",0.01,0.01,0.7,0.33,0);
1342
1343   TPad * pad1a=new TPad("pad1a","pad1a",0.7,0.635,0.99,0.94,0);
1344   TPad * pad2a=new TPad("pad2a","pad2a",0.7,0.33,0.99,0.635,0);
1345   TPad * pad3a=new TPad("pad3a","pad3a",0.7,0.01,0.99,0.33,0);
1346
1347   pad1->Draw();
1348   pad2->Draw();
1349   pad3->Draw();
1350   pad1a->Draw();
1351   pad2a->Draw();
1352   pad3a->Draw();
1353
1354
1355   pad1->SetBottomMargin(0);
1356   pad1->SetBorderSize(0);
1357   pad1->SetRightMargin(0.01);
1358
1359   pad2->SetBottomMargin(0.0);
1360   pad2->SetTopMargin(0);
1361   pad2->SetRightMargin(0.01);
1362   pad2->SetBorderSize(0);
1363
1364   pad3->SetBottomMargin(0.2);
1365   pad3->SetTopMargin(0);
1366   pad3->SetRightMargin(0.01);
1367   pad3->SetBorderSize(0);
1368
1369   pad1a->SetBottomMargin(0);
1370   pad1a->SetBorderSize(0);
1371   pad1a->SetRightMargin(0.01);
1372
1373   pad2a->SetBottomMargin(0.0);
1374   pad2a->SetTopMargin(0);
1375   pad2a->SetRightMargin(0.01);
1376   pad2a->SetBorderSize(0);
1377
1378   pad3a->SetBottomMargin(0.2);
1379   pad3a->SetTopMargin(0);
1380   pad3a->SetRightMargin(0.01);
1381   pad3a->SetBorderSize(0);
1382
1383
1384   // Draw a global picture title
1385   TPaveLabel *title = new TPaveLabel(0.01,0.95,0.99,0.99,
1386                                      Form("%s",period));
1387
1388   title->SetFillColor(0);
1389   title->SetTextFont(52);
1390   title->SetTextColor(4);
1391   title->SetTextSize(0.7);
1392   title->Draw();  
1393
1394
1395   pad1->cd();
1396   graph[kGraphFO2OverAllAC]->Draw("PA");
1397   graph[kGraphFO2OverAllACRej] = GetGraphRej(graph[kGraphFO2OverAllAC] , listRejectedRuns, "Nev_FO2/All AC" ,  meanDummy, 1);
1398
1399   pad1a->cd();
1400   pad1a->cd()->SetGridx(1);
1401   fF02OverAll[0]->Draw();
1402   fF02OverAll[0]->GetXaxis()->SetTitle("F02 / All");
1403
1404
1405   pad2->cd();
1406   graph[kGraphFO2OverAllAC2]->Draw("PA");
1407   graph[kGraphFO2OverAllAC2Rej] = GetGraphRej(graph[kGraphFO2OverAllAC2] , listRejectedRuns, "Nev_FO2/All AC2" ,  meanDummy, 1);
1408
1409   pad2a->cd();
1410   pad2a->cd()->SetGridx(1);
1411   fF02OverAll[1]->Draw();
1412   fF02OverAll[1]->GetXaxis()->SetTitle("F02 / All");
1413
1414   pad3->cd();
1415   graph[kGraphFO2OverAllSMH]->Draw("PA");
1416   graph[kGraphFO2OverAllSMHRej] = GetGraphRej(graph[kGraphFO2OverAllSMH] , listRejectedRuns, "Nev_FO2/All SMH" ,  meanDummy, 1);
1417
1418   pad3a->cd();
1419   pad3a->cd()->SetGridx(1);
1420   fF02OverAll[2]->Draw();
1421   fF02OverAll[2]->GetXaxis()->SetTitle("F02 / All");
1422
1423
1424
1425
1426
1427   c10b->Update();
1428   gSystem->ProcessEvents();
1429   c10b->SaveAs(Form("picturesLHC11hAOD50/c10b_%s.png",c10b->GetName()));
1430
1431
1432
1433   cout << "Files Statistics" << endl;
1434   cout << " Total     [" << ifile         << "]" << endl;
1435   cout << " Available [" << ifileGood     << "]" << endl;  
1436   cout << " Not Empty [" << ifileNotEmpty << "]" <<  endl;
1437
1438   cout << "" << endl;
1439   cout << "All Events" << endl <<  endl;
1440   tableEvts.PrintTable("CSV");
1441
1442   cout << "Empty or missing files" << endl<< endl;  
1443   listEmptyRuns->Sort();
1444   listEmptyRuns->Print();
1445
1446   cout << "Suspicious Runs" << endl << endl;
1447   listRejectedRuns->Sort();
1448   listRejectedRuns->Print();
1449
1450
1451  
1452
1453   
1454
1455 }
1456
1457 TH1F * GetEmpty(const char * name, Int_t nfile) {
1458   TH1F * he01 = new TH1F(TString("hempty")+name, "hempty", nfile, -0.5, nfile-0.5);
1459   for(Int_t ilab = 0; ilab < nfile; ilab++){
1460     he01->GetXaxis()->SetBinLabel(ilab+1, Form("%d", runs[ilab]));
1461   }
1462   he01->SetMinimum(0);
1463   he01->SetMaximum(1);
1464   he01->SetXTitle("run");
1465   he01->SetYTitle(name);
1466   return he01;
1467 }
1468
1469
1470 const char * GetLocalFileName(Int_t run, const char * suffix, const char * path) {
1471   // returns the filename of the local copy of the event_stat file
1472   static TString name;
1473   //  name.Form("%s/event_stat_%s_%d.root", path, suffix, run);
1474   name.Form("%s/event_stat_%d.root", path, run);
1475   return name.Data();
1476
1477 }
1478
1479 double meanMed(double* vec, int np, double nsigmaCut, int &nrej, int *rejList)
1480 {
1481   // compute the mean of the array "vec" rejecting the outliers
1482   // if rejlist array is provided, fill indices of rejected values
1483   // This method works assuming that the "good" points are nearly normally 
1484   // distrubted around the mean (i.e. there is no significant trend)
1485   //
1486   // create copy of the vector
1487   double *vec1 = new Double_t[np];
1488   memcpy(vec1, vec, np*sizeof(double));
1489   //
1490   // get median of the vector as a robust initial estimate of the mean
1491   cout << "Points " << np << endl;
1492   double md = GetMedian(vec1,np);
1493   //
1494   // compute squared differences to median, 
1495   for (int i=0;i<np;i++) {
1496     vec1[i] = TMath::Abs(vec1[i] - md);
1497   }
1498   //
1499   // compute median squared difference for robust estimate of the sigma
1500   double sigmd = GetMedian(vec1,np);
1501   //
1502   printf("Median Value: %+e | Median abs residual: %e\n",md,sigmd);
1503   //
1504   sigmd /= 0.6745; // go to sigma assuming normal distribution
1505   printf("Estimate of sigma: %+e\n",sigmd);
1506   // if the estimate of the sigma is null, do not look for outliers
1507
1508   cout<<"md="<<md<<endl;
1509
1510   if(!sigmd) return md;
1511
1512  
1513
1514
1515   // compute mean rejecting more than nsigmaCut deviations
1516   double mean = 0;
1517   int npok = 0;
1518   for (int i=0;i<np;i++) {
1519     double dev =  TMath::Abs(vec[i]-md)/sigmd;
1520     if ( dev < nsigmaCut ) {
1521       mean += vec[i];
1522       npok++;
1523     }
1524     else {
1525       printf("Reject#%d (Y=%+e) : deviation: %e\n",i,vec[i],dev);
1526       if (rejList) rejList[nrej] = i; 
1527       nrej++;
1528     }
1529   }
1530   //
1531   delete vec1;
1532   return npok ? mean/npok : 0;
1533   //
1534
1535
1536   cout<<"ending meanMed"<<endl;
1537
1538 }
1539
1540 double GetMedian(double* arr, int n)
1541 {
1542   // get median by Wirths algroithm
1543   int i=0,j=0,l=0,m=0;
1544   double x;
1545   l=0 ; m=n-1;
1546   int k = (n&1) ? n/2 : n/2-1;
1547   while (l<m) {
1548     x=arr[k] ;
1549     i=l ;
1550     j=m ;
1551     do {
1552       //      cout << arr[i] << " " << arr[j] << i<<","<<j << " " << arr[k]  << " " << k << " " << x<< endl;
1553       
1554       while (arr[i]<x) i++ ;
1555       while (x<arr[j]) j-- ;
1556       if (i<=j) { // swap elements        
1557         //      cout << "Swapping" << endl;       
1558         double t = arr[i];
1559         arr[i] = arr[j];
1560         arr[j] = t;
1561         i++ ; j-- ;
1562       }      
1563     } while (i<=j) ;
1564     if (j<k) l=i ;
1565     if (k<i) m=j ;
1566   }
1567   return arr[k] ;
1568 }
1569
1570 // double GetMedian(double* arr, int n)
1571 // {
1572 //   // get median by Wirths algroithm
1573 //   int i,j,l,m;
1574 //   double x;
1575 //   l=0 ; m=n-1;
1576 //   int k = (n&1) ? n/2 : n/2-1;
1577 //   while (l<m) {
1578 //     x=arr[k] ;
1579 //     i=l ;
1580 //     j=m ;
1581 //     do {
1582 //       //      cout << i << " " << j << endl;
1583       
1584 //       //skip runs which were not set (value is -1)
1585 //       while (arr[i]<x) i++ ;
1586 //       while (x<arr[j]) j-- ;
1587 //       if(i<=j) { // swap elements
1588 //      // if((TMath::Abs(arr[i]+1)< 0.0001) && ((TMath::Abs(arr[j]+1)<0.0001))){
1589 //      //   i++ ; j-- ;
1590 //      // } else {
1591 //        double t = arr[i];
1592 //        arr[i] = arr[j];
1593 //        arr[j] = t;
1594 //        i++ ; j-- ;
1595 //        //    }
1596 //       }      
1597 //     } while (i<=j) ;
1598 //     if (j<k) l=i ;
1599 //     if (k<i) m=j ;
1600     
1601 //   }
1602 //   return arr[k] ;
1603 // }
1604
1605 TGraphErrors * GetGraphRej(TGraphErrors * gr, TList * rejRunList, const char * reason, Float_t &mean, Bool_t doDraw) {
1606
1607   //Returns a new graph containing only rejected points
1608
1609   const Double_t nsigmaCut = 3;
1610
1611   int *rejList = new Int_t[gr->GetN()];
1612   int nrej = 0;
1613
1614   Double_t * array = new Double_t[gr->GetN()];
1615   Int_t * correspondenceFullArray = new Int_t[gr->GetN()];
1616   Int_t npoint = 0; 
1617   //exclude from the array all the -1
1618   for(Int_t ipoint = 0; ipoint < gr->GetN(); ipoint++){
1619     if (TMath::Abs(gr->GetY()[ipoint]+1)>0.001) { // skip points for which there is no file (==-1)
1620       array[npoint] = gr->GetY()[ipoint];
1621       correspondenceFullArray[npoint] = ipoint;
1622       npoint++;
1623     } 
1624   }
1625
1626   cout<<"calling meanMed"<<endl;
1627   mean = meanMed(array,npoint,nsigmaCut, nrej, rejList);
1628
1629   cout<<"nrej="<<nrej<<"  mean="<<mean<<endl;
1630   
1631  
1632     
1633   TGraphErrors *grrej = new TGraphErrors(nrej);
1634   for (int i=0;i<nrej;i++) {
1635     grrej->SetPoint(i, gr->GetX()[correspondenceFullArray[rejList[i]]], gr->GetY()[correspondenceFullArray[rejList[i]]]);
1636     grrej->SetPointError(i, gr->GetEX()[correspondenceFullArray[rejList[i]]], gr->GetEY()[correspondenceFullArray[rejList[i]]]);
1637     if(!knownProblems.Contains(Form("%d",(int)gr->GetX()[correspondenceFullArray[rejList[i]]])))
1638       rejRunList->Add(new TObjString(Form("[%d], %s",  (int)gr->GetX()[correspondenceFullArray[rejList[i]]], reason)));
1639   }
1640   grrej->SetMarkerColor(kRed);
1641   grrej->SetMarkerStyle(gr->GetMarkerStyle());
1642   
1643
1644
1645   delete rejList;
1646
1647   if(doDraw) {
1648     Float_t minXDraw = gr->GetXaxis()->GetBinLowEdge(1);
1649     Float_t maxXDraw = gr->GetXaxis()->GetBinLowEdge(gr->GetXaxis()->GetNbins());
1650
1651     grrej->Draw("P");
1652     TLine * l = new TLine (minXDraw, mean,maxXDraw, mean);
1653     l->SetLineStyle(kDashed);
1654     l->Draw();
1655
1656
1657   }
1658
1659   return grrej;
1660
1661
1662
1663
1664
1665 void loadlibs()
1666 {
1667   gSystem->Load("libVMC");
1668   gSystem->Load("libTree");
1669   gSystem->Load("libSTEERBase");
1670   gSystem->Load("libESD");
1671   gSystem->Load("libAOD");
1672   gSystem->Load("libANALYSIS");
1673   gSystem->Load("libANALYSISalice");
1674   gSystem->Load("libCORRFW");
1675   gSystem->Load("libMinuit");
1676   gSystem->Load("libPWG2spectra");
1677   gSystem->Load("libPWG0base"); 
1678 }