]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGCF/Correlations/DPhi/FourierDecomposition/AliDhcTask.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Fixed random number initialisation (SetMRPY(1,) does not have effect for AliPythia6)
[u/mrichter/AliRoot.git] / PWGCF / Correlations / DPhi / FourierDecomposition / AliDhcTask.cxx
CommitLineData
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"
8433ff41 10#include "TProfile2D.h"
48a61f36 11#include "TROOT.h"
12#include "TTree.h"
48a61f36 13#include "AliAODEvent.h"
14#include "AliAODInputHandler.h"
15#include "AliAnalysisManager.h"
16#include "AliAnalysisTask.h"
17#include "AliCentrality.h"
18#include "AliDhcTask.h"
19#include "AliESDEvent.h"
20#include "AliESDInputHandler.h"
21#include "AliESDtrackCuts.h"
8433ff41 22#include "AliPool.h"
48a61f36 23#include "AliVParticle.h"
24
25ClassImp(AliDhcTask)
26
27//________________________________________________________________________
28AliDhcTask::AliDhcTask(const char *name)
42036329 29: AliAnalysisTaskSE(name), fVerbosity(0), fEtaMax(1), fZVtxMax(10), fPtMin(0.25), fPtMax(15),
30 fESD(0), fAOD(0), fOutputList(0), fHistPt(0), fHEvt(0), fHTrk(0), fHPtAss(0), fHPtTrg(0),
31 fHCent(0), fHZvtx(0), fNbins(0), fHSs(0), fHMs(0),
32 fIndex(0), fMeanPtTrg(0), fMeanPtAss(0), fMean2PtTrg(0), fMean2PtAss(0),
33 fCentrality(99), fZVertex(99), fEsdTrackCutsTPCOnly(0), fPoolMgr(0)
48a61f36 34{
35 // Constructor
36
37 // Define input and output slots here
38 // Input slot #0 works with a TChain
39 DefineInput(0, TChain::Class());
40 // Output slot #0 id reserved by the base class for AOD
41 // Output slot #1 writes into a TH1 container
42 DefineOutput(1, TList::Class());
43
48a61f36 44 fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.,SPDVertex.,TPCVertex.,Tracks "
45 "AOD:header,tracks,vertices,";
46}
47
48//________________________________________________________________________
49void AliDhcTask::UserCreateOutputObjects()
50{
51 // Create histograms
52 // Called once (per slave on PROOF!)
53
48a61f36 54 fOutputList = new TList();
55 fOutputList->SetOwner(1);
56
8433ff41 57 fEsdTrackCutsTPCOnly = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
58 //fEsdTrackCutsTPCOnly->SetMinNClustersTPC(70);
59 fEsdTrackCutsTPCOnly->SetMinNCrossedRowsTPC(70);
60 fEsdTrackCutsTPCOnly->SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
48a61f36 61
8433ff41 62 BookHistos();
63 InitEventMixer();
48a61f36 64 PostData(1, fOutputList);
65}
66
67//________________________________________________________________________
68void AliDhcTask::BookHistos()
69{
8433ff41 70 // Book histograms.
71
72 Int_t nDeta=20, nPtAssc=12, nPtTrig=12, nCent=12, nDphi=36, nZvtx=8;
48a61f36 73
8433ff41 74 Double_t ptt[] = {0.25, 0.5, 0.75, 1.0, 1.5, 2.0, 2.5, 3.0, 4.0, 5.0, 6.0, 8.0, 15};
75 Double_t pta[] = {0.25, 0.5, 0.75, 1.0, 1.5, 2.0, 2.5, 3.0, 4.0, 5.0, 6.0, 8.0, 15};
76 Double_t cent[] = {0, 1, 2, 3, 4, 5, 10, 20, 30, 40, 50, 60, 90};
77 Double_t zvtx[] = {-10, -6, -4, -2, 0, 2, 4, 6, 10};
48a61f36 78
79 // Event histo
8433ff41 80 fHEvt = new TH2F("fHEvt", "Event-level variables; Zvtx; Cent", 30, -15, 15, 101, 0, 101);
81 fOutputList->Add(fHEvt);
48a61f36 82 // Track histo
83 fHTrk = new TH2F("fHTrk", "Track-level variables",
8433ff41 84 100, 0, TMath::TwoPi(), 100, -fEtaMax, +fEtaMax);
85 fOutputList->Add(fHTrk);
86
48a61f36 87 // Left over from the tutorial :)
8433ff41 88 fHistPt = new TH1F("fHistPt", "P_{T} distribution", 200, 0., fPtMax);
48a61f36 89 fHistPt->GetXaxis()->SetTitle("P_{T} (GeV/c)");
90 fHistPt->GetYaxis()->SetTitle("dN/dP_{T} (c/GeV)");
91 fHistPt->SetMarkerStyle(kFullCircle);
8433ff41 92 fOutputList->Add(fHistPt);
93
94 fHPtAss = new TH1F("fHPtAss","PtAssoc;P_{T} (GeV/c) [GeV/c]",nPtAssc,pta);
95 fOutputList->Add(fHPtAss);
96 fHPtTrg = new TH1F("fHPtTrg","PtTrig;P_{T} (GeV/c) [GeV/c]",nPtTrig,ptt);
97 fOutputList->Add(fHPtTrg);
98 fHCent = new TH1F("fHCent","Cent;bins",nCent,cent);
99 fOutputList->Add(fHCent);
100 fHZvtx = new TH1F("fHZvtx","Zvertex;bins",nZvtx,zvtx);
101 fOutputList->Add(fHZvtx);
102
103 fNbins = nPtTrig*nPtAssc*nCent*nZvtx;
104 fHSs = new TH2*[fNbins];
105 fHMs = new TH2*[fNbins];
106
107 fIndex = new TFormula("GlobIndex","(t-1)*[0]*[1]*[2]+(z-1)*[0]*[1]+(x-1)*[0]+(y-1)+0*[4]");
108 fIndex->SetParameters(nPtTrig,nPtAssc,nZvtx,nCent);
109 fIndex->SetParNames("NTrigBins","NAssocBins", "NZvertexBins", "NCentBins");
110 //fOutputList->Add(fIndex);
111
112 fMeanPtTrg = new TProfile2D*[nPtTrig*nPtAssc];
113 fMeanPtAss = new TProfile2D*[nPtTrig*nPtAssc];
114 fMean2PtTrg = new TProfile2D*[nPtTrig*nPtAssc];
115 fMean2PtAss = new TProfile2D*[nPtTrig*nPtAssc];
116 for (Int_t c=1; c<=nCent; ++c) {
117 TString title(Form("cen=%d (%.1f)",c,fHCent->GetBinCenter(c)));
118 fMeanPtTrg[c-1] = new TProfile2D(Form("hMeanPtTrgCen%d",c),title,nPtTrig,ptt,nPtAssc,pta);
119 fMeanPtAss[c-1] = new TProfile2D(Form("hMeanPtAssCen%d",c),title,nPtTrig,ptt,nPtAssc,pta);
120 fMean2PtTrg[c-1] = new TProfile2D(Form("hMean2PtTrgCen%d",c),title,nPtTrig,ptt,nPtAssc,pta);
121 fMean2PtAss[c-1] = new TProfile2D(Form("hMean2PtAssCen%d",c),title,nPtTrig,ptt,nPtAssc,pta);
122 fOutputList->Add(fMeanPtTrg[c-1]);
123 fOutputList->Add(fMeanPtAss[c-1]);
124 fOutputList->Add(fMean2PtTrg[c-1]);
125 fOutputList->Add(fMean2PtAss[c-1]);
126 }
127
128 Int_t count = 0;
129 for (Int_t c=1; c<=nCent; ++c) {
130 for (Int_t z=1; z<=nZvtx; ++z) {
131 for (Int_t t=1; t<=nPtTrig; ++t) {
132 for (Int_t a=1; a<=nPtAssc; ++a) {
133 fHSs[count] = 0;
134 fHMs[count] = 0;
135 if (a>t) {
136 ++count;
137 continue;
138 }
139 TString title(Form("cen=%d (%.1f), zVtx=%d (%.1f), trig=%d (%.1f), assc=%d (%.1f)",
140 c,fHCent->GetBinCenter(c), z,fHZvtx->GetBinCenter(z),
141 t,fHPtTrg->GetBinCenter(t),a, fHPtAss->GetBinCenter(a)));
142 fHSs[count] = new TH2F(Form("hS%d",count), Form("Signal %s",title.Data()),
143 nDphi,-0.5*TMath::Pi(),1.5*TMath::Pi(),nDeta,-2*fEtaMax,2*fEtaMax);
144 fHMs[count] = new TH2F(Form("hM%d",count), Form("Signal %s",title.Data()),
145 nDphi,-0.5*TMath::Pi(),1.5*TMath::Pi(),nDeta,-2*fEtaMax,2*fEtaMax);
146 fOutputList->Add(fHSs[count]);
147 fOutputList->Add(fHMs[count]);
148 if (fVerbosity>5)
149 cout << count << " " << fIndex->Eval(t,a,z,c) << ": " << title << endl;
150 ++count;
151 }
152 }
153 }
154 }
155
48a61f36 156 return;
157}
158
159//________________________________________________________________________
160void AliDhcTask::InitEventMixer()
161{
162 // The effective pool size in events is set by trackDepth, so more
163 // low-mult events are required to maintain the threshold than
164 // high-mult events. Centrality pools are indep. of data histogram
165 // binning, no need to match.
166
167 Int_t trackDepth = 1000; // # tracks to fill pool
168 Int_t poolsize = 200; // Maximum number of events
169
170 // Centrality pools
8433ff41 171 Int_t nCentBins = 12;
172 Double_t centBins[] = {0,1,2,3,4,5,10,20,30,40,50,60,90.1};
48a61f36 173
174 //Int_t nCentBins = 1;
175 //Double_t centBins[] = {-1,100.1};
176
177 // Z-vertex pools
178 Int_t nZvtxBins = 8;
8433ff41 179 Double_t zvtxbin[] = {-10,-6,-4,-2,0,2,4,6,10};
180
181 fPoolMgr = new AliEvtPoolManager();
182 fPoolMgr->SetTargetTrackDepth(trackDepth);
183 if (fVerbosity>4)
184 fPoolMgr->SetDebug(1);
185 fPoolMgr->InitEventPools(poolsize, nCentBins, centBins, nZvtxBins, zvtxbin);
48a61f36 186
48a61f36 187 return;
188}
189
190//________________________________________________________________________
191void AliDhcTask::UserExec(Option_t *)
192{
193 // Main loop, called for each event.
194
195 static int iEvent = -1; ++iEvent;
196
197 if (fVerbosity>2) {
198 if (iEvent % 10 == 0)
199 cout << iEvent << endl;
200 }
201
202 Int_t dType = -1; // Will be set to kESD or kAOD.
8433ff41 203 MiniEvent* sTracks = 0; // Vector of selected AliMiniTracks.
48a61f36 204 Double_t centCL1 = -1;
205
206 LoadBranches();
207
208 // Get event pointers, check for signs of life
209 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
210 fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
211 if (fESD)
212 dType = kESD;
213 else if (fAOD)
214 dType = kAOD;
215 else {
216 AliError("Neither fESD nor fAOD available");
217 return;
218 }
219
220 // Centrality, vertex, other event variables...
221 if (dType == kESD) {
222 if (!VertexOk(fESD)) {
223 if (fVerbosity > 1)
224 AliInfo(Form("Event REJECTED. z = %.1f", fZVertex));
225 return;
226 }
227 const AliESDVertex* vertex = fESD->GetPrimaryVertex();
228 fZVertex = vertex->GetZ();
229 if(fESD->GetCentrality()) {
230 fCentrality =
231 fESD->GetCentrality()->GetCentralityPercentile("V0M");
232 centCL1 =
233 fESD->GetCentrality()->GetCentralityPercentile("CL1");
234 }
235 }
236 if (dType == kAOD) {
237 const AliAODVertex* vertex = fAOD->GetPrimaryVertex();
238 fZVertex = vertex->GetZ();
239 if (!VertexOk(fAOD)) {
240 if (fVerbosity > 1)
241 Info("Exec()", "Event REJECTED. z = %.1f", fZVertex);
242 return;
243 }
244 const AliCentrality *aodCent = fAOD->GetHeader()->GetCentralityP();
245 if (aodCent) {
246 fCentrality = aodCent->GetCentralityPercentile("V0M");
247 centCL1 = aodCent->GetCentralityPercentile("CL1");
248 }
249 }
250
251 // Fill Event histogram
252 fHEvt->Fill(fZVertex, fCentrality);
253
254 if (fCentrality > 90. || fCentrality < 0) {
255 AliInfo(Form("Event REJECTED. fCentrality = %.1f", fCentrality));
256 return;
257 }
258
259 // Get array of selected tracks
260 if (dType == kESD) {
261 sTracks = GetESDTrax();
262 }
263 if (dType == kAOD) {
264 sTracks = GetAODTrax();
265 }
266
267 // Get pool containing tracks from other events like this one
8433ff41 268 AliEvtPool* pool = fPoolMgr->GetEventPool(fCentrality, fZVertex);
48a61f36 269 if (!pool) {
c71ff6f4 270 AliWarning(Form("No pool found. Centrality %f, ZVertex %f",
271 fCentrality, fZVertex));
272 sTracks->clear();
48a61f36 273 return;
274 }
275
276 if (!pool->IsReady()) {
277 pool->UpdatePool(sTracks);
278 return;
279 }
280
281 if (pool->IsFirstReady()) {
282 // recover events missed before the pool is ready
283 Int_t nEvs = pool->GetCurrentNEvents();
284 if (nEvs>1) {
285 for (Int_t i=0; i<nEvs; ++i) {
286 MiniEvent* evI = pool->GetEvent(i);
287 for (Int_t j=0; j<nEvs; ++j) {
288 MiniEvent* evJ = pool->GetEvent(j);
289 if (i==j) {
290 Correlate(*evI, *evJ, kSameEvt);
291 } else {
292 Correlate(*evI, *evJ, kDiffEvt, 1.0);
293 }
294 }
295 }
296 }
297 } else { /* standard case: same event, then mix*/
298 Correlate(*sTracks, *sTracks, kSameEvt);
299 Int_t nMix = pool->GetCurrentNEvents();
300 for (Int_t jMix=0; jMix<nMix; ++jMix) {
301 MiniEvent* bgTracks = pool->GetEvent(jMix);
302 Correlate(*sTracks, *bgTracks, kDiffEvt, 1.0);
303 }
304 }
305
48a61f36 306 pool->UpdatePool(sTracks);
307 PostData(1, fOutputList);
308 return;
309}
310
311//________________________________________________________________________
312MiniEvent* AliDhcTask::GetESDTrax() const
313{
314 // Loop twice: 1. Count sel. tracks. 2. Fill vector.
315
8433ff41 316 const AliESDVertex *vtxSPD = fESD->GetPrimaryVertexSPD();
317 if (!vtxSPD)
318 return 0;
48a61f36 319
8433ff41 320 Int_t nTrax = fESD->GetNumberOfTracks();
48a61f36 321 if (fVerbosity > 2)
322 AliInfo(Form("%d tracks in event",nTrax));
323
324 // Loop 1.
8433ff41 325 Int_t nSelTrax = 0;
326 TObjArray arr(nTrax);
327 arr.SetOwner(1);
328
48a61f36 329 for (Int_t i = 0; i < nTrax; ++i) {
330 AliESDtrack* esdtrack = fESD->GetTrack(i);
331 if (!esdtrack) {
332 AliError(Form("Couldn't get ESD track %d\n", i));
333 continue;
334 }
335 Bool_t trkOK = fEsdTrackCutsTPCOnly->AcceptTrack(esdtrack);
336 if (!trkOK)
337 continue;
338 Double_t pt = esdtrack->Pt();
339 Bool_t ptOK = pt >= fPtMin && pt < fPtMax;
340 if (!ptOK)
341 continue;
342 Double_t eta = esdtrack->Eta();
8433ff41 343 if (TMath::Abs(eta) > fEtaMax)
48a61f36 344 continue;
8433ff41 345
346 // create a tpc only track
347 AliESDtrack *newtrack = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack->GetID());
348 if(!newtrack)
349 continue;
350 if (newtrack->Pt()<=0) {
351 delete newtrack;
352 continue;
353 }
354
355 AliExternalTrackParam exParam;
356 Bool_t relate = newtrack->RelateToVertexTPC(vtxSPD,fESD->GetMagneticField(),kVeryBig,&exParam);
357 if (!relate) {
358 delete newtrack;
359 continue;
360 }
361
362 // set the constraint parameters to the track
363 newtrack->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
364
365 pt = newtrack->Pt();
366 ptOK = pt >= fPtMin && pt < fPtMax;
367 if (!ptOK) {
368 delete newtrack;
369 continue;
370 }
371 eta = esdtrack->Eta();
372 if (TMath::Abs(eta) > fEtaMax) {
373 delete newtrack;
374 continue;
375 }
376 arr.Add(newtrack);
48a61f36 377 nSelTrax++;
378 }
379
380 MiniEvent* miniEvt = new MiniEvent(0);
381 miniEvt->reserve(nSelTrax);
382
8433ff41 383 // Loop 2.
384 for (Int_t i = 0; i < nSelTrax; ++i) {
385 AliESDtrack* esdtrack = static_cast<AliESDtrack*>(arr.At(i));
48a61f36 386 if (!esdtrack) {
387 AliError(Form("Couldn't get ESD track %d\n", i));
388 continue;
389 }
48a61f36 390 Double_t pt = esdtrack->Pt();
48a61f36 391 Double_t eta = esdtrack->Eta();
48a61f36 392 Double_t phi = esdtrack->Phi();
393 Int_t sign = esdtrack->Charge() > 0 ? 1 : -1;
8433ff41 394 miniEvt->push_back(AliMiniTrack(pt, eta, phi, sign));
48a61f36 395 }
396 return miniEvt;
397}
398
399//________________________________________________________________________
400MiniEvent* AliDhcTask::GetAODTrax() const
401{
402 // Loop twice: 1. Count sel. tracks. 2. Fill vector.
403
404 Int_t nTrax = fAOD->GetNumberOfTracks();
405 Int_t nSelTrax = 0;
406
407 if (fVerbosity > 2)
408 AliInfo(Form("%d tracks in event",nTrax));
409
410 // Loop 1.
411 for (Int_t i = 0; i < nTrax; ++i) {
412 AliAODTrack* aodtrack = fAOD->GetTrack(i);
413 if (!aodtrack) {
414 AliError(Form("Couldn't get AOD track %d\n", i));
415 continue;
416 }
417 // See $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C
418 UInt_t tpcOnly = 1 << 7;
419 Bool_t trkOK = aodtrack->TestFilterBit(tpcOnly);
420 if (!trkOK)
421 continue;
422 Double_t pt = aodtrack->Pt();
423 Bool_t ptOK = pt >= fPtMin && pt < fPtMax;
424 if (!ptOK)
425 continue;
426 Double_t eta = aodtrack->Eta();
8433ff41 427 if (TMath::Abs(eta) > fEtaMax)
48a61f36 428 continue;
429 nSelTrax++;
430 }
431
432 MiniEvent* miniEvt = new MiniEvent(0);
433 miniEvt->reserve(nSelTrax);
434
435 // Loop 2.
436 for (Int_t i = 0; i < nTrax; ++i) {
437 AliAODTrack* aodtrack = fAOD->GetTrack(i);
438 if (!aodtrack) {
439 AliError(Form("Couldn't get AOD track %d\n", i));
440 continue;
441 }
442
443 // See $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C
444 UInt_t tpcOnly = 1 << 7;
445 Bool_t trkOK = aodtrack->TestFilterBit(tpcOnly);
446 if (!trkOK)
447 continue;
448 Double_t pt = aodtrack->Pt();
449 Bool_t ptOK = pt >= fPtMin && pt < fPtMax;
450 if (!ptOK)
451 continue;
452 Double_t eta = aodtrack->Eta();
8433ff41 453 if (TMath::Abs(eta) > fEtaMax)
48a61f36 454 continue;
455
456 Double_t phi = aodtrack->Phi();
457 Int_t sign = aodtrack->Charge() > 0 ? 1 : -1;
8433ff41 458 miniEvt->push_back(AliMiniTrack(pt, eta, phi, sign));
48a61f36 459 }
460 return miniEvt;
461}
462
463//________________________________________________________________________
464Double_t AliDhcTask::DeltaPhi(Double_t phia, Double_t phib,
465 Double_t rangeMin, Double_t rangeMax) const
466{
467 Double_t dphi = -999;
468 Double_t pi = TMath::Pi();
469
470 if (phia < 0) phia += 2*pi;
471 else if (phia > 2*pi) phia -= 2*pi;
472 if (phib < 0) phib += 2*pi;
473 else if (phib > 2*pi) phib -= 2*pi;
474 dphi = phib - phia;
475 if (dphi < rangeMin) dphi += 2*pi;
476 else if (dphi > rangeMax) dphi -= 2*pi;
477
478 return dphi;
479}
480
481//________________________________________________________________________
482Int_t AliDhcTask::Correlate(const MiniEvent &evt1, const MiniEvent &evt2,
8433ff41 483 Int_t pairing, Double_t /*weight*/)
48a61f36 484{
485 // Triggered angular correlations. If pairing is kSameEvt, particles
486 // within evt1 are correlated. If kDiffEvt, correlate triggers from
487 // evt1 with partners from evt2.
8433ff41 488
489 Int_t cbin = fHCent->FindBin(fCentrality);
490 if (fHCent->IsBinOverflow(cbin) ||
491 fHCent->IsBinUnderflow(cbin))
492 return 0;
493
494 Int_t zbin = fHZvtx->FindBin(fZVertex);
495 if (fHZvtx->IsBinOverflow(zbin) ||
496 fHZvtx->IsBinUnderflow(zbin))
497 return 0;
498
48a61f36 499 Int_t iMax = evt1.size();
500 Int_t jMax = evt2.size();
501
8433ff41 502 TH2 **hist = fHMs;
503 if (pairing == kSameEvt) {
504 hist = fHSs;
505 fHCent->AddBinContent(cbin);
506 fHZvtx->AddBinContent(zbin);
507 }
508
509 Int_t nZvtx = fHZvtx->GetNbinsX();
510 Int_t nPtTrig = fHPtTrg->GetNbinsX();
511 Int_t nPtAssc = fHPtAss->GetNbinsX();
512
513 Int_t globIndex = (cbin-1)*nZvtx*nPtTrig*nPtAssc+(zbin-1)*nPtTrig*nPtAssc;
48a61f36 514
515 for (Int_t i=0; i<iMax; ++i) {
516
517 // Trigger particles
8433ff41 518 const AliMiniTrack &a(evt1.at(i));
519
520 Float_t pta = a.Pt();
521 Int_t abin = fHPtTrg->FindBin(pta);
522 if (fHPtTrg->IsBinOverflow(abin) ||
523 fHPtTrg->IsBinUnderflow(abin))
524 continue;
48a61f36 525
526 if (pairing == kSameEvt) {
8433ff41 527 fHistPt->Fill(pta);
48a61f36 528 fHTrk->Fill(a.Phi(),a.Eta());
8433ff41 529 fHPtTrg->AddBinContent(abin);
48a61f36 530 }
531
8433ff41 532 for (Int_t j=0; j<jMax; ++j) {
48a61f36 533 // Associated particles
534 if (pairing == kSameEvt && i==j)
535 continue;
536
8433ff41 537 const AliMiniTrack &b(evt2.at(j));
48a61f36 538
48a61f36 539 Float_t ptb = b.Pt();
540 if (pta < ptb)
541 continue;
542
8433ff41 543 Int_t bbin = fHPtTrg->FindBin(ptb);
544 if (fHPtAss->IsBinOverflow(bbin) ||
545 fHPtAss->IsBinUnderflow(bbin))
546 continue;
547
48a61f36 548 Float_t dphi = DeltaPhi(a.Phi(), b.Phi());
549 Float_t deta = a.Eta() - b.Eta();
550
8433ff41 551 Int_t index = globIndex+(abin-1)*nPtAssc+(bbin-1);
552 hist[index]->Fill(dphi,deta);
553
554 if (pairing == kSameEvt) {
555 fHPtAss->AddBinContent(bbin);
556 fMeanPtTrg[cbin-1]->Fill(pta,ptb,pta);
557 fMeanPtAss[cbin-1]->Fill(pta,ptb,ptb);
558 fMean2PtTrg[cbin-1]->Fill(pta,ptb,pta*pta);
559 fMean2PtAss[cbin-1]->Fill(pta,ptb,ptb*ptb);
560 }
48a61f36 561 }
562 }
563
564 return 0;
565}
566
567//________________________________________________________________________
568void AliDhcTask::Terminate(Option_t *)
569{
570 // Draw result to the screen
571 // Called once at the end of the query
572
8433ff41 573 delete fPoolMgr;
574
575 fHCent->SetEntries(fHCent->Integral());
576 fHZvtx->SetEntries(fHZvtx->Integral());
577 fHPtTrg->SetEntries(fHPtTrg->Integral());
578 fHPtAss->SetEntries(fHPtAss->Integral());
579
48a61f36 580 if (gROOT->IsBatch())
581 return;
582
583 fOutputList = dynamic_cast<TList*> (GetOutputData(1));
584 if (!fOutputList) {
585 AliError("Output list not available");
586 return;
587 }
588
c71ff6f4 589 fHistPt = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistPt"));
48a61f36 590 if (!fHistPt) {
591 AliError("ERROR: fHistPt not available\n");
592 return;
593 }
48a61f36 594 TCanvas *c1 = new TCanvas("AliDhcTask","Pt",10,10,510,510);
595 c1->cd(1)->SetLogy();
596 fHistPt->DrawCopy("E");
597}
598
599//________________________________________________________________________
600Bool_t AliDhcTask::VertexOk(TObject* obj) const
601{
602 // Modified from AliAnalyseLeadingTrackUE::VertexSelection()
603
604 Int_t nContributors = 0;
605 Double_t zVertex = 999;
606 TString name("");
607
608 if (obj->InheritsFrom("AliESDEvent")) {
609 AliESDEvent* esdevt = (AliESDEvent*) obj;
610 const AliESDVertex* vtx = esdevt->GetPrimaryVertex();
611 if (!vtx)
612 return 0;
613 nContributors = vtx->GetNContributors();
614 zVertex = vtx->GetZ();
615 name = vtx->GetName();
616 }
617 else if (obj->InheritsFrom("AliAODEvent")) {
618 AliAODEvent* aodevt = (AliAODEvent*) obj;
619 if (aodevt->GetNumberOfVertices() < 1)
620 return 0;
621 const AliAODVertex* vtx = aodevt->GetPrimaryVertex();
622 nContributors = vtx->GetNContributors();
623 zVertex = vtx->GetZ();
624 name = vtx->GetName();
625 }
626
627 // Reject if TPC-only vertex
8433ff41 628 if (name.CompareTo("TPCVertex")==0)
48a61f36 629 return kFALSE;
630
631 // Check # contributors and range...
632 if( nContributors < 1 || TMath::Abs(zVertex) > fZVtxMax ) {
633 return kFALSE;
634 }
635
636 return kTRUE;
637}