TENDER becomes Tender
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / PiKaPr / TestAOD / HighLevelQA / NumberOfEventsMCvsDATA.C
CommitLineData
22d1b439 1/////////////////////////////////////////////////////////
2// NumberOfEventsMCvsDATA.C //
3// //
4// plots the ratio (# of MC events)/(# of DATA events) //
5// against the run number //
6// //
7// Written by John Groh //
8/////////////////////////////////////////////////////////
9
10// for cout
11#include <iostream>
12using namespace std;
13
14// possible cuts to use - 0 1 2 3
15Double_t CentCutMin[4] = {0, 30, 30, 30};
16Double_t CentCutMax[4] = {5, 40, 40, 40};
17Double_t QvecCutMin[4] = {0, 0, 0, 1.5};
18Double_t QvecCutMax[4] = {100, 100, 0.4, 100};
19Double_t EtaMin[4] = {-0.8,-0.8,-0.8,-0.8};
20Double_t EtaMax[4] = {0.8, 0.8, 0.8, 0.8};
21Double_t Nsigmapid = 3.;
22UInt_t trkbit = 1024;
23
24//--------------------
25// Start of function -
26//--------------------
27// icut selects which cut to use (2 and 3 have too low statistics and cause divide-by-zero errors, so don't use)
28// nSigmaCut is the size of the deviation with which to determine outliers
29void NumberOfEventsMCvsDATA(Int_t icut = 1, const Float_t nSigmaCut = 3)
30{
31 // load libraries (I'm not sure I need all of these...)
32 gSystem->Load("libCore.so");
33 gSystem->Load("libPhysics.so");
34 gSystem->Load("libTree");
35 gSystem->Load("libMatrix");
36 gSystem->Load("libSTEERBase");
37 gSystem->Load("libESD");
38 gSystem->Load("libAOD");
39 gSystem->Load("libANALYSIS");
40 gSystem->Load("libOADB");
41 gSystem->Load("libANALYSISalice");
af472fff 42 gSystem->Load("libTender");
22d1b439 43 gSystem->Load("libCORRFW");
44 gSystem->Load("libPWGTools");
45 gSystem->Load("libPWGLFspectra");
46 gSystem->Load("libProof.so");
47 gSystem->AddIncludePath("-I$ALICE_ROOT/include");
48
49 // get number of runs used
50 Int_t nRunsTemp = 0;
51 Int_t dummy;
52 ifstream runList("output/AODQAChecks/RunList.txt");
53 while (!runList.eof())
54 {
55 runList >> dummy;
56 nRunsTemp++;
57 }
58 runList.close();
59
60 // I know this is silly, but you need a const for array sizes.
61 // Also, it seems to want to read one more line than there is in the file
62 const Int_t nRuns = nRunsTemp - 1;
63
64 // fill an array with the run numbers and another two with the file names
65 runList.open("output/AODQAChecks/RunList.txt");
66 Int_t runs[nRuns];
67 TString dataFilesDATA[nRuns];
68 TString dataFilesMC[nRuns];
69
70 for (Int_t irun=0; irun<nRuns; irun++)
71 {
72 runList >> runs[irun];
73 dataFilesDATA[irun] = Form("output/AODQAChecks/DATA/AnalysisResults%i.root",runs[irun]);
74 dataFilesMC[irun] = Form("output/AODQAChecks/MC/AnalysisResults%i.root",runs[irun]);
75 }
76 runList.close();
77
78 // choose which TDirectories in the .root files to use
79 TString snameDATA;
80 TString snameMC;
81 snameDATA = Form("OutputAODSpectraTask_Data_Cent%.0fto%.0f_QVec%.1fto%.1f_Eta%.1fto%.1f_%.1fSigmaPID_TrBit%d",CentCutMin[icut],CentCutMax[icut],QvecCutMin[icut],QvecCutMax[icut],EtaMin[icut],EtaMax[icut],Nsigmapid,trkbit);
82 snameMC = Form("OutputAODSpectraTask_MC_Cent%.0fto%.0f_QVec%.1fto%.1f_Eta%.1fto%.1f_%.1fSigmaPID_TrBit%d",CentCutMin[icut],CentCutMax[icut],QvecCutMin[icut],QvecCutMax[icut],EtaMin[icut],EtaMax[icut],Nsigmapid,trkbit);
83
84 TFile * fileDATA;
85 TDirectoryFile * DirDATA;
86 AliSpectraAODEventCuts * ecutsDATA;
87 TFile * fileMC;
88 TDirectoryFile * DirMC;
89 AliSpectraAODEventCuts * ecutsMC;
90
91 // histogram to be filled w/ the ratio (# of MC events)/(# of DATA events)
92 TH1F * hNEventsRatio = new TH1F("hNEventsRatio","",nRuns,0,nRuns);
93
94 // loop over all runs
95 for (Int_t irun=0; irun<nRuns; irun++)
96 {
97 Printf("\n--- Processing run %i ---",runs[irun]);
98
99 // open the file
100 fileDATA = TFile::Open(dataFilesDATA[irun].Data());
101 if (!fileDATA)
102 {
103 Printf("\n!!! ERROR: Could not open DATA file %s !!!\n",dataFiles[irun].Data());
104 continue;
105 }
106 fileMC = TFile::Open(dataFilesMC[irun].Data());
107 if (!fileMC)
108 {
109 Printf("\n!!! ERROR: Could not open MC file %s !!!\n",dataFiles[irun].Data());
110 continue;
111 }
112
113 // choose the right directory and event cut objects
114 DirDATA = (TDirectoryFile*)fileDATA->Get(Form("%s",snameDATA.Data()));
115 if (!DirDATA)
116 {
117 Printf("\n!!! ERROR: DirDATA is a null pointer. Skipping to next file !!!\n");
118 continue;
119 }
120 DirMC = (TDirectoryFile*)fileMC->Get(Form("%s",snameMC.Data()));
121 if (!DirMC)
122 {
123 Printf("\n!!! ERROR: DirMC is a null pointer. Skipping to next file !!!\n");
124 continue;
125 }
126 ecutsDATA = (AliSpectraAODEventCuts*)DirDATA->Get("Event Cuts");
127 if (!ecutsDATA)
128 {
129 Printf("\n!!! ERROR: ecutsDATA is a null pointer. Skipping to next file !!!\n");
130 continue;
131 }
132 ecutsMC = (AliSpectraAODEventCuts*)DirMC->Get("Event Cuts");
133 if (!ecutsMC)
134 {
135 Printf("\n!!! ERROR: ecutsMC is a null pointer. Skipping to next file !!!\n");
136 continue;
137 }
138
139 // calculate the ratio and fill the histogram
140 Int_t nEventsDATA = ecutsDATA->NumberOfEvents();
141 Int_t nEventsMC = ecutsMC->NumberOfEvents();
142 Float_t ratio = (Float_t)nEventsMC / (Float_t)nEventsDATA;
143 hNEventsRatio->SetBinContent(irun+1,ratio);
144 Printf("# of MC events: %i\n# of DATA events: %i\nRatio: %.2f",nEventsMC,nEventsDATA,ratio);
145
146 } // end loop over runs
147
148 // draw it
149 TCanvas * cNEventsRatio = new TCanvas("cNEventsRatio","cNEventsRatio");
150 for (Int_t irun=0; irun<nRuns; irun++)
151 hNEventsRatio->GetXaxis()->SetBinLabel(irun+1,Form("%i",runs[irun]));
152 hNEventsRatio->GetYaxis()->SetTitle("(# of MC events)/(# of DATA events)");
153 hNEventsRatio->GetYaxis()->SetTitleOffset(1.2);
154 hNEventsRatio->SetMarkerStyle(21);
155 hNEventsRatio->SetMarkerColor(kBlue);
156 hNEventsRatio->SetStats(kFALSE);
157 hNEventsRatio->DrawCopy("P");
158}
159
160
161
162
163