]>
Commit | Line | Data |
---|---|---|
48a61f36 | 1 | // Dihadron correlations task - simple task to read ESD or AOD input, |
2 | // calculate same- and mixed-event correlations, and fill THnSparse | |
3 | // output. -A. Adare, Apr 2011 | |
4 | ||
5 | #include "TCanvas.h" | |
6 | #include "TChain.h" | |
8433ff41 | 7 | #include "TFormula.h" |
48a61f36 | 8 | #include "TH1.h" |
9 | #include "TH2.h" | |
beb13b6c | 10 | #include "TH3.h" |
d688e049 | 11 | #include "TAxis.h" |
8433ff41 | 12 | #include "TProfile2D.h" |
48a61f36 | 13 | #include "TROOT.h" |
14 | #include "TTree.h" | |
48a61f36 | 15 | #include "AliAODEvent.h" |
16 | #include "AliAODInputHandler.h" | |
17 | #include "AliAnalysisManager.h" | |
18 | #include "AliAnalysisTask.h" | |
19 | #include "AliCentrality.h" | |
20 | #include "AliDhcTask.h" | |
21 | #include "AliESDEvent.h" | |
22 | #include "AliESDInputHandler.h" | |
23 | #include "AliESDtrackCuts.h" | |
beb13b6c | 24 | #include "AliESDMuonTrack.h" |
8433ff41 | 25 | #include "AliPool.h" |
48a61f36 | 26 | #include "AliVParticle.h" |
27 | ||
c64cb1f6 | 28 | using std::cout; |
29 | using std::endl; | |
30 | ||
48a61f36 | 31 | ClassImp(AliDhcTask) |
32 | ||
beb13b6c | 33 | |
34 | //________________________________________________________________________ | |
35 | AliDhcTask::AliDhcTask() | |
36 | : AliAnalysisTaskSE(), fVerbosity(0), fEtaMax(1), fZVtxMax(10), fPtMin(0.25), fPtMax(15), | |
37 | fTrackDepth(1000), fPoolSize(200), fTracksName(), fDoWeights(kFALSE), fFillMuons(kFALSE), | |
b422cf0d | 38 | fPtTACrit(kTRUE), |
beb13b6c | 39 | fEtaTLo(-1.0), fEtaTHi(1.0), fEtaALo(-1.0), fEtaAHi(1.0), |
40 | fESD(0x0), fAOD(0x0), fOutputList(0x0), fHEvt(0x0), fHTrk(0x0), | |
41 | fHPtAss(0x0), fHPtTrg(0x0), fHPtTrgEvt(0x0), | |
42 | fHPtTrgNorm1S(0x0), fHPtTrgNorm1M(0x0), fHPtTrgNorm2S(0x0), fHPtTrgNorm2M(0x0), | |
380fb711 | 43 | fHCent(0x0), fHZvtx(0x0), fNbins(0), fHSs(0x0), fHMs(0x0), fHPts(0x0), |
44 | fHQAT(0x0), fHQAA(0x0), | |
45 | fIndex(0x0), | |
beb13b6c | 46 | fCentrality(99), fZVertex(99), fEsdTPCOnly(0), fPoolMgr(0), |
47 | fCentMethod("V0M"), fNBdeta(20), fNBdphi(36), | |
48 | fBPtT(0x0), fBPtA(0x0), fBCent(0x0), fBZvtx(0x0), | |
49 | fMixBCent(0x0), fMixBZvtx(0x0), | |
50 | fHEffT(0x0), fHEffA(0x0) | |
51 | { | |
52 | ||
53 | } | |
54 | ||
48a61f36 | 55 | //________________________________________________________________________ |
56 | AliDhcTask::AliDhcTask(const char *name) | |
42036329 | 57 | : AliAnalysisTaskSE(name), fVerbosity(0), fEtaMax(1), fZVtxMax(10), fPtMin(0.25), fPtMax(15), |
beb13b6c | 58 | fTrackDepth(1000), fPoolSize(200), fTracksName(), fDoWeights(kFALSE), fFillMuons(kFALSE), |
b422cf0d | 59 | fPtTACrit(kTRUE), |
beb13b6c | 60 | fEtaTLo(-1.0), fEtaTHi(1.0), fEtaALo(-1.0), fEtaAHi(1.0), |
61 | fESD(0x0), fAOD(0x0), fOutputList(0x0), fHEvt(0x0), fHTrk(0x0), | |
62 | fHPtAss(0x0), fHPtTrg(0x0), fHPtTrgEvt(0x0), | |
63 | fHPtTrgNorm1S(0x0), fHPtTrgNorm1M(0x0), fHPtTrgNorm2S(0x0), fHPtTrgNorm2M(0x0), | |
380fb711 | 64 | fHCent(0x0), fHZvtx(0x0), fNbins(0), fHSs(0x0), fHMs(0x0), fHPts(0x0), |
65 | fHQAT(0x0), fHQAA(0x0), | |
66 | fIndex(0x0), | |
8a9d3e12 | 67 | fCentrality(99), fZVertex(99), fEsdTPCOnly(0), fPoolMgr(0), |
d688e049 | 68 | fCentMethod("V0M"), fNBdeta(20), fNBdphi(36), |
69 | fBPtT(0x0), fBPtA(0x0), fBCent(0x0), fBZvtx(0x0), | |
beb13b6c | 70 | fMixBCent(0x0), fMixBZvtx(0x0), |
71 | fHEffT(0x0), fHEffA(0x0) | |
48a61f36 | 72 | { |
73 | // Constructor | |
74 | ||
75 | // Define input and output slots here | |
76 | // Input slot #0 works with a TChain | |
77 | DefineInput(0, TChain::Class()); | |
78 | // Output slot #0 id reserved by the base class for AOD | |
79 | // Output slot #1 writes into a TH1 container | |
80 | DefineOutput(1, TList::Class()); | |
81 | ||
48a61f36 | 82 | fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.,SPDVertex.,TPCVertex.,Tracks " |
83 | "AOD:header,tracks,vertices,"; | |
d688e049 | 84 | |
85 | Double_t ptt[4] = {0.25, 1.0, 2.0, 15.0}; | |
86 | fBPtT = new TAxis(3,ptt); | |
87 | Double_t pta[4] = {0.25, 1.0, 2.0, 15.0}; | |
88 | fBPtA = new TAxis(3,pta); | |
89 | Double_t cent[2] = {-100.0, 100.0}; | |
90 | fBCent = new TAxis(1,cent); | |
91 | Double_t zvtx[2] = {-10, 10}; | |
92 | fBZvtx = new TAxis(1,zvtx); | |
93 | Double_t centmix[2] = {-100.0, 100.0}; | |
94 | fMixBCent = new TAxis(1,centmix); | |
95 | Double_t zvtxmix[9] = {-10,-6,-4,-2,0,2,4,6,10}; | |
96 | fMixBZvtx = new TAxis(8,zvtxmix); | |
48a61f36 | 97 | } |
98 | ||
99 | //________________________________________________________________________ | |
100 | void AliDhcTask::UserCreateOutputObjects() | |
101 | { | |
102 | // Create histograms | |
103 | // Called once (per slave on PROOF!) | |
beb13b6c | 104 | if (fVerbosity > 1) { |
105 | AliInfo("Initialize Dhc Task"); | |
106 | AliInfo(Form(" efficiency correct triggers? %d", fHEffT!=0)); | |
107 | AliInfo(Form(" efficiency correct associates? %d", fHEffA!=0)); | |
108 | AliInfo(Form(" fill muons? %d", fFillMuons)); | |
b422cf0d | 109 | AliInfo(Form(" use pTT > pTA criterion? %d", fPtTACrit)); |
beb13b6c | 110 | AliInfo(Form(" trigger eta range %f .. %f", fEtaTLo, fEtaTHi)); |
111 | AliInfo(Form(" associate eta range %f .. %f", fEtaALo, fEtaAHi)); | |
112 | } | |
48a61f36 | 113 | |
48a61f36 | 114 | fOutputList = new TList(); |
115 | fOutputList->SetOwner(1); | |
116 | ||
8a9d3e12 | 117 | fEsdTPCOnly = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts(); |
118 | //fEsdTPCOnly->SetMinNClustersTPC(70); | |
119 | fEsdTPCOnly->SetMinNCrossedRowsTPC(70); | |
120 | fEsdTPCOnly->SetMinRatioCrossedRowsOverFindableClustersTPC(0.8); | |
48a61f36 | 121 | |
8433ff41 | 122 | BookHistos(); |
123 | InitEventMixer(); | |
48a61f36 | 124 | PostData(1, fOutputList); |
125 | } | |
126 | ||
127 | //________________________________________________________________________ | |
128 | void AliDhcTask::BookHistos() | |
129 | { | |
8433ff41 | 130 | // Book histograms. |
131 | ||
d688e049 | 132 | if (fVerbosity > 1) { |
8a9d3e12 | 133 | AliInfo(Form("Number of pt(t) bins: %d", fBPtT->GetNbins())); |
134 | for (Int_t i=1; i<=fBPtT->GetNbins(); i++) { | |
135 | AliInfo(Form("pt(t) bin %d, %f to %f", i, fBPtT->GetBinLowEdge(i), fBPtT->GetBinUpEdge(i))); | |
136 | } | |
d688e049 | 137 | AliInfo(Form("Number of pt(a) bins: %d", fBPtA->GetNbins())); |
138 | for (Int_t i=1; i<=fBPtA->GetNbins(); i++) { | |
139 | AliInfo(Form("pt(a) bin %d, %f to %f", i, fBPtA->GetBinLowEdge(i), fBPtA->GetBinUpEdge(i))); | |
140 | } | |
141 | } | |
48a61f36 | 142 | |
d688e049 | 143 | Int_t nPtAssc=fBPtA->GetNbins(); |
144 | Int_t nPtTrig=fBPtT->GetNbins(); | |
145 | Int_t nCent=fBCent->GetNbins(); | |
146 | Int_t nZvtx=fBZvtx->GetNbins(); | |
147 | Double_t ptt[nPtTrig+1]; | |
148 | ptt[0] = fBPtT->GetBinLowEdge(1); | |
149 | for (Int_t i=1; i<=nPtTrig; i++) { | |
150 | ptt[i] = fBPtT->GetBinUpEdge(i); | |
151 | } | |
152 | Double_t pta[nPtAssc+1]; | |
153 | pta[0] = fBPtA->GetBinLowEdge(1); | |
154 | for (Int_t i=1; i<=nPtAssc; i++) { | |
155 | pta[i] = fBPtA->GetBinUpEdge(i); | |
156 | } | |
157 | Double_t cent[nCent+1]; | |
158 | cent[0] = fBCent->GetBinLowEdge(1); | |
159 | for (Int_t i=1; i<=nCent; i++) { | |
160 | cent[i] = fBCent->GetBinUpEdge(i); | |
161 | } | |
162 | Double_t zvtx[nZvtx+1]; | |
163 | zvtx[0] = fBZvtx->GetBinLowEdge(1); | |
164 | for (Int_t i=1; i<=nZvtx; i++) { | |
165 | zvtx[i] = fBZvtx->GetBinUpEdge(i); | |
166 | } | |
380fb711 | 167 | |
8a9d3e12 | 168 | fHEvt = new TH2F("fHEvt", "Event-level variables; Zvtx; Cent", nZvtx, zvtx, nCent, cent); |
8433ff41 | 169 | fOutputList->Add(fHEvt); |
380fb711 | 170 | fHTrk = new TH2F("fHTrk", "Track-level variables", |
171 | 100, 0, TMath::TwoPi(), 100, -fEtaMax, +fEtaMax); | |
8433ff41 | 172 | fOutputList->Add(fHTrk); |
173 | ||
8433ff41 | 174 | fHPtAss = new TH1F("fHPtAss","PtAssoc;P_{T} (GeV/c) [GeV/c]",nPtAssc,pta); |
175 | fOutputList->Add(fHPtAss); | |
176 | fHPtTrg = new TH1F("fHPtTrg","PtTrig;P_{T} (GeV/c) [GeV/c]",nPtTrig,ptt); | |
177 | fOutputList->Add(fHPtTrg); | |
beb13b6c | 178 | fHPtTrgEvt = new TH1F("fHPtTrgEvt","PtTrig;P_{T} (GeV/c) [GeV/c]",nPtTrig,ptt); |
179 | fHPtTrgNorm1S = new TH3F("fHPtTrgNorm1S","PtTrig;P_{T} (GeV/c) [GeV/c];centrality;z_{vtx}",nPtTrig,ptt,nCent,cent,nZvtx,zvtx); | |
180 | fOutputList->Add(fHPtTrgNorm1S); | |
181 | fHPtTrgNorm1M = (TH3*) fHPtTrgNorm1S->Clone("fHPtTrgNorm1M"); | |
182 | fOutputList->Add(fHPtTrgNorm1M); | |
183 | fHPtTrgNorm2S = (TH3*) fHPtTrgNorm1S->Clone("fHPtTrgNorm2S"); | |
184 | fOutputList->Add(fHPtTrgNorm2S); | |
185 | fHPtTrgNorm2M = (TH3*) fHPtTrgNorm1S->Clone("fHPtTrgNorm2M"); | |
186 | fOutputList->Add(fHPtTrgNorm2M); | |
380fb711 | 187 | |
8433ff41 | 188 | fHCent = new TH1F("fHCent","Cent;bins",nCent,cent); |
189 | fOutputList->Add(fHCent); | |
190 | fHZvtx = new TH1F("fHZvtx","Zvertex;bins",nZvtx,zvtx); | |
191 | fOutputList->Add(fHZvtx); | |
380fb711 | 192 | |
193 | fHQAT = new TH3F("fHQAT","QA trigger;p_{T} (GeV/c);#eta;#phi", | |
194 | 100,0.0,10.0, | |
195 | 40,fEtaTLo,fEtaTHi, | |
196 | 36,0.0,TMath::TwoPi()); | |
197 | fOutputList->Add(fHQAT); | |
198 | ||
199 | fHQAA = new TH3F("fHQAA","QA associated;p_{T} (GeV/c);#eta;#phi", | |
200 | 100,0.0,10.0, | |
201 | 40,fEtaALo,fEtaAHi, | |
202 | 36,0.0,TMath::TwoPi()); | |
203 | fOutputList->Add(fHQAA); | |
8433ff41 | 204 | |
205 | fNbins = nPtTrig*nPtAssc*nCent*nZvtx; | |
8a9d3e12 | 206 | fHSs = new TH2*[fNbins]; |
207 | fHMs = new TH2*[fNbins]; | |
208 | fHPts = new TH1*[fNbins]; | |
380fb711 | 209 | |
8433ff41 | 210 | fIndex = new TFormula("GlobIndex","(t-1)*[0]*[1]*[2]+(z-1)*[0]*[1]+(x-1)*[0]+(y-1)+0*[4]"); |
211 | fIndex->SetParameters(nPtTrig,nPtAssc,nZvtx,nCent); | |
212 | fIndex->SetParNames("NTrigBins","NAssocBins", "NZvertexBins", "NCentBins"); | |
beb13b6c | 213 | fOutputList->Add(fIndex); |
380fb711 | 214 | |
8433ff41 | 215 | Int_t count = 0; |
216 | for (Int_t c=1; c<=nCent; ++c) { | |
217 | for (Int_t z=1; z<=nZvtx; ++z) { | |
218 | for (Int_t t=1; t<=nPtTrig; ++t) { | |
380fb711 | 219 | for (Int_t a=1; a<=nPtAssc; ++a) { |
220 | fHSs[count] = 0; | |
221 | fHMs[count] = 0; | |
8a9d3e12 | 222 | fHPts[count] = 0; |
380fb711 | 223 | if (a>t) { |
224 | ++count; | |
225 | continue; | |
226 | } | |
8a9d3e12 | 227 | if (t==1 && a==1) { |
228 | TString title(Form("cen=%d (%.1f to %.1f), zVtx=%d (%.1f to %.1f)", | |
380fb711 | 229 | c, fBCent->GetBinLowEdge(c), fBCent->GetBinUpEdge(c), |
8a9d3e12 | 230 | z, fBZvtx->GetBinLowEdge(z), fBZvtx->GetBinUpEdge(z))); |
231 | fHPts[count] = new TH1F(Form("hPt%d",count), Form("Ptdist %s;p_{T} (GeV/c)",title.Data()), 200,0,20); | |
232 | fOutputList->Add(fHPts[count]); | |
233 | } | |
380fb711 | 234 | TString title(Form("cen=%d (%.1f to %.1f), zVtx=%d (%.1f to %.1f), trig=%d (%.1f to %.1f), assoc=%d (%.1f to %.1f)", |
235 | c, fBCent->GetBinLowEdge(c), fBCent->GetBinUpEdge(c), | |
8a9d3e12 | 236 | z, fBZvtx->GetBinLowEdge(z), fBZvtx->GetBinUpEdge(z), |
380fb711 | 237 | t, fBPtT->GetBinLowEdge(t), fBPtT->GetBinUpEdge(t), |
8a9d3e12 | 238 | a, fBPtA->GetBinLowEdge(a), fBPtA->GetBinUpEdge(a))); |
380fb711 | 239 | fHSs[count] = new TH2F(Form("hS%d",count), Form("Signal %s;#Delta #varphi;#Delta #eta",title.Data()), |
240 | fNBdphi,-0.5*TMath::Pi(),1.5*TMath::Pi(),fNBdeta,-2*fEtaMax,2*fEtaMax); | |
241 | fHMs[count] = new TH2F(Form("hM%d",count), Form("Mixed %s;#Delta #varphi;#Delta #eta",title.Data()), | |
242 | fNBdphi,-0.5*TMath::Pi(),1.5*TMath::Pi(),fNBdeta,-2*fEtaMax,2*fEtaMax); | |
243 | fOutputList->Add(fHSs[count]); | |
244 | fOutputList->Add(fHMs[count]); | |
8a9d3e12 | 245 | Int_t tcount = (Int_t)fIndex->Eval(t,a,z,c); |
380fb711 | 246 | if (fVerbosity>5) |
247 | cout << count << " " << tcount << ": " << title << endl; | |
8a9d3e12 | 248 | if (count != tcount) |
249 | AliFatal(Form("Index mismatch: %d %d", count, tcount)); | |
380fb711 | 250 | ++count; |
251 | } | |
8433ff41 | 252 | } |
253 | } | |
254 | } | |
48a61f36 | 255 | } |
256 | ||
beb13b6c | 257 | //________________________________________________________________________ |
b422cf0d | 258 | void AliDhcTask::SetAnaMode(Int_t iAna) |
beb13b6c | 259 | { |
260 | if (iAna==kHH) { | |
261 | SetFillMuons(kFALSE); | |
262 | fEtaTLo = -1.0; | |
263 | fEtaTHi = 1.0; | |
264 | fEtaALo = -1.0; | |
265 | fEtaAHi = 1.0; | |
266 | } else if (iAna==kMuH) { | |
267 | SetFillMuons(kTRUE); | |
268 | fEtaTLo = -5.0; | |
380fb711 | 269 | fEtaTHi = -1.0; |
beb13b6c | 270 | fEtaALo = -1.0; |
271 | fEtaAHi = 1.0; | |
272 | } else if (iAna==kHMu) { | |
273 | SetFillMuons(kTRUE); | |
274 | fEtaTLo = -1.0; | |
275 | fEtaTHi = 1.0; | |
276 | fEtaALo = -5.0; | |
380fb711 | 277 | fEtaAHi = -1.0; |
beb13b6c | 278 | } else if (iAna==kMuMu) { |
279 | SetFillMuons(kTRUE); | |
280 | fEtaTLo = -5.0; | |
380fb711 | 281 | fEtaTHi = -1.0; |
beb13b6c | 282 | fEtaALo = -5.0; |
380fb711 | 283 | fEtaAHi = -1.0; |
beb13b6c | 284 | } else if (iAna==kPSide) { |
285 | SetFillMuons(kFALSE); | |
286 | fEtaTLo = -1.0; | |
287 | fEtaTHi = -0.465; | |
288 | fEtaALo = -1.0; | |
289 | fEtaAHi = -0.465; | |
290 | } else if (iAna==kASide) { | |
291 | SetFillMuons(kFALSE); | |
292 | fEtaTLo = -0.465; | |
293 | fEtaTHi = 1.0; | |
294 | fEtaALo = -0.465; | |
295 | fEtaAHi = 1.0; | |
296 | } else { | |
297 | AliInfo(Form("Unrecognized analysis option: %d", iAna)); | |
298 | } | |
299 | } | |
300 | ||
48a61f36 | 301 | //________________________________________________________________________ |
302 | void AliDhcTask::InitEventMixer() | |
303 | { | |
304 | // The effective pool size in events is set by trackDepth, so more | |
305 | // low-mult events are required to maintain the threshold than | |
306 | // high-mult events. Centrality pools are indep. of data histogram | |
307 | // binning, no need to match. | |
308 | ||
48a61f36 | 309 | // Centrality pools |
d688e049 | 310 | Int_t nCentBins=fMixBCent->GetNbins(); |
311 | Double_t centBins[nCentBins+1]; | |
312 | centBins[0] = fMixBCent->GetBinLowEdge(1); | |
313 | for (Int_t i=1; i<=nCentBins; i++) { | |
314 | centBins[i] = fMixBCent->GetBinUpEdge(i); | |
315 | } | |
48a61f36 | 316 | |
317 | // Z-vertex pools | |
d688e049 | 318 | Int_t nZvtxBins=fMixBZvtx->GetNbins(); |
319 | Double_t zvtxbin[nZvtxBins+1]; | |
320 | zvtxbin[0] = fMixBZvtx->GetBinLowEdge(1); | |
321 | for (Int_t i=1; i<=nZvtxBins; i++) { | |
322 | zvtxbin[i] = fMixBZvtx->GetBinUpEdge(i); | |
323 | } | |
8433ff41 | 324 | |
325 | fPoolMgr = new AliEvtPoolManager(); | |
d688e049 | 326 | fPoolMgr->SetTargetTrackDepth(fTrackDepth); |
8433ff41 | 327 | if (fVerbosity>4) |
328 | fPoolMgr->SetDebug(1); | |
d688e049 | 329 | fPoolMgr->InitEventPools(fPoolSize, nCentBins, centBins, nZvtxBins, zvtxbin); |
48a61f36 | 330 | } |
331 | ||
332 | //________________________________________________________________________ | |
333 | void AliDhcTask::UserExec(Option_t *) | |
334 | { | |
335 | // Main loop, called for each event. | |
48a61f36 | 336 | static int iEvent = -1; ++iEvent; |
337 | ||
338 | if (fVerbosity>2) { | |
339 | if (iEvent % 10 == 0) | |
340 | cout << iEvent << endl; | |
341 | } | |
342 | ||
343 | Int_t dType = -1; // Will be set to kESD or kAOD. | |
f6701c6e | 344 | MiniEvent* sTracks = new MiniEvent(0); // deleted by pool mgr. |
48a61f36 | 345 | |
346 | LoadBranches(); | |
347 | ||
348 | // Get event pointers, check for signs of life | |
349 | fESD = dynamic_cast<AliESDEvent*>(InputEvent()); | |
350 | fAOD = dynamic_cast<AliAODEvent*>(InputEvent()); | |
351 | if (fESD) | |
352 | dType = kESD; | |
353 | else if (fAOD) | |
354 | dType = kAOD; | |
355 | else { | |
356 | AliError("Neither fESD nor fAOD available"); | |
357 | return; | |
358 | } | |
359 | ||
d688e049 | 360 | Bool_t mcgen = 0; |
361 | if (fTracksName.Contains("Gen")) | |
362 | mcgen = 1; | |
363 | ||
48a61f36 | 364 | // Centrality, vertex, other event variables... |
d688e049 | 365 | if (mcgen) { |
366 | fZVertex = 0; | |
367 | TList *list = InputEvent()->GetList(); | |
beb13b6c | 368 | TClonesArray *tcaTracks = 0; |
369 | if (list) | |
370 | tcaTracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracksName)); | |
d688e049 | 371 | if (tcaTracks) |
372 | fCentrality = tcaTracks->GetEntries(); | |
373 | } else { | |
374 | if (dType == kESD) { | |
375 | if (!VertexOk(fESD)) { | |
376 | if (fVerbosity > 1) | |
377 | AliInfo(Form("Event REJECTED (ESD vertex not OK). z = %.1f", fZVertex)); | |
378 | return; | |
379 | } | |
380 | const AliESDVertex* vertex = fESD->GetPrimaryVertex(); | |
381 | fZVertex = vertex->GetZ(); | |
382 | if(fESD->GetCentrality()) { | |
383 | fCentrality = | |
384 | fESD->GetCentrality()->GetCentralityPercentile(fCentMethod); | |
385 | } | |
48a61f36 | 386 | } |
d688e049 | 387 | if (dType == kAOD) { |
388 | const AliAODVertex* vertex = fAOD->GetPrimaryVertex(); | |
389 | fZVertex = vertex->GetZ(); | |
390 | if (!VertexOk(fAOD)) { | |
391 | if (fVerbosity > 1) | |
392 | Info("Exec()", "Event REJECTED (AOD vertex not OK). z = %.1f", fZVertex); | |
393 | return; | |
394 | } | |
395 | const AliCentrality *aodCent = fAOD->GetHeader()->GetCentralityP(); | |
396 | if (aodCent) { | |
397 | fCentrality = aodCent->GetCentralityPercentile(fCentMethod); | |
398 | } | |
48a61f36 | 399 | } |
400 | } | |
d688e049 | 401 | |
48a61f36 | 402 | // Fill Event histogram |
403 | fHEvt->Fill(fZVertex, fCentrality); | |
d688e049 | 404 | if (fCentrality > fBCent->GetXmax() || fCentrality < fBCent->GetXmin()) { |
9370528f | 405 | if (fVerbosity > 1) |
406 | AliInfo(Form("Event REJECTED (centrality out of range). fCentrality = %.1f", fCentrality)); | |
48a61f36 | 407 | return; |
408 | } | |
409 | ||
410 | // Get array of selected tracks | |
411 | if (dType == kESD) { | |
f6701c6e | 412 | GetESDTracks(sTracks); |
48a61f36 | 413 | } |
414 | if (dType == kAOD) { | |
f6701c6e | 415 | GetAODTracks(sTracks); |
48a61f36 | 416 | } |
417 | ||
418 | // Get pool containing tracks from other events like this one | |
8433ff41 | 419 | AliEvtPool* pool = fPoolMgr->GetEventPool(fCentrality, fZVertex); |
48a61f36 | 420 | if (!pool) { |
c71ff6f4 | 421 | AliWarning(Form("No pool found. Centrality %f, ZVertex %f", |
422 | fCentrality, fZVertex)); | |
48a61f36 | 423 | return; |
424 | } | |
425 | ||
426 | if (!pool->IsReady()) { | |
427 | pool->UpdatePool(sTracks); | |
428 | return; | |
429 | } | |
430 | ||
431 | if (pool->IsFirstReady()) { | |
432 | // recover events missed before the pool is ready | |
433 | Int_t nEvs = pool->GetCurrentNEvents(); | |
434 | if (nEvs>1) { | |
435 | for (Int_t i=0; i<nEvs; ++i) { | |
436 | MiniEvent* evI = pool->GetEvent(i); | |
437 | for (Int_t j=0; j<nEvs; ++j) { | |
438 | MiniEvent* evJ = pool->GetEvent(j); | |
439 | if (i==j) { | |
440 | Correlate(*evI, *evJ, kSameEvt); | |
441 | } else { | |
d688e049 | 442 | Correlate(*evI, *evJ, kDiffEvt); |
48a61f36 | 443 | } |
444 | } | |
445 | } | |
446 | } | |
447 | } else { /* standard case: same event, then mix*/ | |
448 | Correlate(*sTracks, *sTracks, kSameEvt); | |
449 | Int_t nMix = pool->GetCurrentNEvents(); | |
450 | for (Int_t jMix=0; jMix<nMix; ++jMix) { | |
451 | MiniEvent* bgTracks = pool->GetEvent(jMix); | |
d688e049 | 452 | Correlate(*sTracks, *bgTracks, kDiffEvt); |
48a61f36 | 453 | } |
454 | } | |
455 | ||
48a61f36 | 456 | pool->UpdatePool(sTracks); |
457 | PostData(1, fOutputList); | |
48a61f36 | 458 | } |
459 | ||
460 | //________________________________________________________________________ | |
f6701c6e | 461 | void AliDhcTask::GetESDTracks(MiniEvent* miniEvt) |
48a61f36 | 462 | { |
463 | // Loop twice: 1. Count sel. tracks. 2. Fill vector. | |
464 | ||
d688e049 | 465 | if (fTracksName.IsNull()) { |
466 | const AliESDVertex *vtxSPD = fESD->GetPrimaryVertexSPD(); | |
467 | if (!vtxSPD) | |
468 | return; | |
469 | ||
470 | Int_t nTrax = fESD->GetNumberOfTracks(); | |
471 | if (fVerbosity > 2) | |
472 | AliInfo(Form("%d tracks in event",nTrax)); | |
473 | ||
474 | // Loop 1. | |
475 | Int_t nSelTrax = 0; | |
476 | TObjArray arr(nTrax); | |
477 | arr.SetOwner(1); | |
478 | ||
479 | for (Int_t i = 0; i < nTrax; ++i) { | |
480 | AliESDtrack* esdtrack = fESD->GetTrack(i); | |
481 | if (!esdtrack) { | |
482 | AliError(Form("Couldn't get ESD track %d\n", i)); | |
483 | continue; | |
484 | } | |
8a9d3e12 | 485 | Bool_t trkOK = fEsdTPCOnly->AcceptTrack(esdtrack); |
d688e049 | 486 | if (!trkOK) |
487 | continue; | |
488 | Double_t pt = esdtrack->Pt(); | |
489 | Bool_t ptOK = pt >= fPtMin && pt < fPtMax; | |
490 | if (!ptOK) | |
491 | continue; | |
492 | Double_t eta = esdtrack->Eta(); | |
493 | if (TMath::Abs(eta) > fEtaMax) | |
494 | continue; | |
495 | ||
496 | // create a tpc only track | |
497 | AliESDtrack *newtrack = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack->GetID()); | |
498 | if(!newtrack) | |
499 | continue; | |
500 | if (newtrack->Pt()<=0) { | |
501 | delete newtrack; | |
502 | continue; | |
503 | } | |
8433ff41 | 504 | |
d688e049 | 505 | AliExternalTrackParam exParam; |
506 | Bool_t relate = newtrack->RelateToVertexTPC(vtxSPD,fESD->GetMagneticField(),kVeryBig,&exParam); | |
507 | if (!relate) { | |
508 | delete newtrack; | |
509 | continue; | |
510 | } | |
8433ff41 | 511 | |
d688e049 | 512 | // set the constraint parameters to the track |
513 | newtrack->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance()); | |
8433ff41 | 514 | |
d688e049 | 515 | pt = newtrack->Pt(); |
516 | ptOK = pt >= fPtMin && pt < fPtMax; | |
517 | if (!ptOK) { | |
518 | delete newtrack; | |
519 | continue; | |
520 | } | |
521 | eta = esdtrack->Eta(); | |
522 | if (TMath::Abs(eta) > fEtaMax) { | |
523 | delete newtrack; | |
524 | continue; | |
525 | } | |
526 | arr.Add(newtrack); | |
527 | nSelTrax++; | |
8433ff41 | 528 | } |
beb13b6c | 529 | |
530 | for(Int_t itrack = 0; itrack < nSelTrax; itrack++) { | |
531 | AliVTrack *esdtrack = static_cast<AliESDtrack*>(arr.At(itrack)); | |
532 | if(!esdtrack) { | |
533 | AliError(Form("ERROR: Could not retrieve esdtrack %d",itrack)); | |
534 | continue; | |
535 | } | |
536 | Double_t pt = esdtrack->Pt(); | |
537 | Double_t eta = esdtrack->Eta(); | |
538 | Double_t phi = esdtrack->Phi(); | |
539 | Int_t sign = esdtrack->Charge() > 0 ? 1 : -1; | |
540 | miniEvt->push_back(AliMiniTrack(pt, eta, phi, sign)); | |
541 | } | |
b673a083 | 542 | } else { |
543 | TList *list = InputEvent()->GetList(); | |
544 | TClonesArray *tcaTracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracksName)); | |
f018e7c8 | 545 | |
b673a083 | 546 | if(!tcaTracks){ |
547 | AliError("Ptr to tcaTracks zero"); | |
548 | return; | |
549 | } | |
f018e7c8 | 550 | |
b673a083 | 551 | const Int_t ntracks = tcaTracks->GetEntries(); |
552 | if (miniEvt) | |
553 | miniEvt->reserve(ntracks); | |
554 | else { | |
555 | AliError("Ptr to miniEvt zero"); | |
556 | return; | |
557 | } | |
48a61f36 | 558 | |
b673a083 | 559 | for (Int_t itrack = 0; itrack < ntracks; itrack++) { |
560 | AliVTrack *vtrack = static_cast<AliVTrack*>(tcaTracks->At(itrack)); | |
561 | if (!vtrack) { | |
562 | AliError(Form("ERROR: Could not retrieve track %d",itrack)); | |
563 | continue; | |
564 | } | |
565 | Double_t pt = vtrack->Pt(); | |
566 | Double_t eta = vtrack->Eta(); | |
567 | Double_t phi = vtrack->Phi(); | |
568 | Int_t sign = vtrack->Charge() > 0 ? 1 : -1; | |
569 | miniEvt->push_back(AliMiniTrack(pt, eta, phi, sign)); | |
48a61f36 | 570 | } |
48a61f36 | 571 | } |
beb13b6c | 572 | |
573 | if (fFillMuons) { | |
beb13b6c | 574 | // count good muons |
575 | Int_t nGoodMuons = 0; | |
576 | for (Int_t iMu = 0; iMu<fESD->GetNumberOfMuonTracks(); iMu++) { | |
577 | AliESDMuonTrack* muonTrack = fESD->GetMuonTrack(iMu); | |
b673a083 | 578 | if (muonTrack) { |
beb13b6c | 579 | if (IsGoodMUONtrack(*muonTrack)) nGoodMuons++; |
580 | } | |
581 | } | |
582 | miniEvt->reserve(miniEvt->size()+nGoodMuons); | |
583 | // fill them into the mini event | |
584 | for (Int_t iMu = 0; iMu<fESD->GetNumberOfMuonTracks(); iMu++) { | |
585 | AliESDMuonTrack* muonTrack = fESD->GetMuonTrack(iMu); | |
b673a083 | 586 | if (muonTrack) { |
587 | if (!IsGoodMUONtrack(*muonTrack)) | |
588 | continue; | |
589 | Double_t ptMu = muonTrack->Pt(); | |
590 | Double_t etaMu = muonTrack->Eta(); | |
591 | Double_t phiMu = muonTrack->Phi(); | |
b422cf0d | 592 | Int_t signMu = muonTrack->Charge() > 0 ? 1 : -1; |
beb13b6c | 593 | miniEvt->push_back(AliMiniTrack(ptMu, etaMu, phiMu, signMu)); |
594 | } | |
595 | } | |
596 | } | |
48a61f36 | 597 | } |
598 | ||
599 | //________________________________________________________________________ | |
f6701c6e | 600 | void AliDhcTask::GetAODTracks(MiniEvent* miniEvt) |
48a61f36 | 601 | { |
602 | // Loop twice: 1. Count sel. tracks. 2. Fill vector. | |
603 | ||
b673a083 | 604 | if (fTracksName.IsNull()) { |
605 | Int_t nTrax = fAOD->GetNumberOfTracks(); | |
606 | Int_t nSelTrax = 0; | |
48a61f36 | 607 | |
b673a083 | 608 | if (fVerbosity > 2) |
609 | AliInfo(Form("%d tracks in event",nTrax)); | |
48a61f36 | 610 | |
b673a083 | 611 | // Loop 1. |
612 | for (Int_t i = 0; i < nTrax; ++i) { | |
613 | AliAODTrack* aodtrack = fAOD->GetTrack(i); | |
614 | if (!aodtrack) { | |
615 | AliError(Form("Couldn't get AOD track %d\n", i)); | |
616 | continue; | |
617 | } | |
618 | // See $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C | |
619 | UInt_t tpcOnly = 1 << 7; | |
620 | Bool_t trkOK = aodtrack->TestFilterBit(tpcOnly); | |
621 | if (!trkOK) | |
622 | continue; | |
623 | Double_t pt = aodtrack->Pt(); | |
624 | Bool_t ptOK = pt >= fPtMin && pt < fPtMax; | |
625 | if (!ptOK) | |
626 | continue; | |
627 | Double_t eta = aodtrack->Eta(); | |
628 | if (TMath::Abs(eta) > fEtaMax) | |
629 | continue; | |
630 | nSelTrax++; | |
48a61f36 | 631 | } |
48a61f36 | 632 | |
b673a083 | 633 | if (miniEvt) |
634 | miniEvt->reserve(nSelTrax); | |
635 | else { | |
636 | AliError("!miniEvt"); | |
637 | return; | |
48a61f36 | 638 | } |
b673a083 | 639 | |
640 | // Loop 2. | |
641 | for (Int_t i = 0; i < nTrax; ++i) { | |
642 | AliAODTrack* aodtrack = fAOD->GetTrack(i); | |
643 | if (!aodtrack) { | |
644 | AliError(Form("Couldn't get AOD track %d\n", i)); | |
645 | continue; | |
646 | } | |
48a61f36 | 647 | |
b673a083 | 648 | // See $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C |
649 | UInt_t tpcOnly = 1 << 7; | |
650 | Bool_t trkOK = aodtrack->TestFilterBit(tpcOnly); | |
651 | if (!trkOK) | |
652 | continue; | |
653 | Double_t pt = aodtrack->Pt(); | |
654 | Bool_t ptOK = pt >= fPtMin && pt < fPtMax; | |
655 | if (!ptOK) | |
656 | continue; | |
657 | Double_t eta = aodtrack->Eta(); | |
658 | if (TMath::Abs(eta) > fEtaMax) | |
659 | continue; | |
660 | ||
661 | Double_t phi = aodtrack->Phi(); | |
662 | Int_t sign = aodtrack->Charge() > 0 ? 1 : -1; | |
663 | miniEvt->push_back(AliMiniTrack(pt, eta, phi, sign)); | |
664 | } | |
665 | } else { | |
666 | TList *list = InputEvent()->GetList(); | |
667 | TClonesArray *tcaTracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracksName)); | |
668 | ||
669 | if (!tcaTracks){ | |
670 | AliError("Ptr to tcaTracks zero"); | |
671 | return; | |
672 | } | |
673 | ||
674 | const Int_t ntracks = tcaTracks->GetEntries(); | |
675 | if (miniEvt) | |
676 | miniEvt->reserve(ntracks); | |
677 | else { | |
678 | AliError("Ptr to miniEvt zero"); | |
679 | return; | |
680 | } | |
681 | ||
682 | for (Int_t itrack = 0; itrack < ntracks; itrack++) { | |
683 | AliVTrack *vtrack = static_cast<AliVTrack*>(tcaTracks->At(itrack)); | |
684 | if (!vtrack) { | |
685 | AliError(Form("ERROR: Could not retrieve vtrack %d",itrack)); | |
686 | continue; | |
687 | } | |
688 | Double_t pt = vtrack->Pt(); | |
689 | Double_t eta = vtrack->Eta(); | |
690 | Double_t phi = vtrack->Phi(); | |
691 | Int_t sign = vtrack->Charge() > 0 ? 1 : -1; | |
692 | miniEvt->push_back(AliMiniTrack(pt, eta, phi, sign)); | |
693 | } | |
694 | } | |
48a61f36 | 695 | |
b673a083 | 696 | if (fFillMuons) { |
697 | // count good muons | |
698 | Int_t nGoodMuons = 0; | |
699 | for (Int_t iMu = 0; iMu<fAOD->GetNumberOfTracks(); iMu++) { | |
700 | AliAODTrack* muonTrack = fAOD->GetTrack(iMu); | |
701 | if(muonTrack) { | |
702 | if (IsGoodMUONtrack(*muonTrack)) | |
703 | nGoodMuons++; | |
704 | } | |
705 | } | |
706 | miniEvt->reserve(miniEvt->size()+nGoodMuons); | |
707 | // fill them into the mini event | |
708 | for (Int_t iMu = 0; iMu<fAOD->GetNumberOfTracks(); iMu++) { | |
709 | AliAODTrack* muonTrack = fAOD->GetTrack(iMu); | |
710 | if (muonTrack) { | |
711 | if (!IsGoodMUONtrack(*muonTrack)) | |
712 | continue; | |
713 | Double_t ptMu = muonTrack->Pt(); | |
714 | Double_t etaMu = muonTrack->Eta(); | |
715 | Double_t phiMu = muonTrack->Phi(); | |
b422cf0d | 716 | Int_t signMu = muonTrack->Charge() > 0 ? 1 : -1; |
b673a083 | 717 | miniEvt->push_back(AliMiniTrack(ptMu, etaMu, phiMu, signMu)); |
718 | } | |
719 | } | |
48a61f36 | 720 | } |
48a61f36 | 721 | } |
722 | ||
723 | //________________________________________________________________________ | |
724 | Double_t AliDhcTask::DeltaPhi(Double_t phia, Double_t phib, | |
725 | Double_t rangeMin, Double_t rangeMax) const | |
726 | { | |
727 | Double_t dphi = -999; | |
728 | Double_t pi = TMath::Pi(); | |
729 | ||
730 | if (phia < 0) phia += 2*pi; | |
731 | else if (phia > 2*pi) phia -= 2*pi; | |
732 | if (phib < 0) phib += 2*pi; | |
733 | else if (phib > 2*pi) phib -= 2*pi; | |
734 | dphi = phib - phia; | |
735 | if (dphi < rangeMin) dphi += 2*pi; | |
736 | else if (dphi > rangeMax) dphi -= 2*pi; | |
737 | ||
738 | return dphi; | |
739 | } | |
740 | ||
741 | //________________________________________________________________________ | |
d688e049 | 742 | Int_t AliDhcTask::Correlate(const MiniEvent &evt1, const MiniEvent &evt2, Int_t pairing) |
48a61f36 | 743 | { |
744 | // Triggered angular correlations. If pairing is kSameEvt, particles | |
745 | // within evt1 are correlated. If kDiffEvt, correlate triggers from | |
746 | // evt1 with partners from evt2. | |
8433ff41 | 747 | |
748 | Int_t cbin = fHCent->FindBin(fCentrality); | |
749 | if (fHCent->IsBinOverflow(cbin) || | |
750 | fHCent->IsBinUnderflow(cbin)) | |
751 | return 0; | |
752 | ||
753 | Int_t zbin = fHZvtx->FindBin(fZVertex); | |
754 | if (fHZvtx->IsBinOverflow(zbin) || | |
755 | fHZvtx->IsBinUnderflow(zbin)) | |
756 | return 0; | |
757 | ||
48a61f36 | 758 | Int_t iMax = evt1.size(); |
759 | Int_t jMax = evt2.size(); | |
760 | ||
8433ff41 | 761 | TH2 **hist = fHMs; |
762 | if (pairing == kSameEvt) { | |
763 | hist = fHSs; | |
b673a083 | 764 | fHCent->Fill(fCentrality); |
765 | fHZvtx->Fill(fZVertex); | |
8433ff41 | 766 | } |
767 | ||
768 | Int_t nZvtx = fHZvtx->GetNbinsX(); | |
769 | Int_t nPtTrig = fHPtTrg->GetNbinsX(); | |
770 | Int_t nPtAssc = fHPtAss->GetNbinsX(); | |
771 | ||
772 | Int_t globIndex = (cbin-1)*nZvtx*nPtTrig*nPtAssc+(zbin-1)*nPtTrig*nPtAssc; | |
48a61f36 | 773 | |
8a9d3e12 | 774 | Int_t ptindex = (Int_t)fIndex->Eval(1,1,zbin,cbin); |
775 | ||
beb13b6c | 776 | fHPtTrgEvt->Reset(); |
383b3bca | 777 | for (Int_t i=0; i<iMax; ++i) { |
778 | const AliMiniTrack &a(evt1.at(i)); | |
779 | Float_t pta = a.Pt(); | |
beb13b6c | 780 | fHPtTrgEvt->Fill(pta); |
383b3bca | 781 | if (pairing == kSameEvt) { |
8a9d3e12 | 782 | fHPts[ptindex]->Fill(pta); |
d688e049 | 783 | } |
784 | } | |
73a051c1 | 785 | |
beb13b6c | 786 | Bool_t bCountTrg; // only count the trigger if an associated particle is found |
48a61f36 | 787 | |
beb13b6c | 788 | for (Int_t i=0; i<iMax; ++i) { |
48a61f36 | 789 | // Trigger particles |
8433ff41 | 790 | const AliMiniTrack &a(evt1.at(i)); |
791 | ||
792 | Float_t pta = a.Pt(); | |
beb13b6c | 793 | Float_t etaa = a.Eta(); |
7e990b20 | 794 | Float_t phia = a.Phi(); |
380fb711 | 795 | |
796 | // brief intermezzo in the trigger particle loop: do associated particle QA here in order to avoid double counting | |
797 | if (pairing == kSameEvt) { | |
798 | if (etaa>fEtaALo && etaa<fEtaAHi) { | |
799 | Int_t bbin = fHPtAss->FindBin(pta); | |
800 | if (!(fHPtAss->IsBinOverflow(bbin) || fHPtAss->IsBinUnderflow(bbin))) { | |
801 | fHQAA->Fill(pta,etaa,phia); // fill every associated particle once | |
802 | } | |
803 | } | |
804 | } | |
805 | ||
806 | // back to triggers | |
8433ff41 | 807 | Int_t abin = fHPtTrg->FindBin(pta); |
beb13b6c | 808 | if (fHPtTrg->IsBinOverflow(abin) || fHPtTrg->IsBinUnderflow(abin)) |
809 | continue; | |
810 | ||
811 | if (etaa<fEtaTLo || etaa>fEtaTHi) | |
8433ff41 | 812 | continue; |
48a61f36 | 813 | |
beb13b6c | 814 | // efficiency weighting |
beb13b6c | 815 | Double_t effWtT = 1.0; |
816 | if (fHEffT) { | |
817 | // trigger particle | |
7e990b20 | 818 | const Int_t nEffDimT = fHEffT->GetNdimensions(); |
819 | Int_t effBinT[nEffDimT]; | |
820 | effBinT[0] = fHEffT->GetAxis(0)->FindBin(etaa); | |
821 | effBinT[1] = fHEffT->GetAxis(1)->FindBin(pta); | |
822 | effBinT[2] = fHEffT->GetAxis(2)->FindBin(fCentrality); | |
823 | effBinT[3] = fHEffT->GetAxis(3)->FindBin(fZVertex); | |
824 | if (nEffDimT>4) { | |
825 | effBinT[4] = fHEffT->GetAxis(4)->FindBin(phia); | |
826 | } | |
827 | effWtT = fHEffT->GetBinContent(effBinT); | |
beb13b6c | 828 | } |
829 | ||
48a61f36 | 830 | if (pairing == kSameEvt) { |
7e990b20 | 831 | fHTrk->Fill(phia,etaa); |
380fb711 | 832 | fHQAT->Fill(pta,etaa,phia); |
833 | fHPtTrg->Fill(pta); | |
beb13b6c | 834 | fHPtTrgNorm1S->Fill(pta,fCentrality,fZVertex,effWtT); |
835 | } else { | |
836 | fHPtTrgNorm1M->Fill(pta,fCentrality,fZVertex,effWtT); | |
48a61f36 | 837 | } |
beb13b6c | 838 | bCountTrg = kFALSE; |
48a61f36 | 839 | |
8433ff41 | 840 | for (Int_t j=0; j<jMax; ++j) { |
48a61f36 | 841 | // Associated particles |
842 | if (pairing == kSameEvt && i==j) | |
beb13b6c | 843 | continue; |
48a61f36 | 844 | |
8433ff41 | 845 | const AliMiniTrack &b(evt2.at(j)); |
48a61f36 | 846 | |
48a61f36 | 847 | Float_t ptb = b.Pt(); |
beb13b6c | 848 | Float_t etab = b.Eta(); |
7e990b20 | 849 | Float_t phib = b.Phi(); |
380fb711 | 850 | if (fPtTACrit&&(pta < ptb)) { |
b422cf0d | 851 | continue; |
852 | } | |
48a61f36 | 853 | |
73a051c1 | 854 | Int_t bbin = fHPtAss->FindBin(ptb); |
beb13b6c | 855 | if (fHPtAss->IsBinOverflow(bbin) || fHPtAss->IsBinUnderflow(bbin)) |
856 | continue; | |
857 | ||
858 | if (etab<fEtaALo || etab>fEtaAHi) | |
859 | continue; | |
8433ff41 | 860 | |
7e990b20 | 861 | Float_t dphi = DeltaPhi(phia, phib); |
beb13b6c | 862 | Float_t deta = etaa - etab; |
48a61f36 | 863 | |
8433ff41 | 864 | Int_t index = globIndex+(abin-1)*nPtAssc+(bbin-1); |
8a9d3e12 | 865 | Double_t weight = 1.0; |
beb13b6c | 866 | // number of trigger weight event-by-event (CMS method) |
867 | if (fDoWeights) { | |
868 | Double_t nTrgs = fHPtTrgEvt->GetBinContent(abin); | |
869 | if (nTrgs==0.0) { | |
870 | AliError(Form("Trying to work with number of triggers weight = %f",nTrgs)); | |
871 | return 0; | |
872 | } | |
873 | weight *= 1./nTrgs; | |
874 | } | |
8a9d3e12 | 875 | |
beb13b6c | 876 | // efficiency weighting |
877 | if (fHEffT) { | |
878 | // trigger particle | |
879 | weight *= effWtT; | |
73a051c1 | 880 | } |
beb13b6c | 881 | if (fHEffA) { |
882 | // associated particle | |
7e990b20 | 883 | const Int_t nEffDimA = fHEffA->GetNdimensions(); |
884 | Int_t effBinA[nEffDimA]; | |
885 | effBinA[0] = fHEffA->GetAxis(0)->FindBin(etab); | |
886 | effBinA[1] = fHEffA->GetAxis(1)->FindBin(ptb); | |
887 | effBinA[2] = fHEffA->GetAxis(2)->FindBin(fCentrality); | |
888 | effBinA[3] = fHEffA->GetAxis(3)->FindBin(fZVertex); | |
889 | if (nEffDimA>4) { | |
890 | effBinA[4] = fHEffA->GetAxis(4)->FindBin(phib); | |
891 | } | |
892 | weight *= fHEffA->GetBinContent(effBinA); | |
beb13b6c | 893 | } |
380fb711 | 894 | if (hist[index]) { // check if this histogram exists, relevant in the fPtTACrit==kFALSE case |
895 | hist[index]->Fill(dphi,deta,weight); | |
896 | bCountTrg = kTRUE; | |
897 | if (pairing == kSameEvt) { | |
898 | fHPtAss->Fill(ptb); // fill every associated particle every time it is used | |
899 | } | |
beb13b6c | 900 | } |
901 | } | |
902 | if (bCountTrg) { | |
903 | if (pairing == kSameEvt) { | |
904 | fHPtTrgNorm2S->Fill(pta,fCentrality,fZVertex,effWtT); | |
905 | } else { | |
906 | fHPtTrgNorm2M->Fill(pta,fCentrality,fZVertex,effWtT); | |
8433ff41 | 907 | } |
48a61f36 | 908 | } |
909 | } | |
910 | ||
beb13b6c | 911 | return 1; |
48a61f36 | 912 | } |
913 | ||
914 | //________________________________________________________________________ | |
915 | void AliDhcTask::Terminate(Option_t *) | |
916 | { | |
917 | // Draw result to the screen | |
918 | // Called once at the end of the query | |
48a61f36 | 919 | } |
920 | ||
921 | //________________________________________________________________________ | |
922 | Bool_t AliDhcTask::VertexOk(TObject* obj) const | |
923 | { | |
924 | // Modified from AliAnalyseLeadingTrackUE::VertexSelection() | |
925 | ||
926 | Int_t nContributors = 0; | |
927 | Double_t zVertex = 999; | |
928 | TString name(""); | |
929 | ||
930 | if (obj->InheritsFrom("AliESDEvent")) { | |
931 | AliESDEvent* esdevt = (AliESDEvent*) obj; | |
932 | const AliESDVertex* vtx = esdevt->GetPrimaryVertex(); | |
933 | if (!vtx) | |
934 | return 0; | |
935 | nContributors = vtx->GetNContributors(); | |
936 | zVertex = vtx->GetZ(); | |
937 | name = vtx->GetName(); | |
938 | } | |
939 | else if (obj->InheritsFrom("AliAODEvent")) { | |
940 | AliAODEvent* aodevt = (AliAODEvent*) obj; | |
941 | if (aodevt->GetNumberOfVertices() < 1) | |
942 | return 0; | |
943 | const AliAODVertex* vtx = aodevt->GetPrimaryVertex(); | |
944 | nContributors = vtx->GetNContributors(); | |
945 | zVertex = vtx->GetZ(); | |
946 | name = vtx->GetName(); | |
947 | } | |
948 | ||
949 | // Reject if TPC-only vertex | |
8433ff41 | 950 | if (name.CompareTo("TPCVertex")==0) |
48a61f36 | 951 | return kFALSE; |
952 | ||
953 | // Check # contributors and range... | |
954 | if( nContributors < 1 || TMath::Abs(zVertex) > fZVtxMax ) { | |
955 | return kFALSE; | |
956 | } | |
957 | ||
958 | return kTRUE; | |
959 | } | |
beb13b6c | 960 | |
961 | //________________________________________________________________________ | |
962 | Bool_t AliDhcTask::IsGoodMUONtrack(AliESDMuonTrack &track) | |
963 | { | |
b673a083 | 964 | // Applying track cuts for MUON tracks |
965 | ||
966 | if (!track.ContainTrackerData()) | |
967 | return kFALSE; | |
968 | if (!track.ContainTriggerData()) | |
969 | return kFALSE; | |
970 | Double_t thetaTrackAbsEnd = TMath::ATan(track.GetRAtAbsorberEnd()/505.) * TMath::RadToDeg(); | |
971 | if ((thetaTrackAbsEnd < 2.) || (thetaTrackAbsEnd > 10.)) | |
972 | return kFALSE; | |
973 | Double_t eta = track.Eta(); | |
974 | if ((eta < -4.) || (eta > -2.5)) | |
975 | return kFALSE; | |
976 | if (track.GetMatchTrigger() < 0.5) | |
977 | return kFALSE; | |
978 | return kTRUE; | |
979 | } | |
980 | ||
981 | //________________________________________________________________________ | |
982 | Bool_t AliDhcTask::IsGoodMUONtrack(AliAODTrack &track) | |
983 | { | |
984 | // Applying track cuts for MUON tracks | |
985 | ||
986 | if (!track.IsMuonTrack()) | |
987 | return kFALSE; | |
988 | Double_t dThetaAbs = TMath::ATan(track.GetRAtAbsorberEnd()/505.)* TMath::RadToDeg(); | |
989 | if ((dThetaAbs<2.) || (dThetaAbs>10.)) | |
990 | return kFALSE; | |
991 | Double_t dEta = track.Eta(); | |
992 | if ((dEta<-4.) || (dEta>2.5)) | |
993 | return kFALSE; | |
994 | if (track.GetMatchTrigger()<0.5) | |
995 | return kFALSE; | |
996 | return kTRUE; | |
997 | } | |
998 | ||
999 | //________________________________________________________________________ | |
1000 | AliDhcTask::~AliDhcTask() | |
1001 | { | |
1002 | //Destructor | |
1003 | if (fPoolMgr) | |
1004 | delete fPoolMgr; | |
beb13b6c | 1005 | } |