Implemented ZDC time cut in phsyics selection and in trigger analysis for MC. Trackle...
[u/mrichter/AliRoot.git] / CORRFW / test / AliCFTaskForUnfolding.C
1 //DEFINITION OF A FEW CONSTANTS
2 const Double_t ptmin  =   0.0 ;
3 const Double_t ptmax  =   1.0 ;
4 const Double_t etamin =  -1.5 ;
5 const Double_t etamax =   1.5 ;
6 const Int_t    PDG    =    211; 
7 const Int_t    minclustersTPC = 40 ;
8 //----------------------------------------------------
9
10 /*
11
12 Macro to prepare the necessary objects to perform spectrum unfolding 
13 (to be used by the class AliCFUnfolding)
14
15 Note that the efficiency calculation (and therfore the container filling) 
16 has to be done using MC values (and not reconstructed values)
17
18 This macro creates the following objects :
19 - (AliCFContainer) container    : used to calculate the efficiency (see AliCFEffGrid)
20 - (THnSparseD)     correlation  : this is the response matrix
21
22 Once you have run this macro, you may launch unfolding using 
23 the example macro CORRFW/test/testUnfolding.C
24
25 */
26
27 Bool_t AliCFTaskForUnfolding()
28 {
29   
30   TBenchmark benchmark;
31   benchmark.Start("AliSingleTrackTask");
32
33   AliLog::SetGlobalDebugLevel(0);
34
35   Load() ; //load the required libraries
36
37   TChain * analysisChain ;
38   printf("\n\nRunning on local file, please check the path\n\n");
39
40   analysisChain = new TChain("esdTree");
41   analysisChain->Add("your_data_path/001/AliESDs.root");
42   analysisChain->Add("your_data_path/002/AliESDs.root");
43
44   
45   Info("AliCFTaskForUnfolding",Form("CHAIN HAS %d ENTRIES",(Int_t)analysisChain->GetEntries()));
46
47   //CONTAINER DEFINITION
48   Info("AliCFTaskForUnfolding","SETUP CONTAINER");
49   //the sensitive variables (2 in this example), their indices
50   UInt_t ipt = 0;
51   UInt_t ieta  = 1;
52   //Setting up the container grid... 
53   UInt_t nstep = 3 ; //number of selection steps MC 
54   const Int_t nvar   = 2 ; //number of variables on the grid:pt,eta
55   const Int_t nbin[nvar] = {20,20} ;
56
57   //arrays for the number of bins in each dimension
58   Int_t iBin[nvar];
59   iBin[0]=nbin[0];
60   iBin[1]=nbin[1];
61
62   //arrays for lower bounds :
63   Double_t *binLim0=new Double_t[nbin[0]+1];
64   Double_t *binLim1=new Double_t[nbin[1]+1];
65
66   //values for bin lower bounds
67   for(Int_t i=0; i<=nbin[0]; i++) binLim0[i]=(Double_t)ptmin  + ( ptmax- ptmin)/nbin[0]*(Double_t)i ; 
68   for(Int_t i=0; i<=nbin[1]; i++) binLim1[i]=(Double_t)etamin + (etamax-etamin)/nbin[1]*(Double_t)i ; 
69
70   //one "container" for MC
71   AliCFContainer* container = new AliCFContainer("container","container for tracks",nstep,nvar,iBin);
72   //setting the bin limits
73   container -> SetBinLimits(ipt,binLim0);
74   container -> SetBinLimits(ieta ,binLim1);
75   container -> SetVarTitle(ipt,"pT");
76   container -> SetVarTitle(ieta,"#eta");
77
78   // Gen-Level kinematic cuts
79   AliCFTrackKineCuts *mcKineCuts = new AliCFTrackKineCuts("mcKineCuts","MC-level kinematic cuts");
80   mcKineCuts->SetPtRange (ptmin ,ptmax );
81   mcKineCuts->SetEtaRange(etamin,etamax);
82
83   //Particle-Level cuts:  
84   AliCFParticleGenCuts* mcGenCuts = new AliCFParticleGenCuts("mcGenCuts","MC particle generation cuts");
85   mcGenCuts->SetRequireIsPrimary();
86   mcGenCuts->SetRequirePdgCode(PDG);
87
88   // Rec-Level kinematic cuts
89   AliCFTrackKineCuts *recKineCuts = new AliCFTrackKineCuts("recKineCuts","rec-level kine cuts");
90   recKineCuts->SetPtRange( ptmin, ptmax);
91   recKineCuts->SetPtRange(etamin,etamax);
92
93   AliCFTrackQualityCuts *recQualityCuts = new AliCFTrackQualityCuts("recQualityCuts","rec-level quality cuts");
94   recQualityCuts->SetMinNClusterTPC(minclustersTPC);
95
96   AliCFTrackIsPrimaryCuts *recIsPrimaryCuts = new AliCFTrackIsPrimaryCuts("recIsPrimaryCuts","rec-level isPrimary cuts");
97   recIsPrimaryCuts->SetMaxNSigmaToVertex(3);
98
99   printf("CREATE MC KINE CUTS\n");
100   TObjArray* mcList = new TObjArray(0) ;
101   mcList->AddLast(mcKineCuts);
102   mcList->AddLast(mcGenCuts);
103
104   printf("CREATE RECONSTRUCTION CUTS\n");
105   TObjArray* recList = new TObjArray(0) ;
106   recList->AddLast(recKineCuts);
107   recList->AddLast(recQualityCuts);
108   recList->AddLast(recIsPrimaryCuts);
109
110   TObjArray* emptyList = new TObjArray(0);
111
112   //CREATE THE INTERFACE TO CORRECTION FRAMEWORK USED IN THE TASK
113   printf("CREATE INTERFACE AND CUTS\n");
114   AliCFManager* man = new AliCFManager() ;
115   man->SetNStepEvent(0);
116   man->SetParticleContainer(container);
117   man->SetParticleCutsList(0,mcList);
118   man->SetParticleCutsList(1,recList);
119   man->SetParticleCutsList(2,emptyList); // this is special step for monte carlo spectrum
120
121   //CREATE THE TASK
122   printf("CREATE TASK\n");
123   // create the task
124   AliCFTaskForUnfolding *task = new AliCFTaskForUnfolding("AliCFTaskForUnfolding");
125   task->SetCFManager(man); //here is set the CF manager
126
127   //create correlation matrix for unfolding
128   Int_t thnDim[2*nvar];
129   for (int k=0; k<nvar; k++) {
130     //first half  : reconstructed 
131     //second half : MC
132     thnDim[k]      = nbin[k];
133     thnDim[k+nvar] = nbin[k];
134   }
135   THnSparseD* correlation = new THnSparseD("correlation","THnSparse with correlations",2*nvar,thnDim);
136   Double_t** binEdges = new Double_t[nvar];
137   for (int k=0; k<nvar; k++) {
138     binEdges[k]=new Double_t[nbin[k]+1];
139     container->GetBinLimits(k,binEdges[k]);
140     correlation->SetBinEdges(k,binEdges[k]);
141     correlation->SetBinEdges(k+nvar,binEdges[k]);
142   }
143   correlation->Sumw2();
144   task->SetCorrelationMatrix(correlation);
145   //---
146
147   //SETUP THE ANALYSIS MANAGER TO READ INPUT CHAIN AND WRITE DESIRED OUTPUTS
148   printf("CREATE ANALYSIS MANAGER\n");
149   // Make the analysis manager
150   AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
151   mgr->SetAnalysisType(AliAnalysisManager::kLocalAnalysis);
152
153   AliMCEventHandler*  mcHandler = new AliMCEventHandler();
154   mgr->SetMCtruthEventHandler(mcHandler);
155  
156   AliInputEventHandler* dataHandler = new AliESDInputHandler();
157   mgr->SetInputEventHandler(dataHandler);
158
159   // Create and connect containers for input/output
160
161   //------ input data ------
162   AliAnalysisDataContainer *cinput0  = mgr->CreateContainer("cchain0",TChain::Class(),AliAnalysisManager::kInputContainer);
163
164   // ----- output data -----
165   
166   //slot 0 : default output tree (by default handled by AliAnalysisTaskSE)
167   AliAnalysisDataContainer *coutput0 = mgr->CreateContainer("ctree0", TTree::Class(),AliAnalysisManager::kOutputContainer,"output.root");
168
169   //now comes user's output objects :
170   
171   // output TH1I for event counting
172   AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("chist0", TH1I::Class(),AliAnalysisManager::kOutputContainer,"output.root");
173   // output Correction Framework Container (for acceptance & efficiency calculations)
174   AliAnalysisDataContainer *coutput2 = mgr->CreateContainer("ccontainer0", AliCFContainer::Class(),AliAnalysisManager::kOutputContainer,"output.root");
175   AliAnalysisDataContainer *coutput3 = mgr->CreateContainer("corr0", THnSparseD::Class(),AliAnalysisManager::kOutputContainer,"output.root");
176
177   cinput0->SetData(analysisChain);
178
179   mgr->AddTask(task);
180   mgr->ConnectInput(task,0,mgr->GetCommonInputContainer());
181   mgr->ConnectOutput(task,0,coutput0);
182   mgr->ConnectOutput(task,1,coutput1);
183   mgr->ConnectOutput(task,2,coutput2);
184   mgr->ConnectOutput(task,3,coutput3);
185  
186   printf("READY TO RUN\n");
187   //RUN !!!
188   if (mgr->InitAnalysis()) {
189     mgr->PrintStatus();
190     mgr->StartAnalysis("local",analysisChain);
191   }
192
193   benchmark.Stop("AliSingleTrackTask");
194   benchmark.Show("AliSingleTrackTask");
195
196   return kTRUE ;
197 }
198
199 void Load() {
200
201   //load the required aliroot libraries
202   gSystem->Load("libANALYSIS") ;
203   gSystem->Load("libANALYSISalice") ;
204   gSystem->Load("libCORRFW.so") ;
205
206   //compile online the task class
207   gSystem->SetIncludePath("-I. -I$ALICE_ROOT/include -I$ROOTSYS/include");
208   gROOT->LoadMacro("./AliCFTaskForUnfolding.cxx+");
209 }