]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/vertexingHF/charmFlow/AliAnalysisTaskSEHFv2.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
possibility to choose a combination of 3 subevents out of TPC,TPC+,TPC-,V0A,V0C for...
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / charmFlow / AliAnalysisTaskSEHFv2.cxx
CommitLineData
a8f6c03f 1/**************************************************************************
2 * Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/////////////////////////////////////////////////////////////
17//
18// AliAnalysisTaskSEHFv2 gives the needed tools for the D
19// mesons v2 analysis with event plane method
20// Authors: Chiara Bianchin, cbianchi@pd.infn.it,
21// Robert Grajcarek, grajcarek@physi.uni-heidelberg.de
22// Giacomo Ortona, ortona@to.infn.it
23// Carlos Eugenio Perez Lara, carlos.eugenio.perez.lara@cern.ch
24//
25/////////////////////////////////////////////////////////////
26
f829c4f7 27/* $Id$ */
28
a8f6c03f 29#include <Riostream.h>
30#include <TClonesArray.h>
31#include <TCanvas.h>
32#include <TList.h>
33#include <TFile.h>
34#include <TString.h>
35#include <TH1F.h>
36#include <TH2F.h>
37#include <TGraphErrors.h>
38#include <TGraph.h>
39#include <TDatabasePDG.h>
40#include <TRandom3.h>
41#include <TVector2.h>
42#include <TArrayF.h>
43
44#include <AliLog.h>
45#include <AliAnalysisDataSlot.h>
46#include <AliAnalysisDataContainer.h>
47#include "AliAnalysisManager.h"
48#include "AliAODHandler.h"
49#include "AliAODEvent.h"
50#include "AliAODVertex.h"
51#include "AliAODTrack.h"
52#include "AliAODMCHeader.h"
53#include "AliAODMCParticle.h"
54#include "AliAODRecoDecayHF3Prong.h"
55#include "AliAODRecoDecayHF.h"
56#include "AliAODRecoDecayHF2Prong.h"
57#include "AliAODRecoDecayHF4Prong.h"
58#include "AliAODRecoCascadeHF.h"
a8f6c03f 59
60#include "AliAnalysisTaskSE.h"
61#include "AliRDHFCutsDplustoKpipi.h"
62#include "AliRDHFCutsD0toKpipipi.h"
63#include "AliRDHFCutsDstoKKpi.h"
64#include "AliRDHFCutsDStartoKpipi.h"
65#include "AliRDHFCutsD0toKpi.h"
66#include "AliRDHFCutsLctopKpi.h"
67
68#include "AliHFMassFitter.h"
69#include "AliEventplane.h"
70#include "AliFlowTrack.h"
71#include "AliFlowVector.h"
72#include "AliFlowTrackCuts.h"
35405405 73#include "AliFlowEvent.h"
a8f6c03f 74
75#include "AliAnalysisTaskSEHFv2.h"
76
77ClassImp(AliAnalysisTaskSEHFv2)
78
79
80//________________________________________________________________________
81AliAnalysisTaskSEHFv2::AliAnalysisTaskSEHFv2():
82AliAnalysisTaskSE(),
83 fhEventsInfo(0),
84 fOutput(0),
85 fRDCuts(0),
590dea9c 86 fLowmasslimit(1.669),
87 fUpmasslimit(2.069),
a8f6c03f 88 fNPtBins(1),
a8f6c03f 89 fNMassBins(200),
35405405 90 fReadMC(kFALSE),
ed352b10 91 fUseAfterBurner(kFALSE),
a8f6c03f 92 fDecChannel(0),
ed352b10 93 fAfterBurner(0),
590dea9c 94 fEventPlaneMeth(kTPCVZERO),
95 fEventPlanesComp(10),
05620b8d 96 fV0EPorder(2),
e68599a9 97 fMinCentr(20),
35405405 98 fMaxCentr(80),
5b559d1e 99 fEtaGap(kFALSE),
100 fSubEvents(0)
a8f6c03f 101{
102 // Default constructor
103}
104
105//________________________________________________________________________
590dea9c 106AliAnalysisTaskSEHFv2::AliAnalysisTaskSEHFv2(const char *name,AliRDHFCuts *rdCuts,Int_t decaychannel):
a8f6c03f 107 AliAnalysisTaskSE(name),
108 fhEventsInfo(0),
109 fOutput(0),
110 fRDCuts(rdCuts),
a8f6c03f 111 fLowmasslimit(0),
112 fUpmasslimit(0),
113 fNPtBins(1),
a8f6c03f 114 fNMassBins(200),
115 fReadMC(kFALSE),
ed352b10 116 fUseAfterBurner(kFALSE),
a8f6c03f 117 fDecChannel(decaychannel),
35405405 118 fAfterBurner(0),
590dea9c 119 fEventPlaneMeth(kTPCVZERO),
120 fEventPlanesComp(10),
05620b8d 121 fV0EPorder(2),
e68599a9 122 fMinCentr(20),
35405405 123 fMaxCentr(80),
5b559d1e 124 fEtaGap(kFALSE),
125 fSubEvents(0)
a8f6c03f 126{
ad34433d 127 // standard constructor
a8f6c03f 128 Int_t pdg=421;
129 switch(fDecChannel){
130 case 0:
131 pdg=411;
132 break;
133 case 1:
134 pdg=421;
135 break;
136 case 2:
137 pdg=413;
138 break;
139 case 3:
140 pdg=431;
141 break;
142 case 4:
143 pdg=421;
144 break;
145 case 5:
146 pdg=4122;
147 break;
148 }
ed352b10 149 fAfterBurner = new AliHFAfterBurner(fDecChannel);
a8f6c03f 150 if(pdg==413) SetMassLimits((Float_t)0.135,(Float_t)0.165);
151 else SetMassLimits((Float_t)0.2,pdg); //check range
152 fNPtBins=fRDCuts->GetNPtBins();
a8f6c03f 153
154 if(fDebug>1)fRDCuts->PrintAll();
155 // Output slot #1 writes into a TH1F container
156 DefineOutput(1,TH1F::Class()); //Info on the number of events etc.
157 // Output slot #2 writes into a TList container
158 DefineOutput(2,TList::Class()); //Main output
159 // Output slot #3 writes into a AliRDHFCuts container (cuts)
160 switch(fDecChannel){
161 case 0:
162 DefineOutput(3,AliRDHFCutsDplustoKpipi::Class()); //Cut object for Dplus
163 break;
164 case 1:
165 DefineOutput(3,AliRDHFCutsD0toKpi::Class()); //Cut object for D0
166 break;
167 case 2:
168 DefineOutput(3,AliRDHFCutsDStartoKpipi::Class()); //Cut object for D*
169 break;
170 }
171 //DefineOutput(4,AliFlowEventSimple::Class());
35405405 172 //DefineOutput(4,TList::Class());
a8f6c03f 173}
174
175//________________________________________________________________________
176AliAnalysisTaskSEHFv2::~AliAnalysisTaskSEHFv2()
177{
178 // Destructor
ad34433d 179 delete fOutput;
180 delete fhEventsInfo;
181 delete fRDCuts;
ad34433d 182 delete fAfterBurner;
a8f6c03f 183}
35405405 184//_________________________________________________________________
a8f6c03f 185void AliAnalysisTaskSEHFv2::SetMassLimits(Float_t range, Int_t pdg){
ad34433d 186 // Set limits for mass spectra plots
a8f6c03f 187 Float_t mass=0;
188 Int_t abspdg=TMath::Abs(pdg);
189 mass=TDatabasePDG::Instance()->GetParticle(abspdg)->Mass();
190 fUpmasslimit = mass+range;
191 fLowmasslimit = mass-range;
192}
193//_________________________________________________________________
194void AliAnalysisTaskSEHFv2::SetMassLimits(Float_t lowlimit, Float_t uplimit){
ad34433d 195 // Set limits for mass spectra plots
a8f6c03f 196 if(uplimit>lowlimit)
197 {
198 fUpmasslimit = uplimit;
199 fLowmasslimit = lowlimit;
200 }
201}
202
203
204//________________________________________________________________________
205void AliAnalysisTaskSEHFv2::LocalInit()
206{
207 // Initialization
208
209 if(fDebug > 1) printf("AnalysisTaskSEHFv2::Init() \n");
210
05620b8d 211 fRDCuts->SetMinCentrality(fMinCentr);
212 fRDCuts->SetMaxCentrality(fMaxCentr);
a8f6c03f 213
214 switch(fDecChannel){
215 case 0:
216 {
217 AliRDHFCutsDplustoKpipi* copycut=new AliRDHFCutsDplustoKpipi(*(static_cast<AliRDHFCutsDplustoKpipi*>(fRDCuts)));
218 // Post the data
219 PostData(3,copycut);
220 }
221 break;
222 case 1:
223 {
224 AliRDHFCutsD0toKpi* copycut=new AliRDHFCutsD0toKpi(*(static_cast<AliRDHFCutsD0toKpi*>(fRDCuts)));
225 // Post the data
226 PostData(3,copycut);
227 }
228 break;
229 case 2:
230 {
231 AliRDHFCutsDStartoKpipi* copycut=new AliRDHFCutsDStartoKpipi(*(static_cast<AliRDHFCutsDStartoKpipi*>(fRDCuts)));
232 // Post the data
233 PostData(3,copycut);
234 }
235 break;
236 default:
237 return;
238 }
239 return;
240}
241//________________________________________________________________________
242void AliAnalysisTaskSEHFv2::UserCreateOutputObjects()
243{
244 // Create the output container
245
246 if(fDebug > 1) printf("AnalysisTaskSEHFv2::UserCreateOutputObjects() \n");
247
65ee3f93 248 fhEventsInfo=new TH1F(GetOutputSlot(1)->GetContainer()->GetName(), "Number of AODs scanned",7,-0.5,6.5);
a8f6c03f 249 fhEventsInfo->GetXaxis()->SetBinLabel(1,"nEventsAnal");
250 fhEventsInfo->GetXaxis()->SetBinLabel(2,"nEvSelected");
251 fhEventsInfo->GetXaxis()->SetBinLabel(3,"nCandidatesSelected");
252 fhEventsInfo->GetXaxis()->SetBinLabel(4,"out of pt bounds");
590dea9c 253 fhEventsInfo->GetXaxis()->SetBinLabel(5,Form("Ev Sel in Centr %.0f-%.0f%s",fRDCuts->GetMinCentrality(),fRDCuts->GetMaxCentrality(),"%"));
254 fhEventsInfo->GetXaxis()->SetBinLabel(6,"mismatch lab");
65ee3f93 255 fhEventsInfo->GetXaxis()->SetBinLabel(7,"non valid TPC EP");
a8f6c03f 256 fhEventsInfo->GetXaxis()->SetNdivisions(1,kFALSE);
257
258
259 // Several histograms are more conveniently managed in a TList
260 fOutput = new TList();
261 fOutput->SetOwner();
262 fOutput->SetName("MainOutput");
263
05620b8d 264 for(Int_t icentr=fMinCentr+5;icentr<=fMaxCentr;icentr=icentr+5){
a8f6c03f 265 TString centrname;centrname.Form("centr%d_%d",icentr-5,icentr);
a8f6c03f 266
590dea9c 267 TH2F* hMPtCand=new TH2F(Form("hMPtCand%s",centrname.Data()),Form("Mass vs pt %s;pt (GeV);M (GeV/c^{2})",centrname.Data()),120,0,24.,fNMassBins,fLowmasslimit,fUpmasslimit);
268 fOutput->Add(hMPtCand);//For <pt> calculation
269
a8f6c03f 270
590dea9c 271 //Candidate distributions
272 for(Int_t i=0;i<fNPtBins;i++){
273 TH2F* hMc2phi=new TH2F(Form("hMc2phi_pt%d%s",i,centrname.Data()),Form("Mass vs cos2#Delta#phi (p_{t} bin %d %s);cos2#Delta#phi;M (GeV/c^{2})",i,centrname.Data()),100,-1.,1.,fNMassBins,fLowmasslimit,fUpmasslimit);
274 fOutput->Add(hMc2phi);//for 2D analysis
275
276 TH2F* hMphi=new TH2F(Form("hMphi_pt%d%s",i,centrname.Data()),Form("Mass vs #Delta#phi %s;#Delta#phi;M (GeV/c^{2})",centrname.Data()),96,0,TMath::Pi(),fNMassBins,fLowmasslimit,fUpmasslimit);
277 fOutput->Add(hMphi);//for phi bins analysis
278
a8f6c03f 279 if (fReadMC){
280 TH2F* hMc2phiS=new TH2F(Form("hMc2phiS_pt%d%s",i,centrname.Data()),Form("Mass vs cos2#Delta#phi (p_{t} bin %d %s);cos2#Delta#phi;M (GeV/c^{2})",i,centrname.Data()),100,-1.,1.,fNMassBins,fLowmasslimit,fUpmasslimit);
281 fOutput->Add(hMc2phiS);
282 TH2F * hMphiS=new TH2F(Form("hMphiS_pt%d%s",i,centrname.Data()),Form("Mass vs #Delta#phi (p_{t} bin %d %s);#Delta#phi;M (GeV/c^{2})",i,centrname.Data()),96,0,2*TMath::Pi(),fNMassBins,fLowmasslimit,fUpmasslimit);
283 fOutput->Add(hMphiS);
284 TH2F* hMc2phiB=new TH2F(Form("hMc2phiB_pt%d%s",i,centrname.Data()),Form("Mass vs cos2#Delta#phi (p_{t} bin %d %s);cos2#Delta#phi;M (GeV/c^{2})",i,centrname.Data()),100,-1.,1.,fNMassBins,fLowmasslimit,fUpmasslimit);
285 fOutput->Add(hMc2phiB);
286 TH2F * hMphiB=new TH2F(Form("hMphiB_pt%d%s",i,centrname.Data()),Form("Mass vs #Delta#phi (p_{t} bin %d %s);#Delta#phi;M (GeV/c^{2})",i,centrname.Data()),96,0,2*TMath::Pi(),fNMassBins,fLowmasslimit,fUpmasslimit);
287 fOutput->Add(hMphiB);
288 if((fDecChannel != AliAnalysisTaskSEHFv2::kDplustoKpipi) &&(fDecChannel != AliAnalysisTaskSEHFv2::kDstartoKpipi)){
289 TH2F* hMc2phiR=new TH2F(Form("hMc2phiR_pt%d%s",i,centrname.Data()),Form("Mass vs cos2#Delta#phi (p_{t} bin %d %s);cos2#Delta#phi;M (GeV/c^{2})",i,centrname.Data()),100,-1.,1.,fNMassBins,fLowmasslimit,fUpmasslimit);
290 fOutput->Add(hMc2phiR);
291 TH2F* hMphiR=new TH2F(Form("hMphiR_pt%d%s",i,centrname.Data()),Form("Mass vs #Delta#phi (p_{t} bin %d %s);#Delta#phi;M (GeV/c^{2})",i,centrname.Data()),96,0,2*TMath::Pi(),fNMassBins,fLowmasslimit,fUpmasslimit);
292 fOutput->Add(hMphiR);
293 }
294 }
295 }
a8f6c03f 296
590dea9c 297
298 //Event Plane
299 TH2F* hEvPlane=new TH2F(Form("hEvPlane%s",centrname.Data()),Form("VZERO/TPC Event plane angle %s;#phi Ev Plane (TPC);#phi Ev Plane (VZERO);Entries",centrname.Data()),200,0.,TMath::Pi(),200,0.,TMath::Pi());
300 fOutput->Add(hEvPlane);
e68599a9 301
65ee3f93 302 TH1F* hEvPlaneA=new TH1F(Form("hEvPlaneA%s",centrname.Data()),Form("Event plane angle %s;#phi Ev Plane;Entries",centrname.Data()),200,0.,TMath::Pi());
303 fOutput->Add(hEvPlaneA);
304
305 TH1F* hEvPlaneB=new TH1F(Form("hEvPlaneB%s",centrname.Data()),Form("Event plane angle %s;#phi Ev Plane;Entries",centrname.Data()),200,0.,TMath::Pi());
306 fOutput->Add(hEvPlaneB);
a8f6c03f 307
a8f6c03f 308 TH1F* hEvPlaneCand=new TH1F(Form("hEvPlaneCand%s",centrname.Data()),Form("Event plane angle - Event plane angle per candidate %s;#phi(Ev Plane Candidate);Entries",centrname.Data()),200,-TMath::Pi(),TMath::Pi());
309 fOutput->Add(hEvPlaneCand);
310
311 TH1F* hEvPlaneReso=new TH1F(Form("hEvPlaneReso%s",centrname.Data()),Form("Event plane angle Resolution %s;cos2(#psi_{A}-#psi_{B});Entries",centrname.Data()),220,-1.1,1.1);
312 fOutput->Add(hEvPlaneReso);
21cafb93 313 if(fEventPlaneMeth>kTPCVZERO){
65ee3f93 314 TH1F* hEvPlaneReso2=new TH1F(Form("hEvPlaneReso2%s",centrname.Data()),Form("Event plane angle Resolution %s;cos2(#psi_{A}-#psi_{B});Entries",centrname.Data()),220,-1.1,1.1);
315 fOutput->Add(hEvPlaneReso2);
316 TH1F* hEvPlaneReso3=new TH1F(Form("hEvPlaneReso3%s",centrname.Data()),Form("Event plane angle Resolution %s;cos2(#psi_{A}-#psi_{B});Entries",centrname.Data()),220,-1.1,1.1);
317 fOutput->Add(hEvPlaneReso3);
318 }
a8f6c03f 319 }
320
a8f6c03f 321 PostData(1,fhEventsInfo);
322 PostData(2,fOutput);
590dea9c 323
a8f6c03f 324 return;
325}
326
327//________________________________________________________________________
328void AliAnalysisTaskSEHFv2::UserExec(Option_t */*option*/)
329{
330 // Execute analysis for current event:
331 // heavy flavor candidates association to MC truth
a8f6c03f 332 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
333 if(fDebug>2) printf("Analysing decay %d\n",fDecChannel);
334 // Post the data already here
335 PostData(1,fhEventsInfo);
336 PostData(2,fOutput);
337
338 TClonesArray *arrayProng =0;
339 Int_t absPdgMom=0;
340 if(!aod && AODEvent() && IsStandardAOD()) {
341 // In case there is an AOD handler writing a standard AOD, use the AOD
342 // event in memory rather than the input (ESD) event.
343 aod = dynamic_cast<AliAODEvent*> (AODEvent());
344 // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
345 // have to taken from the AOD event hold by the AliAODExtension
346 AliAODHandler* aodHandler = (AliAODHandler*)
347 ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
348 if(aodHandler->GetExtensions()) {
349
350 AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
351 AliAODEvent *aodFromExt = ext->GetAOD();
352
353
354 if(fDecChannel==0){
355 absPdgMom=411;
356 arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
357 }
358 if(fDecChannel==1){
359 absPdgMom=421;
360 arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi");
361 }
362 if(fDecChannel==2){
363 absPdgMom=413;
364 arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");
365 }
366 }
471f1b0c 367 } else if(aod){
a8f6c03f 368 if(fDecChannel==0){
369 absPdgMom=411;
370 arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
371 }
372 if(fDecChannel==1){
373 absPdgMom=421;
374 arrayProng=(TClonesArray*)aod->GetList()->FindObject("D0toKpi");
375 }
376 if(fDecChannel==2){
377 absPdgMom=413;
378 arrayProng=(TClonesArray*)aod->GetList()->FindObject("Dstar");
379 }
380 }
381
10287b52 382 if(!aod || !arrayProng) {
a8f6c03f 383 AliError("AliAnalysisTaskSEHFv2::UserExec:Branch not found!\n");
384 return;
385 }
386
387 // fix for temporary bug in ESDfilter
388 // the AODs with null vertex pointer didn't pass the PhysSel
389 if(!aod->GetPrimaryVertex() || TMath::Abs(aod->GetMagneticField())<0.001) return;
390
391 TClonesArray *arrayMC=0;
392 AliAODMCHeader *mcHeader=0;
393
394 // load MC particles
395 if(fReadMC){
396
397 arrayMC = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
398 if(!arrayMC) {
399 AliWarning("AliAnalysisTaskSEHFv2::UserExec:MC particles branch not found!\n");
471f1b0c 400 return;
a8f6c03f 401 }
402
403 // load MC header
404 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
405 if(!mcHeader) {
406 AliError("AliAnalysisTaskSEHFv2::UserExec:MC header branch not found!\n");
407 return;
408 }
409 }
410
411 fhEventsInfo->Fill(0); // count event
412
413 AliAODRecoDecayHF *d=0;
414
415 Int_t nCand = arrayProng->GetEntriesFast();
416 if(fDebug>2) printf("Number of D2H: %d\n",nCand);
417
590dea9c 418 Bool_t isEvSel=fRDCuts->IsEventSelected(aod);
419 if(!isEvSel){
420 if(!fRDCuts->IsEventRejectedDueToCentrality())fhEventsInfo->Fill(4);
421 return;
422 }
a8f6c03f 423
590dea9c 424 fhEventsInfo->Fill(1);
1c451ce0 425 fhEventsInfo->Fill(4);
426
590dea9c 427 AliEventplane *pl=aod->GetEventplane();
428 if(!pl){
429 AliError("AliAnalysisTaskSEHFv2::UserExec:no eventplane! v2 analysis without eventplane not possible!\n");
a8f6c03f 430 return;
431 }
590dea9c 432
433 //Event plane
434 Double_t eventplane=0;
435 Double_t rpangleTPC=0;
436 Double_t rpangleVZERO=0;
437 Double_t planereso=0;
438 Double_t deltaPsi=0;
65ee3f93 439 Double_t rpangleeventB=0;
440 Double_t rpangleeventA=0;
590dea9c 441
442 //For candidate removal from TPC EP
443 TVector2* q=0x0;
35405405 444 TVector2 *qsub1=0x0;
445 TVector2 *qsub2=0x0;
446
a8f6c03f 447 //determine centrality bin
448 Float_t centr=fRDCuts->GetCentrality(aod);
449 Int_t icentr=0;
05620b8d 450 for(Int_t ic=fMinCentr+5;ic<=fMaxCentr;ic=ic+5){
a8f6c03f 451 if(ic>centr){
452 icentr=ic;
453 break;
454 }
455 }
456 TString centrbinname=Form("centr%d_%d",icentr-5,icentr);
590dea9c 457
ed352b10 458 if(fReadMC){
ed352b10 459 TRandom3 *g = new TRandom3(0);
65ee3f93 460 eventplane=g->Rndm()*TMath::Pi();
ed352b10 461 delete g;g=0x0;
590dea9c 462 if(fUseAfterBurner)fAfterBurner->SetEventPlane((Double_t)eventplane);
a8f6c03f 463 }else{
590dea9c 464 if(fEventPlaneMeth!=kTPC){
5b559d1e 465
590dea9c 466 //VZERO EP and resolution
65ee3f93 467 rpangleVZERO=GetPhi0Pi(pl->GetEventplane("V0",aod,fV0EPorder));
21cafb93 468 if(fEventPlaneMeth>kTPCVZERO){
5b559d1e 469 //Using V0A/C for VZERO resolution
470 if(fEventPlaneMeth==kVZEROA){
471 rpangleVZERO=GetPhi0Pi(pl->GetEventplane("V0A",aod,fV0EPorder));
472 if(fSubEvents!=kSingleV0Side)rpangleeventB=GetPhi0Pi(pl->GetEventplane("V0C",aod,fV0EPorder));
473 }
474 else if(fEventPlaneMeth==kVZEROC){
475 rpangleVZERO=GetPhi0Pi(pl->GetEventplane("V0C",aod,fV0EPorder));
476 if(fSubEvents!=kSingleV0Side)rpangleeventB=GetPhi0Pi(pl->GetEventplane("V0A",aod,fV0EPorder));
477 }
478 // deltaPsi =rpangleeventA-rpangleeventB;
479 eventplane=rpangleVZERO;
ed352b10 480 }
590dea9c 481 }
5b559d1e 482
65ee3f93 483 // TPC event plane
484 rpangleTPC = pl->GetEventplane("Q");
5b559d1e 485 if(rpangleTPC<0){
486 fhEventsInfo->Fill(6);
487 return;
488 }
489 rpangleeventA = rpangleTPC;
490 if(fSubEvents>kFullTPC||fEtaGap||fEventPlaneMeth==kVZERO){
491 qsub1 = pl->GetQsub1();
492 qsub2 = pl->GetQsub2();
493 if(!qsub1 || !qsub2){
494 fhEventsInfo->Fill(6);
495 return;
496 }
497 if(fSubEvents==kPosTPC||fSubEvents==kSingleV0Side||fEventPlaneMeth==kVZERO){
65ee3f93 498 rpangleeventA = qsub1->Phi()/2.;
5b559d1e 499 if(fSubEvents==kSingleV0Side||fEventPlaneMeth==kVZERO)rpangleeventB = qsub2->Phi()/2.;
590dea9c 500 }
5b559d1e 501 if(fSubEvents==kNegTPC)rpangleeventA = qsub2->Phi()/2.;
502 }
503
504 if(fEventPlaneMeth<=kTPCVZERO){
505 q = pl->GetQVector();
65ee3f93 506 deltaPsi = pl->GetQsubRes();
5b559d1e 507 }else{
508 deltaPsi=eventplane-rpangleeventA;
a8f6c03f 509 }
a8f6c03f 510 }
65ee3f93 511
512 if(TMath::Abs(deltaPsi)>TMath::Pi()/2.){
513 if(deltaPsi>0.) deltaPsi-=TMath::Pi();
514 else deltaPsi +=TMath::Pi();
515 } // difference of subevents reaction plane angle cannot be bigger than phi/2
516 planereso = TMath::Cos(2.*deltaPsi); // reaction plane resolution
590dea9c 517 if(TMath::Abs(rpangleTPC-rpangleVZERO)>fEventPlanesComp)return;
65ee3f93 518
21cafb93 519 if(fDebug>2)printf("Filling EP-related histograms\n");
590dea9c 520 //Filling EP-related histograms
590dea9c 521 ((TH2F*)fOutput->FindObject(Form("hEvPlane%s",centrbinname.Data())))->Fill(rpangleTPC,rpangleVZERO); // reaction plane angle without autocorrelations removal
522 ((TH1F*)fOutput->FindObject(Form("hEvPlaneReso%s",centrbinname.Data())))->Fill(planereso); //RP resolution
21cafb93 523 if(fEventPlaneMeth>kTPCVZERO||fEtaGap){
65ee3f93 524 ((TH1F*)fOutput->FindObject(Form("hEvPlaneA%s",centrbinname.Data())))->Fill(rpangleeventA); //Angle of first subevent
525 ((TH1F*)fOutput->FindObject(Form("hEvPlaneB%s",centrbinname.Data())))->Fill(rpangleeventB); //Angle of second subevent
526 }
21cafb93 527 if(fEventPlaneMeth>kTPCVZERO){
5b559d1e 528 Double_t deltaSub=rpangleeventA-rpangleeventB;
65ee3f93 529 if(TMath::Abs(deltaSub)>TMath::Pi()/2.){// difference of subevents reaction plane angle cannot be bigger than phi/2
530 if(deltaSub>0.) deltaSub-=TMath::Pi();
531 else deltaSub +=TMath::Pi();
532 }
533 ((TH1F*)fOutput->FindObject(Form("hEvPlaneReso2%s",centrbinname.Data())))->Fill(TMath::Cos(2.*deltaSub)); //RP resolution
5b559d1e 534 deltaSub =eventplane-rpangleeventB;
65ee3f93 535 if(TMath::Abs(deltaSub)>TMath::Pi()/2.){// difference of subevents reaction plane angle cannot be bigger than phi/2
536 if(deltaSub>0.) deltaSub-=TMath::Pi();
537 else deltaSub +=TMath::Pi();
538 }
539 ((TH1F*)fOutput->FindObject(Form("hEvPlaneReso3%s",centrbinname.Data())))->Fill(TMath::Cos(2.*deltaSub)); //RP resolution
540 }
a8f6c03f 541
21cafb93 542 if(fDebug>2)printf("Loop on D candidates\n");
590dea9c 543 //Loop on D candidates
a8f6c03f 544 for (Int_t iCand = 0; iCand < nCand; iCand++) {
a8f6c03f 545 d=(AliAODRecoDecayHF*)arrayProng->UncheckedAt(iCand);
a8f6c03f 546 Bool_t isSelBit=kTRUE;
547 if(fDecChannel==0) isSelBit=d->HasSelectionBit(AliRDHFCuts::kDplusCuts);
548 if(fDecChannel==1) isSelBit=d->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts);
549 if(fDecChannel==2) isSelBit=d->HasSelectionBit(AliRDHFCuts::kD0fromDstarCuts);
590dea9c 550 if(!isSelBit)continue;
551 Int_t ptbin=fRDCuts->PtBin(d->Pt());
552 if(ptbin<0) {
553 fhEventsInfo->Fill(3);
554 continue;
555 }
556 Bool_t isFidAcc = fRDCuts->IsInFiducialAcceptance(d->Pt(),d->Y(absPdgMom));
557 if(!isFidAcc)continue;
558 Int_t isSelected= fRDCuts->IsSelected(d,AliRDHFCuts::kCandidate,aod);
559 if(!isSelected)continue;
560
561 fhEventsInfo->Fill(2); // candidate selected
562 if(fDebug>3) printf("+++++++Is Selected\n");
a8f6c03f 563
590dea9c 564 Float_t* invMass=0x0;
565 Int_t nmasses;
566 CalculateInvMasses(d,invMass,nmasses);
a8f6c03f 567
21cafb93 568 if(fEventPlaneMeth<=kTPCVZERO){
590dea9c 569 eventplane = GetEventPlaneForCandidate(d,q,pl,qsub1,qsub2); // remove autocorrelations
570 ((TH1F*)fOutput->FindObject(Form("hEvPlaneCand%s",centrbinname.Data())))->Fill(rpangleTPC-eventplane);
571 }
a8f6c03f 572
590dea9c 573 Double_t phi=d->Phi();
574 if(fReadMC&&fUseAfterBurner)phi=fAfterBurner->GetNewAngle(d,arrayMC);
575 Float_t deltaphi=GetPhi0Pi(phi-eventplane);
576
577 //fill the histograms with the appropriate method
578 if(fDecChannel==0)FillDplus(d,arrayMC,ptbin,deltaphi,invMass,isSelected,icentr);
579 else if(fDecChannel==1)FillD02p(d,arrayMC,ptbin,deltaphi,invMass,isSelected,icentr);
580 else if(fDecChannel==2)FillDstar(d,arrayMC,ptbin,deltaphi,invMass,isSelected,icentr);
581
471f1b0c 582 delete [] invMass;
a8f6c03f 583 }
584
585 PostData(1,fhEventsInfo);
586 PostData(2,fOutput);
590dea9c 587
a8f6c03f 588 return;
589}
590
591//***************************************************************************
592
593// Methods used in the UserExec
594
595void AliAnalysisTaskSEHFv2::CalculateInvMasses(AliAODRecoDecayHF* d,Float_t*& masses,Int_t& nmasses){
596 //Calculates all the possible invariant masses for each candidate
597 //NB: the implementation for each candidate is responsibility of the corresponding developer
598
599 if(fDecChannel==0){
600 //D+ -- Giacomo
601 nmasses=1;
602 masses=new Float_t[nmasses];
603 Int_t pdgdaughters[3] = {211,321,211};
604 masses[0]=d->InvMass(3,(UInt_t*)pdgdaughters);
605 }
606 if(fDecChannel==1){
607 //D0 (Kpi) -- Chiara
608 const Int_t ndght=2;
609 nmasses=2;
610 masses=new Float_t[nmasses];
611 Int_t pdgdaughtersD0[ndght]={211,321};//pi,K
612 masses[0]=d->InvMass(ndght,(UInt_t*)pdgdaughtersD0); //D0
613 Int_t pdgdaughtersD0bar[ndght]={321,211};//K,pi
614 masses[1]=d->InvMass(ndght,(UInt_t*)pdgdaughtersD0bar); //D0bar
615 }
616 if(fDecChannel==2){
617 //D* -- Robert,Yifei, Alessandro
618 nmasses=1;
619 masses=new Float_t[nmasses];
620 masses[0]=((AliAODRecoCascadeHF*)d)->DeltaInvMass();
621 }
622}
623
624//******************************************************************************
625
590dea9c 626//Methods to fill the histograms, one for each channel
a8f6c03f 627//NB: the implementation for each candidate is responsibility of the corresponding developer
628
590dea9c 629//******************************************************************************
ad34433d 630void AliAnalysisTaskSEHFv2::FillDplus(AliAODRecoDecayHF* d,TClonesArray *arrayMC,Int_t ptbin,Float_t deltaphi, const Float_t* masses,Int_t isSel,Int_t icentr){
a8f6c03f 631 //D+ channel
632 if(!isSel){
633 if(fDebug>3)AliWarning("Candidate not selected\n");
634 return;
635 }
636 if(!masses){
637 if(fDebug>3)AliWarning("Masses not calculated\n");
638 return;
639 }
640
a8f6c03f 641 ((TH2F*)fOutput->FindObject(Form("hMc2phi_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[0]);
590dea9c 642 ((TH2F*)fOutput->FindObject(Form("hMphi_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[0]);
e68599a9 643 ((TH2F*)fOutput->FindObject(Form("hMPtCandcentr%d_%d",icentr-5,icentr)))->Fill(d->Pt(),masses[0]);
a8f6c03f 644 Int_t pdgdaughters[3] = {211,321,211};
645
646 if(fReadMC){
647 Int_t lab=-1;
ed352b10 648 if(fUseAfterBurner){
649 Bool_t isSignal=fAfterBurner->GetIsSignal();
650 if(isSignal)lab=10;
651 }else {
652 lab = d->MatchToMC(411,arrayMC,3,pdgdaughters);
653 }
a8f6c03f 654 if(lab>=0){ //signal
a8f6c03f 655 ((TH2F*)fOutput->FindObject(Form("hMc2phiS_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[0]);
656 ((TH2F*)fOutput->FindObject(Form("hMphiS_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[0]);
657 } else{ //background
a8f6c03f 658 ((TH2F*)fOutput->FindObject(Form("hMc2phiB_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[0]);
659 ((TH2F*)fOutput->FindObject(Form("hMphiB_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[0]);
660 }
661 }
662}
663
590dea9c 664//******************************************************************************
ad34433d 665void AliAnalysisTaskSEHFv2::FillD02p(AliAODRecoDecayHF* d,TClonesArray *arrayMC,Int_t ptbin,Float_t deltaphi, const Float_t* masses,Int_t isSel,Int_t icentr){
a8f6c03f 666
667 //D0->Kpi channel
668
669 //mass histograms
670 if(!masses){
671 if(fDebug>3)AliWarning("Masses not calculated\n");
672 return;
673 }
a8f6c03f 674 if(isSel==1 || isSel==3) {
a8f6c03f 675 ((TH2F*)fOutput->FindObject(Form("hMc2phi_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[0]);
590dea9c 676 ((TH2F*)fOutput->FindObject(Form("hMphi_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[0]);
e68599a9 677 ((TH2F*)fOutput->FindObject(Form("hMPtCandcentr%d_%d",icentr-5,icentr)))->Fill(d->Pt(),masses[0]);
a8f6c03f 678 }
679 if(isSel>=2) {
a8f6c03f 680 ((TH2F*)fOutput->FindObject(Form("hMc2phi_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[1]);
590dea9c 681 ((TH2F*)fOutput->FindObject(Form("hMphi_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[1]);
e68599a9 682 ((TH2F*)fOutput->FindObject(Form("hMPtCandcentr%d_%d",icentr-5,icentr)))->Fill(d->Pt(),masses[1]);
a8f6c03f 683 }
684
a8f6c03f 685 //MC histograms
686 if(fReadMC){
687
688 Int_t matchtoMC=-1;
689
690 //D0
691 Int_t pdgdaughters[2];
692 pdgdaughters[0]=211;//pi
693 pdgdaughters[1]=321;//K
694 Int_t nprongs=2;
695 Int_t absPdgMom=421;
696
697 matchtoMC = d->MatchToMC(absPdgMom,arrayMC,nprongs,pdgdaughters);
698
699 Int_t prongPdgPlus=421,prongPdgMinus=(-1)*421;
700 if((isSel==1 || isSel==3)){ //D0
701 if(matchtoMC>=0){
702 AliAODMCParticle *dMC = (AliAODMCParticle*)arrayMC->At(matchtoMC);
703 Int_t pdgMC = dMC->GetPdgCode();
704
705 if(pdgMC==prongPdgPlus) {
a8f6c03f 706 ((TH2F*)fOutput->FindObject(Form("hMc2phiS_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[0]);
707 ((TH2F*)fOutput->FindObject(Form("hMphiS_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[0]);
708 }
709 else {
a8f6c03f 710 ((TH2F*)fOutput->FindObject(Form("hMc2phiR_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[0]);
590dea9c 711 ((TH2F*)fOutput->FindObject(Form("hMphiR_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[0]);
a8f6c03f 712 }
713 } else {
a8f6c03f 714 ((TH2F*)fOutput->FindObject(Form("hMc2phiB_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[0]);
590dea9c 715 ((TH2F*)fOutput->FindObject(Form("hMphi_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[0]);
a8f6c03f 716 }
717 }
718 if(isSel>=2){ //D0bar
719 if(matchtoMC>=0){
720 AliAODMCParticle *dMC = (AliAODMCParticle*)arrayMC->At(matchtoMC);
721 Int_t pdgMC = dMC->GetPdgCode();
722
723 if(pdgMC==prongPdgMinus) {
a8f6c03f 724 ((TH2F*)fOutput->FindObject(Form("hMc2phiS_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[1]);
725 ((TH2F*)fOutput->FindObject(Form("hMphiS_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[1]);
726 }
727 else {
a8f6c03f 728 ((TH2F*)fOutput->FindObject(Form("hMc2phiR_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[1]);
590dea9c 729 ((TH2F*)fOutput->FindObject(Form("hMphi_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[1]);
a8f6c03f 730 }
731 } else {
a8f6c03f 732 ((TH2F*)fOutput->FindObject(Form("hMc2phiB_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[1]);
733 ((TH2F*)fOutput->FindObject(Form("hMphiB_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[1]);
734 }
735 }
736 }
737}
590dea9c 738//******************************************************************************
ad34433d 739void AliAnalysisTaskSEHFv2::FillDstar(AliAODRecoDecayHF* d,TClonesArray *arrayMC,Int_t ptbin,Float_t deltaphi, const Float_t* masses,Int_t isSel,Int_t icentr){
a8f6c03f 740 //D* channel
741 if(!isSel){
742 if(fDebug>3)AliWarning("Candidate not selected\n");
743 return;
744 }
745 if(!masses){
746 if(fDebug>3)AliWarning("Masses not calculated\n");
747 return;
748 }
590dea9c 749
a8f6c03f 750 ((TH2F*)fOutput->FindObject(Form("hMc2phi_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[0]);
590dea9c 751 ((TH2F*)fOutput->FindObject(Form("hMphi_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[0]);
e68599a9 752 ((TH2F*)fOutput->FindObject(Form("hMPtCandcentr%d_%d",icentr-5,icentr)))->Fill(d->Pt(),masses[0]);
a8f6c03f 753 Int_t pdgDgDStartoD0pi[2]={421,211};
754 Int_t pdgDgD0toKpi[2]={321,211};
755
756 if(fReadMC){
757 Int_t lab=-1;
758 lab = ((AliAODRecoCascadeHF*)d)->MatchToMC(413,421,pdgDgDStartoD0pi,pdgDgD0toKpi,arrayMC);
759 if(lab>=0){ //signal
a8f6c03f 760 ((TH2F*)fOutput->FindObject(Form("hMc2phiS_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[0]);
761 ((TH2F*)fOutput->FindObject(Form("hMphiS_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[0]);
762 } else{ //background
a8f6c03f 763 ((TH2F*)fOutput->FindObject(Form("hMc2phiB_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(TMath::Cos(2*deltaphi),masses[0]);
764 ((TH2F*)fOutput->FindObject(Form("hMphiB_pt%dcentr%d_%d",ptbin,icentr-5,icentr)))->Fill(deltaphi,masses[0]);
765 }
766 }
a8f6c03f 767}
768
a8f6c03f 769//________________________________________________________________________
590dea9c 770void AliAnalysisTaskSEHFv2::SetEventPlaneMethod(Int_t method){
1c451ce0 771 if(method>kVZEROC||method<0){
590dea9c 772 AliWarning("No EP method associated to the selection, setting to TPC EP\n");
773 method=kTPCVZERO;
a8f6c03f 774 }
590dea9c 775 fEventPlaneMeth=method;
a8f6c03f 776}
590dea9c 777
a8f6c03f 778//________________________________________________________________________
779Float_t AliAnalysisTaskSEHFv2::GetPhi0Pi(Float_t phi){
ad34433d 780 // Sets the phi angle in the range 0-pi
a8f6c03f 781 Float_t result=phi;
782 while(result<0){
783 result=result+TMath::Pi();
784 }
785 while(result>TMath::Pi()){
786 result=result-TMath::Pi();
787 }
590dea9c 788 return result;
a8f6c03f 789}
790
791//________________________________________________________________________
ad34433d 792Float_t AliAnalysisTaskSEHFv2::GetEventPlaneForCandidate(AliAODRecoDecayHF* d, const TVector2* q,AliEventplane *pl, const TVector2* qsub1, const TVector2* qsub2){
a8f6c03f 793 // remove autocorrelations
794
35405405 795 TArrayF* qx = 0x0;
796 TArrayF* qy = 0x0;
797 TVector2 qcopy;
798 if(!fEtaGap){
799 qx = pl->GetQContributionXArray();
800 qy = pl->GetQContributionYArray();
801 qcopy = *q;
802 }
803 else {
804 if(d->Eta()>0.){
805 qx = pl->GetQContributionXArraysub1();
806 qy = pl->GetQContributionYArraysub1();
807 qcopy = *qsub1;
808 }
809 else{
810 qx = pl->GetQContributionXArraysub2();
811 qy = pl->GetQContributionYArraysub2();
812 qcopy = *qsub2;
813 }
814 }
815
a8f6c03f 816 if(fDecChannel==2){
817 //D* -- Yifei, Alessandro,Robert
818 AliAODRecoDecayHF2Prong* theD0particle = ((AliAODRecoCascadeHF*)d)->Get2Prong();
819 AliAODTrack *track0 = (AliAODTrack*)theD0particle->GetDaughter(0);
820 AliAODTrack *track1 = (AliAODTrack*)theD0particle->GetDaughter(1);
821 AliAODTrack *track2 = ((AliAODRecoCascadeHF*)d)->GetBachelor();
822 // reduce global q vector
823
824 TVector2 q0;
825 if((track0->GetID()) < qx->fN){
826 q0.Set(qx->At(track0->GetID()),qy->At(track0->GetID()));}
827
828 TVector2 q1;
829 if((track1->GetID()) < qx->fN){
830 q1.Set(qx->At(track1->GetID()),qy->At(track1->GetID()));}
831
832 TVector2 q2;
833 if((track2->GetID()) < qx->fN){
834 q2.Set(qx->At(track2->GetID()),qy->At(track2->GetID()));}
835
836 qcopy = qcopy -(q0+q1+q2);
837
838 }
839
840 // reduce Q vector for D+ and D0
841
842 if(fDecChannel==1){
843 //D0 -- Chiara
844 AliAODTrack *track0 = (AliAODTrack*)d->GetDaughter(0);
845 AliAODTrack *track1 = (AliAODTrack*)d->GetDaughter(1);
846
847 TVector2 q0;
848 if((track0->GetID()) < qx->fN){
849 q0.Set(qx->At(track0->GetID()),qy->At(track0->GetID()));}
850
851 TVector2 q1;
852 if((track1->GetID()) < qx->fN){
853 q1.Set(qx->At(track1->GetID()),qy->At(track1->GetID()));}
854
855 qcopy = qcopy -(q0+q1);
856 }
857
858 if(fDecChannel==0){
859 //D+ -- Giacomo
860 AliAODTrack *track0 = (AliAODTrack*)d->GetDaughter(0);
861 AliAODTrack *track1 = (AliAODTrack*)d->GetDaughter(1);
862 AliAODTrack *track2 = (AliAODTrack*)d->GetDaughter(2);
863
864 TVector2 q0;
865 if((track0->GetID()) < qx->fN){
866 q0.Set(qx->At(track0->GetID()),qy->At(track0->GetID()));}
867
868 TVector2 q1;
869 if((track1->GetID()) < qx->fN){
870 q1.Set(qx->At(track1->GetID()),qy->At(track1->GetID()));}
871
872 TVector2 q2;
873 if((track2->GetID()) < qx->fN){
874 q2.Set(qx->At(track2->GetID()),qy->At(track2->GetID()));}
875
876 qcopy = qcopy -(q0+q1+q2);
877
878 }
879
880 return qcopy.Phi()/2.;
881
882}
590dea9c 883// //________________________________________________________________________
884// Float_t AliAnalysisTaskSEHFv2::GetEventPlaneFromV0(AliAODEvent *aodEvent){
885// // Compute event plane for VZERO - Obsolete, used for 2010 data
886
887// Int_t centr=fRDCuts->GetCentrality(aodEvent);
888// centr=centr-centr%10;
889// //temporary fix
890// if(centr<20)centr=20;
891// if(centr>70)centr=70;
892// //end temporary fix
893// Int_t binx=0;
894// Int_t iParHist=(centr-20)/10;
895
896// TString name;name.Form("parhist%d_%d",centr,centr+10);
897
898// if(fDebug>15)printf("EPfromV0 centr %d, iparhist %d (%p-%p)\n",centr,iParHist,fParHist->FindObject(name.Data()),fParHist->At(iParHist));
899
900// Int_t runnumber=aodEvent->GetRunNumber();
901// if(fParHist->At(iParHist)){
902// for(Int_t i=1;i<=((TH2D*)fParHist->At(iParHist))->GetNbinsX()&&binx<=0;i++){
903// Int_t run=atoi(((TH2D*)fParHist->At(iParHist))->GetXaxis()->GetBinLabel(i));
904// if(run>=runnumber)binx=i;
905// }
906// }else{
907// fhEventsInfo->Fill(7);
908// }
a8f6c03f 909
590dea9c 910// AliFlowTrackCuts* cutsRP = AliFlowTrackCuts::GetStandardVZEROOnlyTrackCuts();
911// cutsRP->SetEvent(aodEvent, MCEvent());//, 0x0);
912// cutsRP->SetName( Form("rp_cuts") );
913// AliFlowTrackCuts* dummy = new AliFlowTrackCuts("null_cuts");
914// dummy->SetParamType(AliFlowTrackCuts::kGlobal);
915// dummy->SetPtRange(+1,-1); // select nothing QUICK
916// dummy->SetEtaRange(+1,-1); // select nothing VZERO
917// dummy->SetEvent(aodEvent,MCEvent());
918
919// //////////////// construct the flow event container ////////////
920// AliFlowEvent flowEvent(cutsRP,dummy);
921// flowEvent.SetReferenceMultiplicity( 64 );
922// for(Int_t i=0;i<64&&binx>0;i++){
923// AliFlowTrack *flowTrack=flowEvent.GetTrack(i);
924// Double_t inte=((TH2D*)fParHist->At(iParHist))->Integral(binx,binx,i+1,i+1);
925// if(inte>0)flowTrack->SetWeight(flowTrack->Weight()/inte);
926// }
927// if(fDebug>15)printf("EPfromV0 flow tracks weights done\n");
ed352b10 928
590dea9c 929// AliFlowVector qvec=flowEvent.GetQ(fV0EPorder);
930// Double_t angleEP=(1./(Double_t)fV0EPorder)*qvec.Phi();
931// if(fDebug>15)printf("EPfromV0 phi %f\n",angleEP);
932// return angleEP;
933// }
a8f6c03f 934//________________________________________________________________________
935void AliAnalysisTaskSEHFv2::Terminate(Option_t */*option*/)
936{
a8f6c03f 937 // Terminate analysis
938 //
939 if(fDebug > 1) printf("AnalysisTaskSEHFv2: Terminate() \n");
590dea9c 940
a8f6c03f 941 return;
942}
943