AliCFMuonResTask1.C :
[u/mrichter/AliRoot.git] / CORRFW / test / muons / AliCFMuonResTask1.C
CommitLineData
9a0fd15b 1// DEFINITION OF A FEW CONSTANTS
2
3const Double_t nevtmin= 1;
4const Double_t nevtmax = 15000;
5
6// Muons
7const Double_t ymin = -4.0 ;
8const Double_t ymax = -2.5 ;
9
10const Double_t phimin = -180;
11const Double_t phimax = 180;
12
13// Resonance
14const Int_t PDG = 443;
15
16const Double_t ptmin = 0.0 ;
17const Double_t ptmax = 30 ;
18const Double_t pmin = 0.0 ;
19const Double_t pmax = 700 ;
20const Int_t charge = 0 ;
21const Double_t mmin = 0.1 ;
22const Double_t mmax = 6 ;
23const Double_t mymin = -5 ;
24const Double_t mymax = -1.5 ;
25
26//----------------------------------------------------
27
28Bool_t AliCFMuonResTask1(
29 const Bool_t useGrid = 0,
30 const Bool_t readAOD = 0,
31 const char * kTagXMLFile="wn.xml" // XML file containing tags
32 )
33{
34
35 TBenchmark benchmark;
36 benchmark.Start("AliMuonResTask1");
37
38 AliLog::SetGlobalDebugLevel(0);
39
40 Load() ; // load the required libraries
41
42 TChain * analysisChain ;
43
44///// INPUT
45
46 if (useGrid) { // data located on AliEn
47 TGrid::Connect("alien://") ; // Create an AliRunTagCuts and an AliEventTagCuts Object
48 // and impose some selection criteria
49 AliRunTagCuts *runCuts = new AliRunTagCuts();
50 AliEventTagCuts *eventCuts = new AliEventTagCuts();
51 AliLHCTagCuts *lhcCuts = new AliLHCTagCuts();
52 AliDetectorTagCuts *detCuts = new AliDetectorTagCuts();
53 eventCuts->SetMultiplicityRange(0,2000);
54
55 // Create an AliTagAnalysis Object and chain the tags
56 AliTagAnalysis *tagAna = new AliTagAnalysis();
57 if (readAOD) tagAna->SetType("AOD"); // for aliroot > v4-05
58 else tagAna->SetType("ESD"); // for aliroot > v4-05
59 TAlienCollection *coll = TAlienCollection::Open(kTagXMLFile);
60 TGridResult *tagResult = coll->GetGridResult("",0,0);
61 tagResult->Print();
62 tagAna->ChainGridTags(tagResult);
63
64 // Create a new esd chain and assign the chain that is returned by querying the tags
65 analysisChain = tagAna->QueryTags(runCuts,lhcCuts,detCuts,eventCuts);
66 }
67
68 else {// local data
69 // here put your input data path
70 printf("\n\nRunning on local file, please check the path\n\n");
71
72 if (readAOD) {
73 analysisChain = new TChain("aodTree");
74 analysisChain->Add("AliAOD.root");
75 }
76 else {
77 analysisChain = new TChain("esdTree");
78 analysisChain->Add("/scratch/lopez/PDC08jpsi/run2-300/AliESDs.root");
0dbc8694 79 analysisChain->Add("/scratch/lopez/PDC08jpsi/run3-1000/AliESDs.root");
80 analysisChain->Add("/scratch/lopez/PDC08jpsi/run4-1000/AliESDs.root");
81 analysisChain->Add("/scratch/lopez/PDC08jpsi/run5-1000/AliESDs.root");
82 analysisChain->Add("/scratch/lopez/PDC08jpsi/run7-1000/AliESDs.root");
83 analysisChain->Add("/scratch/lopez/PDC08jpsi/run8-1000/AliESDs.root");
84 analysisChain->Add("/scratch/lopez/PDC08jpsi/run9-1000/AliESDs.root");
85 analysisChain->Add("/scratch/lopez/PDC08jpsi/run10-1000/AliESDs.root");
86 analysisChain->Add("/scratch/lopez/PDC08jpsi/run11-1000/AliESDs.root");
87 analysisChain->Add("/scratch/lopez/PDC08jpsi/run12-1000/AliESDs.root");
88 analysisChain->Add("/scratch/lopez/PDC08jpsi/run13-1000/AliESDs.root");
9a0fd15b 89 }
90 }
91
92///// END INPUT
93
94
95 Info("AliCFMuonResTask1",Form("CHAIN HAS %d ENTRIES",(Int_t)analysisChain->GetEntries()));
96
97 // CONTAINER DEFINITION
98 Info("AliCFMuonResTask1","SETUP CONTAINER");
99
100 // The sensitive variables (9 in this example), their indices
101 UInt_t nevt = 0;
102 UInt_t y1 = 1;
103 UInt_t phi1 = 2;
104 UInt_t y2 = 3;
105 UInt_t phi2 = 4;
106 UInt_t imass = 5;
107 UInt_t y = 6;
108 UInt_t pt = 7;
109 UInt_t p = 8;
110
111
112 // Setting up the container grid
113 UInt_t nstep = 2 ; //number of selection steps : MC and ESD
114 const Int_t nvar = 9 ; //number of variables on the grid
115 const Int_t nbin1 = nevtmax ;
116 const Int_t nbin2 = 100 ;
117 const Int_t nbin3 = 360 ;
118 const Int_t nbin4 = 100 ;
119 const Int_t nbin5 = 360 ;
120 const Int_t nbin6 = 100 ;
121 const Int_t nbin7 = 100 ;
122 const Int_t nbin8 = 100 ;
123 const Int_t nbin9 = 100 ;
124
125 // arrays for the number of bins in each dimension
126 Int_t iBin[nvar];
127 iBin[0]=nbin1;
128 iBin[1]=nbin2;
129 iBin[2]=nbin3;
130 iBin[3]=nbin4;
131 iBin[4]=nbin5;
132 iBin[5]=nbin6;
133 iBin[6]=nbin7;
134 iBin[7]=nbin8;
135 iBin[8]=nbin9;
136
137 // arrays for lower bounds :
138 Double_t *binLim1=new Double_t[nbin1+1];
139 Double_t *binLim2=new Double_t[nbin2+1];
140 Double_t *binLim3=new Double_t[nbin3+1];
141 Double_t *binLim4=new Double_t[nbin4+1];
142 Double_t *binLim5=new Double_t[nbin5+1];
143 Double_t *binLim6=new Double_t[nbin6+1];
144 Double_t *binLim7=new Double_t[nbin7+1];
145 Double_t *binLim8=new Double_t[nbin8+1];
146 Double_t *binLim9=new Double_t[nbin9+1];
147
148 // values for bin lower bounds
149 for(Int_t i=0; i<=nbin1; i++) binLim1[i]=(Double_t)nevtmin + (nevtmax-nevtmin) /nbin1*(Double_t)i ;
150 for(Int_t i=0; i<=nbin2; i++) binLim2[i]=(Double_t)ymin + (ymax-ymin) /nbin2*(Double_t)i ;
151 for(Int_t i=0; i<=nbin3; i++) binLim3[i]=(Double_t)phimin + (phimax-phimin) /nbin3*(Double_t)i ;
152 for(Int_t i=0; i<=nbin4; i++) binLim4[i]=(Double_t)ymin + (ymax-ymin) /nbin4*(Double_t)i ;
153 for(Int_t i=0; i<=nbin5; i++) binLim5[i]=(Double_t)phimin + (phimax-phimin) /nbin5*(Double_t)i ;
154 for(Int_t i=0; i<=nbin6; i++) binLim6[i]=(Double_t)mmin + (mmax-mmin) /nbin6*(Double_t)i ;
155 for(Int_t i=0; i<=nbin7; i++) binLim7[i]=(Double_t)mymin + (mymax-mymin) /nbin7*(Double_t)i ;
156 for(Int_t i=0; i<=nbin8; i++) binLim8[i]=(Double_t)ptmin + (ptmax-ptmin)/nbin8*(Double_t)i ;
157 for(Int_t i=0; i<=nbin9; i++) binLim9[i]=(Double_t)pmin + (pmax-pmin)/nbin9*(Double_t)i ;
158
159 // one container of 2 steps (MC and ESD) with 9 variables
160 AliCFContainer* container = new AliCFContainer("container","container for tracks",nstep,nvar,iBin);
161 // setting the bin limits
162 container -> SetBinLimits(nevt,binLim1);
163 container -> SetBinLimits(y1,binLim2);
164 container -> SetBinLimits(phi1,binLim3);
165 container -> SetBinLimits(y2,binLim4);
166 container -> SetBinLimits(phi2,binLim5);
167 container -> SetBinLimits(imass,binLim6);
168 container -> SetBinLimits(y,binLim7);
169 container -> SetBinLimits(pt,binLim8);
170 container -> SetBinLimits(p,binLim9);
171
172 // Set list
173 TList* qaList = new TList();
174
175 //CREATE THE CUTS
176 // Choice of the Resonance
177 AliCFParticleGenCuts* mcGenCuts = new AliCFParticleGenCuts("mcGenCuts","MC particle generation cuts");
178 mcGenCuts->SetRequirePdgCode(PDG);
179 mcGenCuts->SetQAOn(qaList);
180
181 // Set a pt range of the resonance
182 AliCFTrackKineCuts *mcKineCuts = new AliCFTrackKineCuts("mcKineCuts","MC-level kinematic cuts");
183 mcKineCuts->SetChargeMC(charge);
184 mcKineCuts->SetPtRange(ptmin,ptmax);
185 mcKineCuts->SetQAOn(qaList);
186
187 // Create and fill the list associated
188 TObjArray* mcList = new TObjArray(0) ;
189 mcList->AddLast(mcKineCuts);
190 mcList->AddLast(mcGenCuts);
191
192 // kinematic cuts on muons rapidity
193 AliCFTrackKineCuts *recKineCuts = new AliCFTrackKineCuts("recKineCuts","rec-level kine cuts");
194 recKineCuts->SetRapidityRange(ymin,ymax);
195 recKineCuts->SetQAOn(qaList);
196 TObjArray* recList = new TObjArray(0) ;
197 recList->AddLast(recKineCuts);
198
199 // CREATE THE INTERFACE TO CORRECTION FRAMEWORK USED IN THE TASK
0dbc8694 200
201 printf("CREATE INTERFACE AND CUTS\n");
9a0fd15b 202 AliCFManager* man = new AliCFManager() ;
203 man->SetParticleContainer (container);
9a0fd15b 204
0dbc8694 205 man->SetParticleCutsList(AliCFManager::kPartGenCuts,mcList);
206 man->SetParticleCutsList(AliCFManager::kPartAccCuts,recList);
9a0fd15b 207
208 //CREATE THE TASK
209 printf("CREATE TASK\n");
210 // create the task
211 AliCFMuonResTask1 *task = new AliCFMuonResTask1("AliMuonResTask1");
212 task->SetCFManager(man); //here is set the CF manager
213 task->SetQAList(qaList);
214 if (readAOD) task->SetReadAODData() ;
215
216 //SETUP THE ANALYSIS MANAGER TO READ INPUT CHAIN AND WRITE DESIRED OUTPUTS
217 printf("CREATE ANALYSIS MANAGER\n");
218 // Make the analysis manager
219 AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
220
221 if (useGrid) mgr->SetAnalysisType(AliAnalysisManager::kGridAnalysis);
222 else mgr->SetAnalysisType(AliAnalysisManager::kLocalAnalysis);
223
224
225 AliMCEventHandler* mcHandler = new AliMCEventHandler();
226 mgr->SetMCtruthEventHandler(mcHandler);
227
228 AliInputEventHandler* dataHandler ;
229
230 if (readAOD) dataHandler = new AliAODInputHandler();
231 else dataHandler = new AliESDInputHandler();
232 mgr->SetInputEventHandler(dataHandler);
233
234 // Create and connect containers for input/output
235
236 // input data
237 AliAnalysisDataContainer *cinput0 = mgr->CreateContainer("cchain0",TChain::Class(),AliAnalysisManager::kInputContainer);
238
239 // output data
240 Char_t file[256];
241 sprintf(file,"CFMuonResTask1.root");
242 printf("Analysis output in %s \n",file);
243
244 // output TH1I for event counting
245 AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("chist0", TH1I::Class(),AliAnalysisManager::kOutputContainer,file);
246 // output Correction Framework Container (for acceptance & efficiency calculations)
247 AliAnalysisDataContainer *coutput2 = mgr->CreateContainer("ccontainer0", AliCFContainer::Class(),AliAnalysisManager::kOutputContainer,file);
248
249 cinput0->SetData(analysisChain);
250 mgr->AddTask(task);
0dbc8694 251
252 mgr->ConnectInput(task,0,mgr->GetCommonInputContainer());
253
9a0fd15b 254 mgr->ConnectOutput(task,1,coutput1);
255 mgr->ConnectOutput(task,2,coutput2);
256
257 printf("READY TO RUN\n");
258 //RUN !!!
259 if (mgr->InitAnalysis()) {
260 mgr->PrintStatus();
261 mgr->StartAnalysis("local",analysisChain);
262 }
263
264 benchmark.Stop("AliMuonResTask1");
265 benchmark.Show("AliMuonResTask1");
266
267 return kTRUE ;
268}
269
270void Load() {
271
272 //load the required aliroot libraries
273 gSystem->Load("libANALYSIS") ;
9a0fd15b 274 gSystem->Load("libANALYSISalice") ;
275
0dbc8694 276// gSystem->Load("libCORRFW.so") ;
9a0fd15b 277 gSystem->Load("$ALICE_ROOT/lib/tgt_linux/libCORRFW.so") ;
278
279 //compile online the task class
280 gSystem->SetIncludePath("-I. -I$ALICE_ROOT/include -I$ROOTSYS/include");
281 gROOT->LoadMacro("./AliCFMuonResTask1.cxx+");
282}