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