]>
Commit | Line | Data |
---|---|---|
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 | } |