]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskEfficiencyBF.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'master_patch'
[u/mrichter/AliRoot.git] / PWGCF / EBYE / BalanceFunctions / AliAnalysisTaskEfficiencyBF.cxx
CommitLineData
84164254 1#include "TChain.h"
2#include "TTree.h"
3#include "TH1F.h"
4#include "TH2F.h"
40b41983 5#include "TH3D.h"
84164254 6#include "TCanvas.h"
7#include "TParticle.h"
8#include "TObjArray.h"
9
10#include "AliAnalysisTask.h"
11#include "AliAnalysisManager.h"
12
13#include "AliStack.h"
14#include "AliMCEvent.h"
15#include "AliESDEvent.h"
16#include "AliESDInputHandler.h"
17#include "AliESDtrackCuts.h"
18#include "AliCentrality.h"
6d6eb2df 19#include "AliLog.h"
84164254 20
21#include "AliAnalysisTaskEfficiencyBF.h"
22
23// ---------------------------------------------------------------------
24//
25// Task for calculating the efficiency of the Balance Function
26// for single particles and pairs
27//
28// Authors: Panos Christakoglou, Michael Weber
29//
30// ---------------------------------------------------------------------
31
32ClassImp(AliAnalysisTaskEfficiencyBF)
33
34//________________________________________________________________________
35AliAnalysisTaskEfficiencyBF::AliAnalysisTaskEfficiencyBF(const char *name)
36 : AliAnalysisTaskSE(name), fESD(0), fQAList(0), fOutputList(0),
37 fHistEventStats(0), fHistCentrality(0), fHistNMult(0),
40b41983 38 fHistGeneratedEtaPtPhiPlus(0), fHistFindableEtaPtPhiPlus(0),
39 fHistReconstructedEtaPtPhiPlus(0), fHistSurvivedEtaPtPhiPlus(0),
40 fHistGeneratedEtaPtPhiMinus(0), fHistFindableEtaPtPhiMinus(0),
41 fHistReconstructedEtaPtPhiMinus(0), fHistSurvivedEtaPtPhiMinus(0),
84164254 42 fHistGeneratedEtaPtPlusControl(0), fHistFindableEtaPtPlusControl(0),
43 fHistReconstructedEtaPtPlusControl(0), fHistSurvivedEtaPtPlusControl(0),
44 fHistGeneratedEtaPtMinusControl(0), fHistFindableEtaPtMinusControl(0),
45 fHistReconstructedEtaPtMinusControl(0), fHistSurvivedEtaPtMinusControl(0),
46 fHistGeneratedEtaPtPlusPlus(0), fHistFindableEtaPtPlusPlus(0),
47 fHistReconstructedEtaPtPlusPlus(0), fHistSurvivedEtaPtPlusPlus(0),
48 fHistGeneratedEtaPtMinusMinus(0), fHistFindableEtaPtMinusMinus(0),
49 fHistReconstructedEtaPtMinusMinus(0), fHistSurvivedEtaPtMinusMinus(0),
50 fHistGeneratedEtaPtPlusMinus(0), fHistFindableEtaPtPlusMinus(0),
51 fHistReconstructedEtaPtPlusMinus(0), fHistSurvivedEtaPtPlusMinus(0),
40b41983 52 fHistGeneratedPhiEtaPlusPlus(0), fHistFindablePhiEtaPlusPlus(0),
53 fHistReconstructedPhiEtaPlusPlus(0), fHistSurvivedPhiEtaPlusPlus(0),
54 fHistGeneratedPhiEtaMinusMinus(0), fHistFindablePhiEtaMinusMinus(0),
55 fHistReconstructedPhiEtaMinusMinus(0), fHistSurvivedPhiEtaMinusMinus(0),
56 fHistGeneratedPhiEtaPlusMinus(0), fHistFindablePhiEtaPlusMinus(0),
57 fHistReconstructedPhiEtaPlusMinus(0), fHistSurvivedPhiEtaPlusMinus(0),
84164254 58 fESDtrackCuts(0), fAnalysisMode(0),
59 fCentralityEstimator("V0M"), fCentralityPercentileMin(0.0), fCentralityPercentileMax(5.0),
60 fVxMax(3.0), fVyMax(3.0), fVzMax(10.),
61 fMinNumberOfTPCClusters(80), fMaxChi2PerTPCCluster(4.0), fMaxDCAxy(3.0), fMaxDCAz(3.0),
d9906391 62 fMinPt(0.3), fMaxPt(1.5), fMinEta(-0.8), fMaxEta(0.8),fEtaRangeMin(0.0), fEtaRangeMax(1.6), fPtRangeMin(0.3), fPtRangeMax(1.5), fPhiRangeMin(0.0),fPhiRangeMax(360.),fdPhiRangeMax(180.), fEtaBin(100),fdEtaBin(64),fPtBin(49),fPhiBin(100),fdPhiBin(90) {
84164254 63 // Define input and output slots here
64 // Input slot #0 works with a TChain
65 DefineInput(0, TChain::Class());
66 // Output slot #0 id reserved by the base class for AOD
67 // Output slot #1 writes into a TH1 container
68 DefineOutput(1, TList::Class());
69 DefineOutput(2, TList::Class());
70}
71
72//________________________________________________________________________
73void AliAnalysisTaskEfficiencyBF::UserCreateOutputObjects() {
74 // Create histograms
75 // Called once
76
211b716d 77 // global switch disabling the reference
78 // (to avoid "Replacing existing TH1" if several wagons are created in train)
79 Bool_t oldStatus = TH1::AddDirectoryStatus();
80 TH1::AddDirectory(kFALSE);
81
84164254 82 fQAList = new TList();
83 fQAList->SetName("QAList");
84 fQAList->SetOwner();
85
86 fOutputList = new TList();
87 fOutputList->SetName("OutputList");
88 fOutputList->SetOwner();
89
90 //Event stats.
91 TString gCutName[4] = {"Total","Offline trigger",
92 "Vertex","Analyzed"};
93 fHistEventStats = new TH1F("fHistEventStats",
94 "Event statistics;;N_{events}",
95 4,0.5,4.5);
96 for(Int_t i = 1; i <= 4; i++)
97 fHistEventStats->GetXaxis()->SetBinLabel(i,gCutName[i-1].Data());
98 fQAList->Add(fHistEventStats);
99
100 //ESD analysis
101 fHistCentrality = new TH1F("fHistCentrality",";Centrality bin;Events",
102 20,0.5,20.5);
103 fQAList->Add(fHistCentrality);
104
105 //multiplicity (good MC tracks)
106 TString histName;
107 histName = "fHistNMult";
108 fHistNMult = new TH1F(histName.Data(),
109 ";N_{mult.}",
110 200,0,20000);
111 fQAList->Add(fHistNMult);
112
113 //eta vs pt for MC positives
40b41983 114 fHistGeneratedEtaPtPhiPlus = new TH3D("fHistGeneratedEtaPtPhiPlus",
115 "Generated positive primaries;#eta;p_{T} (GeV/c);#phi",
622c9e7f 116 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax,fPhiBin,fPhiRangeMin,fPhiRangeMax);
40b41983 117 fOutputList->Add(fHistGeneratedEtaPtPhiPlus);
118 fHistFindableEtaPtPhiPlus = new TH3D("fHistFindableEtaPtPhiPlus",
119 "Findable positive primaries;#eta;p_{T} (GeV/c);#phi",
622c9e7f 120 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax,fPhiBin,fPhiRangeMin,fPhiRangeMax);
40b41983 121 fOutputList->Add(fHistFindableEtaPtPhiPlus);
122 fHistReconstructedEtaPtPhiPlus = new TH3D("fHistReconstructedEtaPtPhiPlus",
123 "Reconstructed positive primaries;#eta;p_{T} (GeV/c);#phi",
622c9e7f 124 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax,fPhiBin,fPhiRangeMin,fPhiRangeMax);
40b41983 125 fOutputList->Add(fHistReconstructedEtaPtPhiPlus);
126 fHistSurvivedEtaPtPhiPlus = new TH3D("fHistSurvivedEtaPtPhiPlus",
127 "Survived positive primaries;#eta;p_{T} (GeV/c);#phi",
622c9e7f 128 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax,fPhiBin,fPhiRangeMin,fPhiRangeMax);
40b41983 129 fOutputList->Add(fHistSurvivedEtaPtPhiPlus);
84164254 130
131 //eta vs pt for MC negatives
40b41983 132 fHistGeneratedEtaPtPhiMinus = new TH3D("fHistGeneratedEtaPtPhiMinus",
133 "Generated positive primaries;#eta;p_{T} (GeV/c);#phi",
622c9e7f 134 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax,fPhiBin,fPhiRangeMin,fPhiRangeMax);
40b41983 135 fOutputList->Add(fHistGeneratedEtaPtPhiMinus);
136 fHistFindableEtaPtPhiMinus = new TH3D("fHistFindableEtaPtPhiMinus",
137 "Findable positive primaries;#eta;p_{T} (GeV/c);#phi",
622c9e7f 138 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax,fPhiBin,fPhiRangeMin,fPhiRangeMax);
40b41983 139 fOutputList->Add(fHistFindableEtaPtPhiMinus);
140 fHistReconstructedEtaPtPhiMinus = new TH3D("fHistReconstructedEtaPtPhiMinus",
141 "Reconstructed positive primaries;#eta;p_{T} (GeV/c);#phi",
622c9e7f 142 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax,fPhiBin,fPhiRangeMin,fPhiRangeMax);
40b41983 143 fOutputList->Add(fHistReconstructedEtaPtPhiMinus);
144 fHistSurvivedEtaPtPhiMinus = new TH3D("fHistSurvivedEtaPtPhiMinus",
145 "Survived positive primaries;#eta;p_{T} (GeV/c);#phi",
622c9e7f 146 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax,fPhiBin,fPhiRangeMin,fPhiRangeMax);
40b41983 147 fOutputList->Add(fHistSurvivedEtaPtPhiMinus);
84164254 148
149 //eta vs pt for MC positives (control)
150 fHistGeneratedEtaPtPlusControl = new TH2F("fHistGeneratedEtaPtPlusControl",
151 "Generated positive primaries;#eta;p_{T} (GeV/c)",
622c9e7f 152 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 153 fOutputList->Add(fHistGeneratedEtaPtPlusControl);
154 fHistFindableEtaPtPlusControl = new TH2F("fHistFindableEtaPtPlusControl",
155 "Findable positive primaries;#eta;p_{T} (GeV/c)",
622c9e7f 156 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 157 fOutputList->Add(fHistFindableEtaPtPlusControl);
158 fHistReconstructedEtaPtPlusControl = new TH2F("fHistReconstructedEtaPtPlusControl",
159 "Reconstructed positive primaries;#eta;p_{T} (GeV/c)",
622c9e7f 160 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 161 fOutputList->Add(fHistReconstructedEtaPtPlusControl);
162 fHistSurvivedEtaPtPlusControl = new TH2F("fHistSurvivedEtaPtPlusControl",
163 "Survived positive primaries;#eta;p_{T} (GeV/c)",
622c9e7f 164 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 165 fOutputList->Add(fHistSurvivedEtaPtPlusControl);
166
167 //eta vs pt for MC negatives (control)
168 fHistGeneratedEtaPtMinusControl = new TH2F("fHistGeneratedEtaPtMinusControl",
169 "Generated positive primaries;#eta;p_{T} (GeV/c)",
622c9e7f 170 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 171 fOutputList->Add(fHistGeneratedEtaPtMinusControl);
172 fHistFindableEtaPtMinusControl = new TH2F("fHistFindableEtaPtMinusControl",
173 "Findable positive primaries;#eta;p_{T} (GeV/c)",
622c9e7f 174 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 175 fOutputList->Add(fHistFindableEtaPtMinusControl);
176 fHistReconstructedEtaPtMinusControl = new TH2F("fHistReconstructedEtaPtMinusControl",
177 "Reconstructed positive primaries;#eta;p_{T} (GeV/c)",
622c9e7f 178 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 179 fOutputList->Add(fHistReconstructedEtaPtMinusControl);
180 fHistSurvivedEtaPtMinusControl = new TH2F("fHistSurvivedEtaPtMinusControl",
181 "Survived positive primaries;#eta;p_{T} (GeV/c)",
622c9e7f 182 fEtaBin,fMinEta,fMaxEta,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 183 fOutputList->Add(fHistSurvivedEtaPtMinusControl);
184
185 //eta vs pt for MC ++
186 fHistGeneratedEtaPtPlusPlus = new TH2F("fHistGeneratedEtaPtPlusPlus",
187 "Generated ++ primaries;#Delta#eta;p_{T} (GeV/c)",
d9906391 188 fdEtaBin,fEtaRangeMin,fEtaRangeMax,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 189 fOutputList->Add(fHistGeneratedEtaPtPlusPlus);
190 fHistFindableEtaPtPlusPlus = new TH2F("fHistFindableEtaPtPlusPlus",
191 "Findable ++ primaries;#Delta#eta;p_{T} (GeV/c)",
d9906391 192 fdEtaBin,fEtaRangeMin,fEtaRangeMax,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 193 fOutputList->Add(fHistFindableEtaPtPlusPlus);
194 fHistReconstructedEtaPtPlusPlus = new TH2F("fHistReconstructedEtaPtPlusPlus",
195 "Reconstructed ++ primaries;#Delta#eta;p_{T} (GeV/c)",
d9906391 196 fdEtaBin,fEtaRangeMin,fEtaRangeMax,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 197 fOutputList->Add(fHistReconstructedEtaPtPlusPlus);
198 fHistSurvivedEtaPtPlusPlus = new TH2F("fHistSurvivedEtaPtPlusPlus",
199 "Survived ++ primaries;#Delta#eta;p_{T} (GeV/c)",
d9906391 200 fdEtaBin,fEtaRangeMin,fEtaRangeMax,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 201 fOutputList->Add(fHistSurvivedEtaPtPlusPlus);
202
203 //eta vs pt for MC --
204 fHistGeneratedEtaPtMinusMinus = new TH2F("fHistGeneratedEtaPtMinusMinus",
205 "Generated -- primaries;#Delta#eta;p_{T} (GeV/c)",
d9906391 206 fdEtaBin,fEtaRangeMin,fEtaRangeMax,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 207 fOutputList->Add(fHistGeneratedEtaPtMinusMinus);
208 fHistFindableEtaPtMinusMinus = new TH2F("fHistFindableEtaPtMinusMinus",
209 "Findable -- primaries;#Delta#eta;p_{T} (GeV/c)",
d9906391 210 fdEtaBin,fEtaRangeMin,fEtaRangeMax,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 211 fOutputList->Add(fHistFindableEtaPtMinusMinus);
212 fHistReconstructedEtaPtMinusMinus = new TH2F("fHistReconstructedEtaPtMinusMinus",
213 "Reconstructed -- primaries;#Delta#eta;p_{T} (GeV/c)",
d9906391 214 fdEtaBin,fEtaRangeMin,fEtaRangeMax,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 215 fOutputList->Add(fHistReconstructedEtaPtMinusMinus);
216 fHistSurvivedEtaPtMinusMinus = new TH2F("fHistSurvivedEtaPtMinusMinus",
217 "Survived -- primaries;#Delta#eta;p_{T} (GeV/c)",
d9906391 218 fdEtaBin,fEtaRangeMin,fEtaRangeMax,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 219 fOutputList->Add(fHistSurvivedEtaPtMinusMinus);
220
221 //eta vs pt for MC +-
222 fHistGeneratedEtaPtPlusMinus = new TH2F("fHistGeneratedEtaPtPlusMinus",
223 "Generated +- primaries;#Delta#eta;p_{T} (GeV/c)",
d9906391 224 fdEtaBin,fEtaRangeMin,fEtaRangeMax,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 225 fOutputList->Add(fHistGeneratedEtaPtPlusMinus);
226 fHistFindableEtaPtPlusMinus = new TH2F("fHistFindableEtaPtPlusMinus",
227 "Findable +- primaries;#Delta#eta;p_{T} (GeV/c)",
d9906391 228 fdEtaBin,fEtaRangeMin,fEtaRangeMax,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 229 fOutputList->Add(fHistFindableEtaPtPlusMinus);
230 fHistReconstructedEtaPtPlusMinus = new TH2F("fHistReconstructedEtaPtPlusMinus",
231 "Reconstructed +- primaries;#Delta#eta;p_{T} (GeV/c)",
d9906391 232 fdEtaBin,fEtaRangeMin,fEtaRangeMax,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 233 fOutputList->Add(fHistReconstructedEtaPtPlusMinus);
234 fHistSurvivedEtaPtPlusMinus = new TH2F("fHistSurvivedEtaPtPlusMinus",
235 "Survived +- primaries;#Delta#eta;p_{T} (GeV/c)",
d9906391 236 fdEtaBin,fEtaRangeMin,fEtaRangeMax,fPtBin,fPtRangeMin,fPtRangeMax);
84164254 237 fOutputList->Add(fHistSurvivedEtaPtPlusMinus);
238
40b41983 239 //=============================//
240 //phi vs eta for MC ++
241 fHistGeneratedPhiEtaPlusPlus = new TH2F("fHistGeneratedPhiEtaPlusPlus",
242 "Generated ++ primaries;#Delta#phi",
d9906391 243 fdPhiBin,fPhiRangeMin,fdPhiRangeMax,fdEtaBin,fEtaRangeMin,fEtaRangeMax);
40b41983 244 fOutputList->Add(fHistGeneratedPhiEtaPlusPlus);
245 fHistFindablePhiEtaPlusPlus = new TH2F("fHistFindablePhiEtaPlusPlus",
246 "Findable ++ primaries;#Delta#phi;#Delta#eta",
d9906391 247 fdPhiBin,fPhiRangeMin,fdPhiRangeMax,fdEtaBin,fEtaRangeMin,fEtaRangeMax);
40b41983 248 fOutputList->Add(fHistFindablePhiEtaPlusPlus);
249 fHistReconstructedPhiEtaPlusPlus = new TH2F("fHistReconstructedPhiEtaPlusPlus",
250 "Reconstructed ++ primaries;#Delta#phi;#Delta#eta",
d9906391 251 fdPhiBin,fPhiRangeMin,fdPhiRangeMax,fdEtaBin,fEtaRangeMin,fEtaRangeMax);
40b41983 252 fOutputList->Add(fHistReconstructedPhiEtaPlusPlus);
253 fHistSurvivedPhiEtaPlusPlus = new TH2F("fHistSurvivedPhiEtaPlusPlus",
254 "Survived ++ primaries;#Delta#phi;#Delta#eta",
d9906391 255 fdPhiBin,fPhiRangeMin,fdPhiRangeMax,fdEtaBin,fEtaRangeMin,fEtaRangeMax);
40b41983 256 fOutputList->Add(fHistSurvivedPhiEtaPlusPlus);
257
258 //phi vs eta for MC --
259 fHistGeneratedPhiEtaMinusMinus = new TH2F("fHistGeneratedPhiEtaMinusMinus",
260 "Generated -- primaries;#Delta#phi;#Delta#eta",
d9906391 261 fdPhiBin,fPhiRangeMin,fdPhiRangeMax,fdEtaBin,fEtaRangeMin,fEtaRangeMax);
40b41983 262 fOutputList->Add(fHistGeneratedPhiEtaMinusMinus);
263 fHistFindablePhiEtaMinusMinus = new TH2F("fHistFindablePhiEtaMinusMinus",
264 "Findable -- primaries;#Delta#phi;#Delta#eta",
d9906391 265 fdPhiBin,fPhiRangeMin,fdPhiRangeMax,fdEtaBin,fEtaRangeMin,fEtaRangeMax);
40b41983 266 fOutputList->Add(fHistFindablePhiEtaMinusMinus);
267 fHistReconstructedPhiEtaMinusMinus = new TH2F("fHistReconstructedPhiEtaMinusMinus",
268 "Reconstructed -- primaries;#Delta#phi;#Delta#eta",
d9906391 269 fdPhiBin,fPhiRangeMin,fdPhiRangeMax,fdEtaBin,fEtaRangeMin,fEtaRangeMax);
40b41983 270 fOutputList->Add(fHistReconstructedPhiEtaMinusMinus);
271 fHistSurvivedPhiEtaMinusMinus = new TH2F("fHistSurvivedPhiEtaMinusMinus",
272 "Survived -- primaries;#Delta#phi;#Delta#eta",
d9906391 273 fdPhiBin,fPhiRangeMin,fdPhiRangeMax,fdEtaBin,fEtaRangeMin,fEtaRangeMax);
40b41983 274 fOutputList->Add(fHistSurvivedPhiEtaMinusMinus);
275
276 //phi vs eta for MC +-
277 fHistGeneratedPhiEtaPlusMinus = new TH2F("fHistGeneratedPhiEtaPlusMinus",
278 "Generated +- primaries;#Delta#phi;#Delta#eta",
d9906391 279 fdPhiBin,fPhiRangeMin,fdPhiRangeMax,fdEtaBin,fEtaRangeMin,fEtaRangeMax);
40b41983 280 fOutputList->Add(fHistGeneratedPhiEtaPlusMinus);
281 fHistFindablePhiEtaPlusMinus = new TH2F("fHistFindablePhiEtaPlusMinus",
282 "Findable +- primaries;#Delta#phi;#Delta#eta",
d9906391 283 fdPhiBin,fPhiRangeMin,fdPhiRangeMax,fdEtaBin,fEtaRangeMin,fEtaRangeMax);
40b41983 284 fOutputList->Add(fHistFindablePhiEtaPlusMinus);
285 fHistReconstructedPhiEtaPlusMinus = new TH2F("fHistReconstructedPhiEtaPlusMinus",
286 "Reconstructed +- primaries;#Delta#phi;#Delta#eta",
d9906391 287 fdPhiBin,fPhiRangeMin,fdPhiRangeMax,fdEtaBin,fEtaRangeMin,fEtaRangeMax);
40b41983 288 fOutputList->Add(fHistReconstructedPhiEtaPlusMinus);
289 fHistSurvivedPhiEtaPlusMinus = new TH2F("fHistSurvivedPhiEtaPlusMinus",
290 "Survived +- primaries;#Delta#phi;#Delta#eta",
d9906391 291 fdPhiBin,fPhiRangeMin,fdPhiRangeMax,fdEtaBin,fEtaRangeMin,fEtaRangeMax);
40b41983 292 fOutputList->Add(fHistSurvivedPhiEtaPlusMinus);
293 //=============================//
294
84164254 295 fQAList->Print();
296 fOutputList->Print();
297
298 PostData(1, fQAList);
299 PostData(2, fOutputList);
211b716d 300
301 TH1::AddDirectory(oldStatus);
302
84164254 303}
304
305//________________________________________________________________________
306void AliAnalysisTaskEfficiencyBF::UserExec(Option_t *) {
307 // Main loop
308 // Called for each event
309
310 // Post output data.
311 //ESD analysis
312 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
313 if (!fESD) {
314 printf("ERROR: fESD not available\n");
315 return;
316 }
317
318 AliMCEvent* mcEvent = MCEvent();
319 if (!mcEvent) {
6d6eb2df 320 AliError("ERROR: Could not retrieve MC event");
84164254 321 return;
322 }
323 AliStack* stack = mcEvent->Stack();
324 if (!stack) {
6d6eb2df 325 AliError("ERROR: Could not retrieve MC stack");
84164254 326 return;
327 }
328
329 // arrays for 2 particle histograms
330 Int_t nMCLabelCounter = 0;
331 const Int_t maxMCLabelCounter = 20000;
332
84164254 333 Double_t eta[maxMCLabelCounter];
40b41983 334 Double_t pt[maxMCLabelCounter];
335 Double_t phi[maxMCLabelCounter];
84164254 336 Int_t level[maxMCLabelCounter];
337 Int_t charge[maxMCLabelCounter];
338
84164254 339 //AliInfo(Form("%d %d",mcEvent->GetNumberOfTracks(),fESD->GetNumberOfTracks()));
340 fHistEventStats->Fill(1); //all events
341
342 //Centrality stuff
343 AliCentrality *centrality = fESD->GetCentrality();
344 Int_t nCentrality = 0;
345 nCentrality = (Int_t)(centrality->GetCentralityPercentile(fCentralityEstimator.Data())/10.);
346
347 //Printf("Centrality: %lf",centrality->GetCentralityPercentile(fCentralityEstimator.Data()));
348
349 if(centrality->IsEventInCentralityClass(fCentralityPercentileMin,
350 fCentralityPercentileMax,
351 fCentralityEstimator.Data())) {
352 fHistEventStats->Fill(2); //triggered + centrality
353 fHistCentrality->Fill(nCentrality+1);
354
355 //Printf("Centrality selection: %lf - %lf",fCentralityPercentileMin,fCentralityPercentileMax);
356
357 if(fAnalysisMode.CompareTo("TPC") == 0 ) {
358 const AliESDVertex *vertex = fESD->GetPrimaryVertexTPC();
359 if(vertex) {
360 if(vertex->GetNContributors() > 0) {
361 if(vertex->GetZRes() != 0) {
362 fHistEventStats->Fill(3); //events with a proper vertex
e690d4d0 363 if(TMath::Abs(vertex->GetX()) < fVxMax) {
364 if(TMath::Abs(vertex->GetY()) < fVyMax) {
365 if(TMath::Abs(vertex->GetZ()) < fVzMax) {
84164254 366 fHistEventStats->Fill(4); //analyzed events
367
368 Int_t nMCParticles = mcEvent->GetNumberOfTracks();
369 TArrayI labelMCArray(nMCParticles);
370
371 for (Int_t iTracks = 0; iTracks < mcEvent->GetNumberOfTracks(); iTracks++) {
372 AliMCParticle *mcTrack = (AliMCParticle*) mcEvent->GetTrack(iTracks);
373 if (!mcTrack) {
6d6eb2df 374 AliError(Form("ERROR: Could not receive track %d (mc loop)", iTracks));
84164254 375 continue;
376 }
377
378 //exclude particles generated out of the acceptance
379 Double_t vz = mcTrack->Zv();
380 if (TMath::Abs(vz) > 50.) continue;
40b41983 381
84164254 382 //acceptance
622c9e7f 383 if(TMath::Abs(mcTrack->Eta()) > fMaxEta)
40b41983 384 continue;
622c9e7f 385 if((mcTrack->Pt() > fMaxPt)||(mcTrack->Pt() < fMinPt))
84164254 386 continue;
40b41983 387 if((mcTrack->Phi() > fPhiRangeMax)||(mcTrack->Phi() < fPhiRangeMin))
84164254 388 continue;
389
390 TParticle* particle = mcTrack->Particle();
391 if(!particle) continue;
392 if(!stack->IsPhysicalPrimary(iTracks)) continue;
393
394 if(iTracks <= stack->GetNprimary()) {
395 Short_t gMCCharge = mcTrack->Charge();
77e61e9e 396 Double_t phiRad = particle->Phi();
397 Double_t phiDeg = phiRad*TMath::RadToDeg();
398
84164254 399 if(gMCCharge > 0)
40b41983 400 fHistGeneratedEtaPtPhiPlus->Fill(particle->Eta(),
401 particle->Pt(),
77e61e9e 402 phiDeg);
84164254 403 else if(gMCCharge < 0)
40b41983 404 fHistGeneratedEtaPtPhiMinus->Fill(particle->Eta(),
405 particle->Pt(),
77e61e9e 406 phiDeg);
84164254 407
408
409 // findable tracks --> DOES NOT WORK????
410 // Loop over Track References
411 Bool_t labelTPC = kTRUE;
412 AliTrackReference* trackRef = 0;
413
414 for (Int_t iTrackRef = 0; iTrackRef < mcTrack->GetNumberOfTrackReferences(); iTrackRef++) {
415 trackRef = mcTrack->GetTrackReference(iTrackRef);
416 if(trackRef) {
417 Int_t detectorId = trackRef->DetectorId();
418 if (detectorId == AliTrackReference::kTPC) {
419 labelTPC = kTRUE;
420 break;
421 }
422 }
423 }//loop over track references
424
425 if(labelTPC) {
426 labelMCArray.AddAt(iTracks,nMCLabelCounter);
427
428 if(nMCLabelCounter >= maxMCLabelCounter){
429 AliWarning(Form("MC Label Counter > Limit (%d) --> stop loop here",maxMCLabelCounter));
430 break;
431 }
432
40b41983 433 //fill the arrays for 2 particle analysis
84164254 434 eta[nMCLabelCounter] = particle->Eta();
435 pt[nMCLabelCounter] = particle->Pt();
77e61e9e 436 phi[nMCLabelCounter] = particle->Phi()*TMath::RadToDeg();
84164254 437 charge[nMCLabelCounter] = gMCCharge;
77e61e9e 438
84164254 439 // findable = generated in this case!
40b41983 440
84164254 441 level[nMCLabelCounter] = 1;
84164254 442 nMCLabelCounter += 1;
40b41983 443
444
84164254 445 if(gMCCharge > 0)
40b41983 446 fHistFindableEtaPtPhiPlus->Fill(particle->Eta(),
447 particle->Pt(),
77e61e9e 448 phiDeg);
84164254 449 else if(gMCCharge < 0)
40b41983 450 fHistFindableEtaPtPhiMinus->Fill(particle->Eta(),
451 particle->Pt(),
77e61e9e 452 phiDeg);
84164254 453 }
454 }//primaries
455 }//loop over MC particles
456
457 fHistNMult->Fill(nMCLabelCounter);
458
459 // not used so far
460 //Float_t dcaXY = 0.0, dcaZ = 0.0;
461
462 //ESD track loop
463 Int_t nGoodTracks = fESD->GetNumberOfTracks();
464
465 TArrayI labelArray(nGoodTracks);
466 Int_t labelCounter = 0;
467 for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) {
468 AliESDtrack* track = fESD->GetTrack(iTracks);
469 //AliESDtrack* track = fESDtrackCuts->GetTPCOnlyTrack(fESD,iTracks);
470 if(!track) continue;
471
472 AliESDtrack *tpcOnlyTrack = new AliESDtrack();
473
474 if (!track->FillTPCOnlyTrack(*tpcOnlyTrack)) {
475 delete tpcOnlyTrack;
476 continue;
477 }
478
479 Int_t label = TMath::Abs(track->GetTPCLabel());
480 if(IsLabelUsed(labelArray,label)) continue;
481 labelArray.AddAt(label,labelCounter);
482 labelCounter += 1;
483
84164254 484 Int_t mcGoods = nMCLabelCounter;
485 for (Int_t k = 0; k < mcGoods; k++) {
486 Int_t mcLabel = labelMCArray.At(k);
84164254 487
488 if (mcLabel != TMath::Abs(label)) continue;
489 if(mcLabel != label) continue;
490 if(label > stack->GetNtrack()) continue;
491
492 TParticle *particle = stack->Particle(label);
493 if(!particle) continue;
494
495 //acceptance
622c9e7f 496 if(TMath::Abs(particle->Eta()) > fMaxEta)
40b41983 497 continue;
622c9e7f 498 if((particle->Pt() > fMaxPt)||(particle->Pt() < fMinPt))
84164254 499 continue;
40b41983 500 if((particle->Phi() > fPhiRangeMax)||(particle->Phi() < fPhiRangeMin))
84164254 501 continue;
40b41983 502
84164254 503 if(!stack->IsPhysicalPrimary(label)) continue;
504
505 if(label <= stack->GetNprimary()) {
506
507 // reconstructed
508 level[k] = 2;
509
510 Short_t gCharge = track->Charge();
77e61e9e 511 Double_t phiRad = particle->Phi();
512 Double_t phiDeg = phiRad*TMath::RadToDeg();
513
84164254 514 if(gCharge > 0)
40b41983 515 fHistReconstructedEtaPtPhiPlus->Fill(particle->Eta(),
516 particle->Pt(),
77e61e9e 517 phiDeg);
84164254 518 else if(gCharge < 0)
40b41983 519 fHistReconstructedEtaPtPhiMinus->Fill(particle->Eta(),
520 particle->Pt(),
77e61e9e 521 phiDeg);
84164254 522
523 // track cuts + analysis kinematic cuts
524 if(fESDtrackCuts->AcceptTrack(track) && TMath::Abs(track->Eta()) < fMaxEta && track->Pt() > fMinPt && track->Pt() < fMaxPt ){
525
526 // survived
527 level[k] = 3;
528
529 if(gCharge > 0)
40b41983 530 fHistSurvivedEtaPtPhiPlus->Fill(particle->Eta(),
531 particle->Pt(),
77e61e9e 532 phiDeg);
84164254 533 else if(gCharge < 0)
40b41983 534 fHistSurvivedEtaPtPhiMinus->Fill(particle->Eta(),
535 particle->Pt(),
77e61e9e 536 phiDeg);
40b41983 537
84164254 538 }//track cuts
539 }//primary particles
540 }//findable track loop
541 }//ESD track loop
40b41983 542
84164254 543 labelMCArray.Reset();
544 labelArray.Reset();
40b41983 545
546
84164254 547 }//Vz cut
548 }//Vy cut
549 }//Vx cut
550 }//Vz resolution
551 }//number of contributors
552 }//valid vertex
553 }//TPC analysis mode
554 }//centrality
555
556
557
558 // Here comes the 2 particle analysis
559 // loop over all good MC particles
560 for (Int_t i = 0; i < nMCLabelCounter ; i++) {
40b41983 561
84164254 562 // control 1D histograms (charge might be different?)
563 if(charge[i] > 0){
564 if(level[i] > 0) fHistGeneratedEtaPtPlusControl->Fill(eta[i],pt[i]);
565 if(level[i] > 1) fHistReconstructedEtaPtPlusControl->Fill(eta[i],pt[i]);
566 if(level[i] > 2) fHistSurvivedEtaPtPlusControl->Fill(eta[i],pt[i]);
567 }
568 else if(charge[i] < 0){
569 if(level[i] > 0) fHistGeneratedEtaPtMinusControl->Fill(eta[i],pt[i]);
570 if(level[i] > 1) fHistReconstructedEtaPtMinusControl->Fill(eta[i],pt[i]);
571 if(level[i] > 2) fHistSurvivedEtaPtMinusControl->Fill(eta[i],pt[i]);
572 }
40b41983 573
574
84164254 575 for (Int_t j = i+1; j < nMCLabelCounter ; j++) {
576
577 if(charge[i] > 0 && charge[j] > 0 ){
40b41983 578 if(level[i] > 0 && level[j] > 0) {
579 fHistGeneratedEtaPtPlusPlus->Fill(TMath::Abs(eta[i]-eta[j]),pt[i]);
77e61e9e 580 if (TMath::Abs(phi[i]-phi[j]) < 180)
40b41983 581 fHistGeneratedPhiEtaPlusPlus->Fill(TMath::Abs(phi[i]-phi[j]),TMath::Abs(eta[i]-eta[j]));
582 }
583 if(level[i] > 1 && level[j] > 1) {
584 fHistReconstructedEtaPtPlusPlus->Fill(TMath::Abs(eta[i]-eta[j]),pt[i]);
77e61e9e 585 if (TMath::Abs(phi[i]-phi[j]) < 180)
40b41983 586 fHistReconstructedPhiEtaPlusPlus->Fill(TMath::Abs(phi[i]-phi[j]),TMath::Abs(eta[i]-eta[j]));
587 }
588 if(level[i] > 2 && level[j] > 2) {
589 fHistSurvivedEtaPtPlusPlus->Fill(TMath::Abs(eta[i]-eta[j]),pt[i]);
77e61e9e 590 if (TMath::Abs(phi[i]-phi[j]) < 180)
40b41983 591 fHistSurvivedPhiEtaPlusPlus->Fill(TMath::Abs(phi[i]-phi[j]),TMath::Abs(eta[i]-eta[j]));
592 }
84164254 593 }
594
595 else if(charge[i] < 0 && charge[j] < 0 ){
40b41983 596 if(level[i] > 0 && level[j] > 0) {
597 fHistGeneratedEtaPtMinusMinus->Fill(TMath::Abs(eta[i]-eta[j]),pt[i]);
77e61e9e 598 if (TMath::Abs(phi[i]-phi[j]) < 180)
40b41983 599 fHistGeneratedPhiEtaMinusMinus->Fill(TMath::Abs(phi[i]-phi[j]),TMath::Abs(eta[i]-eta[j]));
600 }
601 if(level[i] > 1 && level[j] > 1) {
602 fHistReconstructedEtaPtMinusMinus->Fill(TMath::Abs(eta[i]-eta[j]),pt[i]);
603 fHistReconstructedPhiEtaMinusMinus->Fill(TMath::Abs(phi[i]-phi[j]),TMath::Abs(eta[i]-eta[j]));
604 }
605 if(level[i] > 2 && level[j] > 2) {
77e61e9e 606 fHistSurvivedEtaPtMinusMinus->Fill(TMath::Abs(eta[i]-eta[j]),pt[i]);
607 if (TMath::Abs(phi[i]-phi[j]) < 180)
40b41983 608 fHistSurvivedPhiEtaMinusMinus->Fill(TMath::Abs(phi[i]-phi[j]),TMath::Abs(eta[i]-eta[j]));
609 }
84164254 610 }
40b41983 611
84164254 612 else if((charge[i] > 0 && charge[j] < 0)||(charge[i] < 0 && charge[j] > 0)){
40b41983 613 if(level[i] > 0 && level[j] > 0) {
614 fHistGeneratedEtaPtPlusMinus->Fill(TMath::Abs(eta[i]-eta[j]),pt[i]);
77e61e9e 615 if (TMath::Abs(phi[i]-phi[j]) < 180)
40b41983 616 fHistGeneratedPhiEtaPlusMinus->Fill(TMath::Abs(phi[i]-phi[j]),TMath::Abs(eta[i]-eta[j]));
617 }
618 if(level[i] > 1 && level[j] > 1) {
619 fHistReconstructedEtaPtPlusMinus->Fill(TMath::Abs(eta[i]-eta[j]),pt[i]);
620 fHistReconstructedPhiEtaPlusMinus->Fill(TMath::Abs(phi[i]-phi[j]),TMath::Abs(eta[i]-eta[j]));
621 }
622 if(level[i] > 2 && level[j] > 2) {
623 fHistSurvivedEtaPtPlusMinus->Fill(TMath::Abs(eta[i]-eta[j]),pt[i]);
77e61e9e 624 if (TMath::Abs(phi[i]-phi[j]) < 180)
40b41983 625 fHistSurvivedPhiEtaPlusMinus->Fill(TMath::Abs(phi[i]-phi[j]),TMath::Abs(eta[i]-eta[j]));
626 }
84164254 627 }
628 }
40b41983 629 }
630
631 //Checking Entries for generated and survived
632 /*TH1F* GeneratedEtaPt = (TH1F*)fHistGeneratedEtaPtPlusMinus->ProjectionX("GeneratedEtaPt",5,15);
633 Int_t xGeneratedPt = fHistGeneratedEtaPtPlusMinus->GetNbinsX();
634 for (Int_t h=1; h < xGeneratedPt+1; h++){
635 Double_t binEntriesGenerated = GeneratedEtaPt->GetBinContent(h);
636 Printf("binEntriesGenerated: %lf - xGeneratedPt: %d",binEntriesGenerated,h);
637 }
638 TH1F* GeneratedPhiEta = (TH1F*)fHistGeneratedPhiEtaPlusMinus->ProjectionY("GeneratedPhiEta",5,15);
639 Int_t yGeneratedPhi = fHistGeneratedPhiEtaPlusMinus->GetNbinsY();
640 for (Int_t h=1; h < yGeneratedPhi+1; h++){
641 Double_t binEntriesGenerated = GeneratedPhiEta->GetBinContent(h);
642 Printf("binEntriesGenerated: %lf - yGeneratedPhi: %d",binEntriesGenerated,h);
643 }*/
644
645 /*TH1F* SurvivedEtaPt = (TH1F*)fHistSurvivedEtaPtPlusMinus->ProjectionX("SurvivedEtaPt",5,15);
646 Int_t xSurvivedPt = fHistSurvivedEtaPtPlusMinus->GetNbinsX();
647 for (Int_t h=1; h < xSurvivedPt+1; h++){
648 Double_t binEntriesSurvived = SurvivedEtaPt->GetBinContent(h);
649 Printf("binEntriesSurvived: %lf - xSurvivedPt: %d",binEntriesSurvived,h);
650 }
651 TH1F* SurvivedPhiEta = (TH1F*)fHistSurvivedPhiEtaPlusMinus->ProjectionY("SurvivedPhiEta",5,15);
652 Int_t ySurvivedPhi = fHistSurvivedPhiEtaPlusMinus->GetNbinsY();
653 for (Int_t h=1; h < ySurvivedPhi+1; h++){
654 Double_t binEntriesSurvived = SurvivedPhiEta->GetBinContent(h);
655 Printf("binEntriesSurvived: %lf - ySurvivedPhi: %d",binEntriesSurvived,h);
656 }*/
84164254 657}
40b41983 658 //________________________________________________________________________
84164254 659void AliAnalysisTaskEfficiencyBF::Terminate(Option_t *) {
660 // Draw result to the screen
661 // Called once at the end of the query
662
663 fOutputList = dynamic_cast<TList*> (GetOutputData(1));
664 if (!fOutputList) {
665 printf("ERROR: Output list not available\n");
666 return;
667 }
668}
669
670//____________________________________________________________________//
671Bool_t AliAnalysisTaskEfficiencyBF::IsLabelUsed(TArrayI labelArray, Int_t label) {
672 //Checks if the label is used already
673 Bool_t status = kFALSE;
674 for(Int_t i = 0; i < labelArray.GetSize(); i++) {
675 if(labelArray.At(i) == label)
676 status = kTRUE;
677 }
678
679 return status;
680}