AliPhysics  32954cd (32954cd)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
trendingCluster.C
Go to the documentation of this file.
1 
27 
28 /* $Id: $ */
29 //--------------------------------------------------
30 
31 //
32 // Author: Yaxian Mao
33 // Modified M. Germain
34 //
35 
36 #include <TFile.h>
37 #include <TH1F.h>
38 #include <TH2F.h>
39 #include <TH3D.h>
40 #include <Riostream.h>
41 #include <TCanvas.h>
42 #include <TGraphErrors.h>
43 #include <TGrid.h>
44 #include <TFileMerger.h>
45 #include <TMultiGraph.h>
46 #include <TROOT.h>
47 #include <TString.h>
48 #include <TStyle.h>
49 #include <TLegend.h>
50 #include <TGridCollection.h>
51 #include <TROOT.h>
52 #include <TGridResult.h>
53 #include <TClonesArray.h>
54 #include <TObjString.h>
55 #include <stdio.h>
56 
57 #include <string>
58 #include <fstream>
59 #include <iostream>
60 #include <sstream> // Required for stringstreams
61 
62 using namespace std;
63 
64 void trendingCluster(TString fCalorimeter = "EMCAL", TString period = "LHC11h", TString pass = "pass1_HLT", const Int_t n = 10, TString fTrigger = "MB"){
65 
66  FILE * pFile;
67 
68  TString file = "";
69  if (fTrigger=="EMC") file = "/scratch/alicehp2/germain/QA/"+period+"/"+ pass + "/runlist.txt" ;
70  else file = "/scratch/alicehp2/germain/QA/"+period+"/"+ pass + "/runlistMB.txt" ;
71 
72  //cout<<file<<endl;
73  pFile = fopen(file.Data(), "r"); //open the text file where include the run list and correct run index
74 
75  cout<<file<<endl;
76  cout << " fcalo: " << fCalorimeter << "; period: " << period << "; pass: " << pass << " trigger "<<fTrigger<< endl;
77  char outfilename [100] ;
78 
79  Int_t index[500];
80  Int_t p;
81  Int_t q;
82  Int_t ncols;
83  Int_t nlines = 0 ;
84  Int_t RunId[500] ;
85  Double_t x[500] ;
86  Double_t xrun[500] ;
87  //TString label;
88 
89  while (1){
90  ncols = fscanf(pFile,"%d %d ",&p,&q);
91  if (ncols< 0) break;
92  x[nlines]=p;
93  index[nlines]=p;
94  RunId[nlines]=q;
95  xrun[nlines]=1.*q;
96  nlines++;
97  }
98  fclose(pFile);
99  const Int_t nRun = nlines ;
100 
101  Double_t xe[nRun] ;
102  Double_t nEvent[nRun] ;
103  // Double_t TimeMean[nRun] ;
104  // Double_t TimeRMS[nRun] ;
105  Double_t CellMean[nRun] ;
106  Double_t CellRMS[nRun] ;
107  Double_t ClusterMean[nRun] ;
108  Double_t ClusterRMS[nRun] ;
109  Double_t EtotalMean[nRun] ;//total energy deposited per event
110  Double_t EtotalRMS[nRun] ;
111 
112  Double_t CellPerClusterMean[nRun] ;
113  Double_t CellPerClusterRMS[nRun] ;
114  Double_t ECell1Mean[nRun] ;//total energy deposited per event without 1 cell clusters
115  Double_t ECell1RMS[nRun] ;
116 
117  TFile * f ;
118  TDirectoryFile * dir;
119  TList * outputList;
120 
121  TH1F * fhNEvents;
122  TH1F * fhE;
123  TH1F * fhNClusters;
124  TH1F * fhNCells;
125  TH2F * fhIM ;
126  TH3D * fhNCellsPerCluster ;
127  TH2F * fhTimeAmp ;
128 
129  TH1F * NCells[n];
130  TH1F * NClusters[n];
131  TH2F * NCellsPerCluster[n];
132  TH1F * E[n];
133 
134  TGraphErrors * AverNcellsSM[n];
135  TGraphErrors * AverNclustersSM[n];
136  TGraphErrors * AverNcellsPerClusterSM[n];
137  TGraphErrors * AverESM[n];
138  TGraphErrors * AverTimeSM[n];
139  TGraphErrors * AverMggSM[n];
140  TGraphErrors * AverEcell1SM[n];
141 
142  Double_t CellMeanSM[n][nRun] ;
143  Double_t CellRMSSM[n][nRun] ;
144  Double_t ClusterMeanSM[n][nRun] ;
145  Double_t ClusterRMSSM[n][nRun] ;
146  Double_t EtotalMeanSM[n][nRun] ;//total energy deposited per event
147  Double_t EtotalRMSSM[n][nRun] ;
148  Double_t CellPerClusterMeanSM[n][nRun] ;
149  Double_t CellPerClusterRMSSM[n][nRun] ;
150  Double_t ECell1MeanSM[n][nRun] ;//total energy deposited per event without 1 cell clusters
151  Double_t ECell1RMSSM[n][nRun] ;
152 
153 
154  TString namefile = "/scratch/alicehp2/germain/QA/"+period+"/"+pass+"/"+ fCalorimeter + period + pass + fTrigger+"data.txt";
155 
156  fstream QAData(namefile, ios::out); //write the QA check values at the end
157 
158  cout << " namefile " << namefile << endl;
159  cout << " nRun " << nRun << " index(nRun) " << index[nRun-1] << endl;
160 
161  for(Int_t i = 0 ; i < nRun ; i++){
162 
163  xe[i] = 0.5 ;
164 
165  TString name = "/scratch/alicehp2/germain/QA/"+period +"/"+ pass + "/";
166 
167  name += RunId[i] ;
168  name += ".root";
169  f = TFile::Open(name.Data(),"read") ;
170  if (!f) continue;
171 
172  // cout << " i = " << i << " file opend Output" << RunId[i]<< endl;
173 
174  if(fTrigger=="EMC"){ dir = (TDirectoryFile *)f->Get("CaloQA_EMC7");
175  outputList = (TList*)dir->Get("CaloQA_EMC7");
176  }
177  else{
178  dir = (TDirectoryFile *)f->Get("CaloQA_default");
179  outputList = (TList*)dir->Get("CaloQA_default");
180  }
181 
182  //define the averages for checking Histograms
183  //
184 
185  fhNEvents =(TH1F *)outputList->FindObject("hNEvents");
186  nEvent[i]=fhNEvents->GetEntries();
187  cout << " Run " << RunId[i]<< " nevent " << nEvent[i]<< endl;
188  if( nEvent[i] == 0) continue ;
189  if( nEvent[i] < 2) continue ;
190  fhE = (TH1F *)outputList->FindObject(fCalorimeter+"_hE");
191 
192  Double_t energy = 0. ;
193 
194  for(Int_t ibin = fhE->FindBin(0.3) ; ibin <fhE->FindBin(50.) ; ibin++){ //Starting from 0.3eV
195  energy+=fhE->GetBinCenter(ibin)*fhE->GetBinContent(ibin);
196  }
197  EtotalMean[i]=energy/fhE->Integral(fhE->FindBin(0.3), fhE->FindBin(50.)) ;
198  EtotalRMS[i]=fhE->GetMeanError();
199 
200 
201 
202  //for single module check
203  for(Int_t ism = 0 ; ism < n ; ism++){
204  TString nameNCell = Form("%s_hNCells_Mod%d",fCalorimeter.Data(),ism);
205 
206  TString nameNCluster = Form("%s_hNClusters_Mod%d",fCalorimeter.Data(),ism);
207  TString nameNCellPerCluster = Form("%s_hNCellsPerCluster_Mod%d",fCalorimeter.Data(),ism);
208  TString nameE = Form("%s_hE_Mod%d",fCalorimeter.Data(),ism);
209 
210  NCells[ism] = (TH1F*)outputList->FindObject(nameNCell.Data());
211  NClusters[ism] = (TH1F*)outputList->FindObject(nameNCluster.Data());
212  NCellsPerCluster[ism] = (TH2F*)outputList->FindObject(nameNCellPerCluster.Data());
213  E[ism] = (TH1F*)outputList->FindObject(nameE.Data());
214  CellMeanSM[ism][i]=NCells[ism]->GetMean();
215  CellRMSSM[ism][i]=NCells[ism]->GetMeanError();
216  ClusterMeanSM[ism][i]=NClusters[ism]->GetMean();
217  ClusterRMSSM[ism][i]=NClusters[ism]->GetMeanError();
218  CellPerClusterMeanSM[ism][i]=NCellsPerCluster[ism]->GetMean(2);
219  CellPerClusterRMSSM[ism][i]=NCellsPerCluster[ism]->GetMeanError(2);
220  // cout<<"SM = "<<ism<<" Mean : = "<<CellPerClusterMeanSM[ism][i]<<endl ;
221  ECell1MeanSM[ism][i] =NCellsPerCluster[ism]->ProjectionX("",2,300,"")->Integral(5,100)/(nEvent[i]);
222  ECell1RMSSM[ism][i] =NCellsPerCluster[ism]->ProjectionX("",2,300,"")->GetMeanError();
223  Double_t energySM = 0. ;
224  for(Int_t ibin = E[ism]->FindBin(0.3) ; ibin <E[ism]->FindBin(50.) ; ibin++){ //starting from 0.3GeV
225  energySM+=E[ism]->GetBinCenter(ibin)*(E[ism]->GetBinContent(ibin));
226  }
227  EtotalMeanSM[ism][i]=energySM/(E[ism]->Integral(E[ism]->FindBin(0.3),E[ism]->FindBin(50.)));
228 
229  EtotalRMSSM[ism][i]=E[ism]->GetMeanError();
230 
231  if(ism==0) {
232  fhNCells = (TH1F*)NCells[ism]->Clone("NCells");
233  fhNClusters = (TH1F*)NClusters[ism]->Clone("NClusters");
234  }
235  else {
236  fhNCells->Add(NCells[ism],1);
237  fhNClusters->Add(NClusters[ism],1);
238  }
239  } //per SM loop
240  ClusterMean[i]=fhNClusters->GetMean();
241  ClusterRMS[i]=fhNClusters->GetMeanError();
242  CellMean[i]=fhNCells->GetMean();
243  CellRMS[i]=fhNCells->GetMeanError();
244  outputList->Clear() ;
245  dir->Close();
246  f->Close();
247  f=NULL;
248  dir=NULL;
249  outputList=NULL;
250 
251  //if you want to write all the QA check values at the end, otherwise just comment out below
252  // becareful with different detectors as the check output are different since different modules/SM for different detetcor
253 
254  QAData <<i+1<<" "<< RunId[i] <<" "<< nEvent[i]
255  // <<" "<< RunId[i]<<" "<<CellMean[i]
256  // <<" "<< CellMeanSM[0][i] <<" "<< CellMeanSM[1][i]
257  // <<" "<<CellMeanSM[2][i]
258  // <<" "<< RunId[i] <<" "<<ClusterMean[i]
259  // <<" "<< ClusterMeanSM[0][i]<<" "<< ClusterMeanSM[1][i]
260  // <<" "<< ClusterMeanSM[2][i]
261  // <<" "<< RunId[i]<<" "<< CellPerClusterMean[i]
262  // <<" "<< CellPerClusterMeanSM[0][i]<<" "<< CellPerClusterMeanSM[1][i]
263  // <<" "<< CellPerClusterMeanSM[2][i]
264  // <<" "<< RunId[i]<<" "<< EtotalMean[i]
265  // <<" "<< EtotalMeanSM[0][i] <<" "<< EtotalMeanSM[1][i]
266  // <<" "<< EtotalMeanSM[2][i]
267 
268  <<"\n" ;
269 
270  } // end loop on nrun
271 
272 
273  QAData.close();
274 
275  TString base = "/scratch/alicehp2/germain/QA/";
276  base += period ;
277  base += "/";
278  base += pass ;
279  base += "/";
280 
281  base += fTrigger;
282  TString ClusterAverages ; ClusterAverages = base + "ClusterAverages.gif";
283  TString Entries ; Entries = base + "Nentries.gif";
284  TString ClusterAveragesEnergy ; ClusterAveragesEnergy = base + "ClusterAveragesEnergy.gif";
285  TString ClusterAveragesEntries ; ClusterAveragesEntries = base + "ClusterAveragesEntries.gif";
286  TString ClusterAveragesCells ; ClusterAveragesCells = base + "ClusterAveragescells.gif";
287 
288 
289  cout << "c11 nEvents" << endl;
290  //just for the canvas defination
291 
292  cout << " index(0)" << index[nRun-1] << " index(20) " << index[20] <<endl;
293  TH1F * dummy = new TH1F("dummy", "dummy", nRun, 0., nRun+0.5);
294  // TH1F * dummy = new TH1F("dummy", "dummy", index[nRun-1], 0., index[nRun-1]+0.5);
295  dummy->SetTitle("") ;
296  dummy->SetStats(kFALSE) ;
297  dummy->SetAxisRange(0., nRun, "X") ;
298 
299  for(Int_t i = 0 ; i < nRun ; i++){
300  TString label = " ";
301  label+=RunId[i];
302  cout <<" run "<< RunId[i] << " label " <<label << endl;
303 
304  dummy->GetXaxis()->SetBinLabel(i+1,label.Data());
305  dummy->GetXaxis()->LabelsOption("v");
306  }
307 
308 
309 
310  //number of events passes physics selection for each run
311  TCanvas * c11 = new TCanvas("nEvents", "nEvents", 1000, 500);
312  c11->SetFillColor(0);
313  c11->SetBorderSize(0);
314  c11->SetFrameBorderMode(0);
315  gStyle->SetOptStat(0);
316  c11->SetLogy();
317  c11->SetGrid();
318  // dummy->GetXaxis()->SetTitle("RUN");
319  dummy->GetXaxis()->SetTitleOffset(0.05);
320  dummy->GetYaxis()->SetTitle("N_{events}");
321  dummy->SetMinimum(1.) ; //should addjust based on the statistics
322  dummy->SetMaximum(1.e6) ; //should addjust based on the statistics
323  dummy->Draw();
324  TGraph * nEvents = new TGraph(nRun, x, nEvent);
325  nEvents->SetMarkerStyle(20);
326  nEvents->SetMarkerColor(1);
327  nEvents->SetLineColor(2);
328  nEvents->Draw("same PL") ;
329  c11->Update();
330 
331 
332  if (fTrigger=="MB")sprintf(outfilename,"nEventMB.gif");
333  if (fTrigger=="EMC")sprintf(outfilename,"nEventEMC.gif");
334 
335 
336  c11->SaveAs(Entries);
337 
338 
339  cout << "c1 Aver NCell" << endl;
340 
341  TCanvas * c1 = new TCanvas("AverNCell", "AverNCell", 600, 600);
342  // c1->SetLogy();
343  c1->SetFillColor(0);
344  c1->SetBorderSize(0);
345  c1->SetFrameBorderMode(0);
346  gStyle->SetOptStat(0);
347  TH1F * h1 = (TH1F*)dummy->Clone("");
348  h1->GetXaxis()->SetTitle("RUN Index");
349  h1->GetYaxis()->SetTitle("<N_{cells}>");
350  h1->SetMinimum(0.) ;
351  h1->SetMaximum(10.) ;
352  if(fCalorimeter=="EMCAL") h1->SetMaximum(5.) ;
353  h1->Draw();
354 
355  TGraphErrors * AverNcells = new TGraphErrors(nRun, x, CellMean, xe, CellRMS);
356 
357  AverNcells->SetMarkerColor(1);
358  AverNcells->SetMarkerStyle(20);
359  AverNcells->Draw("same P") ;
360  for(Int_t ism = 0 ; ism < n ; ism++){
361  AverNcellsSM[ism] = new TGraphErrors(nRun, x, CellMeanSM[ism], xe, CellRMSSM[ism]);
362  AverNcellsSM[ism]->SetMarkerColor(ism+2);
363  AverNcellsSM[ism]->SetMarkerStyle(21+ism);
364  AverNcellsSM[ism]->Draw("same P");
365  }
366 
367 
368 
369  TLegend * l1 = new TLegend(0.4, 0.6, 0.75, 0.85);
370  l1->SetFillColor(0);
371  l1->SetBorderSize(0);
372  l1->SetTextSize(0.02);
373  l1->AddEntry(AverNcells, "<# of cells>", "") ;
374  l1->AddEntry(AverNcells, Form("det = %s",fCalorimeter.Data()), "") ;
375  l1->AddEntry(AverNcells,"average", "p");
376  for(Int_t ism = 0 ; ism < n ; ism++){
377  TString projname = Form("SM_ %d",ism);
378  l1->AddEntry(AverNcellsSM[ism],projname.Data(), "p");
379  }
380  l1->Draw("same");
381  c1->Update();
382 
383 
384  TCanvas * c200 = new TCanvas("ClusterAverages", "ClusterAverages", 1000, 500);
385  c200->SetFillColor(0);
386  c200->SetBorderSize(0);
387  c200->SetFrameBorderMode(0);
388  c200->SetGrid();
389 
390  gPad->SetLeftMargin(0.08);
391  gPad->SetRightMargin(0.02);
392  gPad->SetGrid();
393 
394  TH1F * h2 = (TH1F*)dummy->Clone("");
395 
396  h2->GetYaxis()->SetTitle("<N_{clusters}>/event");
397 
398  if (fTrigger=="EMC") h2->SetMaximum(70) ;
399  else h2->SetMaximum(50) ;
400  h2->SetMinimum(10) ; // for Pb Pb
401  h2->Draw();
402 
403  TGraphErrors * AverNclusters = new TGraphErrors(nRun, x, ClusterMean, xe, ClusterRMS);
404  AverNclusters->SetMarkerStyle(20);
405  AverNclusters->SetMarkerColor(1);
406  AverNclusters->Draw("same P") ;
407 
408  for(Int_t ism = 0 ; ism < n ; ism++){
409  AverNclustersSM[ism] = new TGraphErrors(nRun, x, ClusterMeanSM[ism], xe, ClusterRMSSM[ism]);
410  AverNclustersSM[ism]->SetMarkerColor(ism+2);
411  AverNclustersSM[ism]->SetMarkerStyle(21+ism);
412 
413  AverNclustersSM[ism]->Draw("same P");
414  }
415 
416 
417  TLegend * l200 = new TLegend(0.4, 0.6, 0.75, 0.85);
418  l200->SetFillColor(0);
419  l200->SetBorderSize(0);
420  l200->SetTextSize(0.02);
421  l200->AddEntry(AverNclusters, "<# of clusters>", "") ;
422  l200->AddEntry(AverNclusters, Form("det = %s",fCalorimeter.Data()), "") ;
423  l200->AddEntry(AverNclusters,"average", "p");
424 
425  for(Int_t ism = 0 ; ism < n ; ism++){
426  TString projname = Form("SM_ %d",ism);
427  l200->AddEntry(AverNclustersSM[ism],projname.Data(), "p");
428  }
429 
430  l200->Draw("same");
431  c200->Update();
432  c200->SaveAs(ClusterAveragesEntries);
433 
434  TCanvas * c201 = new TCanvas("ClusterAveragesEnergy", "ClusterAveragesEnergy", 1000, 500);
435  c201->SetFillColor(0);
436  c201->SetBorderSize(0);
437  c201->SetFrameBorderMode(0);
438  c201->SetGrid();
439 
440 
441  gPad->SetLeftMargin(0.08);
442  gPad->SetRightMargin(0.02);
443  gPad->SetGrid();
444 
445 
446  TH1F * h3 = (TH1F*)dummy->Clone("");
447 
448  h3->GetYaxis()->SetTitle("<E> (GeV)");
449  h3->SetMinimum(0.2) ;
450  h3->SetMaximum(1.) ;
451  h3->Draw();
452 
453  TGraphErrors * AverE = new TGraphErrors(nRun, x, EtotalMean, xe, EtotalRMS);
454  AverE->SetMarkerStyle(20);
455  AverE->SetMarkerColor(1);
456  AverE->Draw("same P");
457 
458  for(Int_t ism = 0 ; ism < n ; ism++){
459  AverESM[ism] = new TGraphErrors(nRun, x, EtotalMeanSM[ism], xe, EtotalRMSSM[ism]);
460  AverESM[ism]->SetMarkerColor(ism+2);
461  AverESM[ism]->SetMarkerStyle(21+ism);
462  AverESM[ism]->Draw("same P");
463 
464  }
465 
466  TLegend * l3 = new TLegend(0.4, 0.6, 0.75, 0.85);
467  l3->SetFillColor(0);
468  l3->SetBorderSize(0);
469  l3->SetTextSize(0.02);
470  l3->AddEntry(AverE, "<E>", "") ;
471  l3->AddEntry(AverE, Form("det = %s",fCalorimeter.Data()), "") ;
472  l3->AddEntry(AverE,"average", "p");
473 
474  for(Int_t ism = 0 ; ism < n ; ism++){
475  TString projname = Form("SM_ %d",ism);
476  l3->AddEntry(AverESM[ism],projname.Data(), "p");
477  }
478 
479  l3->Draw("same");
480 
481  c201->SaveAs(ClusterAveragesEnergy);
482 
483 
484  TCanvas * c202 = new TCanvas("ClusterAveragesCells", "ClusterAveragesCells", 1000, 500);
485  c202->SetFillColor(0);
486  c202->SetBorderSize(0);
487  c202->SetFrameBorderMode(0);
488  c202->SetGrid();
489 
490  gPad->SetLeftMargin(0.08);
491  gPad->SetRightMargin(0.02);
492 
493  gPad->SetGrid();
494 
495 
496  TH1F * h4 = (TH1F*)dummy->Clone("");
497 
498  h4->GetYaxis()->SetTitle("<N_{CellsPerCluster}>");
499  h4->SetMinimum(1.5) ;
500  h4->SetMaximum(5.5) ;
501  h4->Draw();
502 
503 
504  TGraphErrors * AverCellPerCluster = new TGraphErrors(nRun, x, CellPerClusterMean, xe, CellPerClusterRMS);
505  AverCellPerCluster->SetMarkerStyle(20);
506  AverCellPerCluster->SetMarkerColor(1);
507 
508  for(Int_t ism = 0 ; ism < n ; ism++){
509  AverNcellsPerClusterSM[ism] = new TGraphErrors(nRun, x, CellPerClusterMeanSM[ism], xe, CellPerClusterRMSSM[ism]);
510  AverNcellsPerClusterSM[ism]->SetMarkerColor(ism+2);
511  AverNcellsPerClusterSM[ism]->SetMarkerStyle(21+ism);
512  AverNcellsPerClusterSM[ism]->Draw("same P");
513 
514  }
515 
516  TLegend * l4 = new TLegend(0.4, 0.6, 0.75, 0.85);
517  l4->SetFillColor(0);
518  l4->SetBorderSize(0);
519  l4->SetTextSize(0.02);
520  l4->AddEntry(AverCellPerCluster, "<# of cells per cluster>", "") ;
521  l4->AddEntry(AverCellPerCluster, Form("det = %s",fCalorimeter.Data()), "") ;
522  l4->AddEntry(AverCellPerCluster,"average", "p");
523 
524  for(Int_t ism = 0 ; ism < n ; ism++){
525  TString projname = Form("SM_ %d",ism);
526  l4->AddEntry(AverNcellsPerClusterSM[ism],projname.Data(), "p");
527  }
528 
529  l4->Draw("same");
530 
531  c202->SaveAs(ClusterAveragesCells);
532 
533 }
Double_t
double Double_t
Definition: External.C:58
TH2F
Definition: External.C:236
energy
energy
Definition: HFPtSpectrum.C:44
Int_t
int Int_t
Definition: External.C:63
TH3D
Definition: External.C:252
nEvents
Double_t nEvents
plot quality messages
Definition: DrawAnaCaloTrackQA.C:95
file
TFile * file
TList with histograms for a given trigger.
Definition: DrawAnaCaloTrackQA.C:89
dir
TDirectoryFile * dir
Definition: DrawAnaCaloTrackQA.C:87