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