Fix for leak
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliAnalysisTaskCombinHF.cxx
CommitLineData
811e9cbd 1/**************************************************************************
2 * Copyright(c) 1998-2018, 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/* $Id: $ */
17
18//*************************************************************************
19// Class AliAnalysisTaskCombinHF
20// AliAnalysisTaskSE to build D meson candidates by combining tracks
21// background is computed LS and track rotations is
22// Authors: F. Prino, A. Rossi
23/////////////////////////////////////////////////////////////
24
25#include <TList.h>
26#include <TH1F.h>
ebc460de 27#include <TH2F.h>
811e9cbd 28#include <TH3F.h>
29#include <THnSparse.h>
30
31#include "AliAnalysisManager.h"
32#include "AliInputEventHandler.h"
33#include "AliPIDResponse.h"
34#include "AliAODHandler.h"
35#include "AliAODEvent.h"
36#include "AliAODMCParticle.h"
37#include "AliAODMCHeader.h"
38#include "AliAODVertex.h"
39#include "AliAODTrack.h"
ebc460de 40#include "AliVertexingHFUtils.h"
811e9cbd 41#include "AliAnalysisTaskCombinHF.h"
42
43ClassImp(AliAnalysisTaskCombinHF)
44
45
46//________________________________________________________________________
47AliAnalysisTaskCombinHF::AliAnalysisTaskCombinHF():
48 AliAnalysisTaskSE(),
49 fOutput(0x0),
50 fHistNEvents(0x0),
51 fHistTrackStatus(0x0),
ebc460de 52 fHistCheckOrigin(0x0),
53 fHistCheckOriginSel(0x0),
54 fHistCheckDecChan(0x0),
55 fHistCheckDecChanAcc(0x0),
56 fPtVsYGen(0x0),
57 fPtVsYGenLimAcc(0x0),
58 fPtVsYGenAcc(0x0),
59 fPtVsYReco(0x0),
811e9cbd 60 fMassVsPtVsY(0x0),
61 fMassVsPtVsYRot(0x0),
62 fMassVsPtVsYLSpp(0x0),
63 fMassVsPtVsYLSmm(0x0),
64 fMassVsPtVsYSig(0x0),
65 fMassVsPtVsYRefl(0x0),
66 fMassVsPtVsYBkg(0x0),
67 fNSelected(0x0),
68 fNormRotated(0x0),
69 fDeltaMass(0x0),
70 fDeltaMassFullAnalysis(0x0),
71 fFilterMask(BIT(4)),
72 fTrackCutsAll(0x0),
73 fTrackCutsPion(0x0),
74 fTrackCutsKaon(0x0),
75 fPidHF(new AliAODPidHF()),
76 fAnalysisCuts(0x0),
77 fMinMass(1.750),
78 fMaxMass(2.150),
ebc460de 79 fEtaAccCut(0.9),
80 fPtAccCut(0.1),
811e9cbd 81 fNRotations(9),
82 fMinAngleForRot(5*TMath::Pi()/6),
83 fMaxAngleForRot(7*TMath::Pi()/6),
84 fMinAngleForRot3(2*TMath::Pi()/6),
85 fMaxAngleForRot3(4*TMath::Pi()/6),
86 fCounter(0x0),
87 fMeson(kDzero),
88 fReadMC(kFALSE),
ebc460de 89 fPromptFeeddown(kPrompt),
90 fGoUpToQuark(kTRUE),
811e9cbd 91 fFullAnalysis(0)
92{
93 // default constructor
94}
95
96//________________________________________________________________________
97AliAnalysisTaskCombinHF::AliAnalysisTaskCombinHF(Int_t meson, AliRDHFCuts* analysiscuts):
98 AliAnalysisTaskSE("DmesonCombin"),
99 fOutput(0x0),
100 fHistNEvents(0x0),
101 fHistTrackStatus(0x0),
ebc460de 102 fHistCheckOrigin(0x0),
103 fHistCheckOriginSel(0x0),
104 fHistCheckDecChan(0x0),
105 fHistCheckDecChanAcc(0x0),
106 fPtVsYGen(0x0),
107 fPtVsYGenLimAcc(0x0),
108 fPtVsYGenAcc(0x0),
109 fPtVsYReco(0x0),
811e9cbd 110 fMassVsPtVsY(0x0),
111 fMassVsPtVsYRot(0x0),
112 fMassVsPtVsYLSpp(0x0),
113 fMassVsPtVsYLSmm(0x0),
114 fMassVsPtVsYSig(0x0),
115 fMassVsPtVsYRefl(0x0),
116 fMassVsPtVsYBkg(0x0),
117 fNSelected(0x0),
118 fNormRotated(0x0),
119 fDeltaMass(0x0),
120 fDeltaMassFullAnalysis(0x0),
121 fFilterMask(BIT(4)),
122 fTrackCutsAll(0x0),
123 fTrackCutsPion(0x0),
124 fTrackCutsKaon(0x0),
125 fPidHF(new AliAODPidHF()),
126 fAnalysisCuts(analysiscuts),
127 fMinMass(1.750),
128 fMaxMass(2.150),
ebc460de 129 fEtaAccCut(0.9),
130 fPtAccCut(0.1),
811e9cbd 131 fNRotations(9),
132 fMinAngleForRot(5*TMath::Pi()/6),
133 fMaxAngleForRot(7*TMath::Pi()/6),
134 fMinAngleForRot3(2*TMath::Pi()/6),
135 fMaxAngleForRot3(4*TMath::Pi()/6),
136 fCounter(0x0),
137 fMeson(meson),
138 fReadMC(kFALSE),
ebc460de 139 fPromptFeeddown(1),
140 fGoUpToQuark(kTRUE),
811e9cbd 141 fFullAnalysis(0)
142
143{
144 // standard constructor
145
146 DefineOutput(1,TList::Class()); //My private output
147 DefineOutput(2,AliNormalizationCounter::Class());
148 }
149
150//________________________________________________________________________
151AliAnalysisTaskCombinHF::~AliAnalysisTaskCombinHF()
152{
153 //
154 // Destructor
155 //
156 delete fOutput;
157 delete fHistNEvents;
158 delete fHistTrackStatus;
ebc460de 159 delete fHistCheckOrigin;
160 delete fHistCheckOriginSel;
161 delete fHistCheckDecChan;
162 delete fHistCheckDecChanAcc;
163 delete fPtVsYGen;
164 delete fPtVsYGenLimAcc;
165 delete fPtVsYGenAcc;
166 delete fPtVsYReco;
811e9cbd 167 delete fMassVsPtVsY;
168 delete fMassVsPtVsYLSpp;
169 delete fMassVsPtVsYLSmm;
170 delete fMassVsPtVsYRot;
171 delete fMassVsPtVsYSig;
172 delete fMassVsPtVsYRefl;
173 delete fMassVsPtVsYBkg;
174 delete fNSelected;
175 delete fNormRotated;
176 delete fDeltaMass;
177 delete fCounter;
178 delete fTrackCutsAll;
179 delete fTrackCutsPion;
180 delete fTrackCutsKaon;
181 delete fPidHF;
182 delete fAnalysisCuts;
183}
184
185//________________________________________________________________________
186void AliAnalysisTaskCombinHF::UserCreateOutputObjects()
187{
188 // Create the output container
189 //
190 if(fDebug > 1) printf("AnalysisTaskCombinHF::UserCreateOutputObjects() \n");
191
192 fOutput = new TList();
193 fOutput->SetOwner();
194 fOutput->SetName("OutputHistos");
195
196 fHistNEvents = new TH1F("hNEvents", "number of events ",8,-0.5,7.5);
197 fHistNEvents->GetXaxis()->SetBinLabel(1,"nEventsAnal");
198 fHistNEvents->GetXaxis()->SetBinLabel(2,"n. passing IsEvSelected");
199 fHistNEvents->GetXaxis()->SetBinLabel(3,"n. rejected due to trigger");
200 fHistNEvents->GetXaxis()->SetBinLabel(4,"n. rejected due to not reco vertex");
201 fHistNEvents->GetXaxis()->SetBinLabel(5,"n. rejected for contr vertex");
202 fHistNEvents->GetXaxis()->SetBinLabel(6,"n. rejected for vertex out of accept");
203 fHistNEvents->GetXaxis()->SetBinLabel(7,"n. rejected for pileup events");
204 fHistNEvents->GetXaxis()->SetBinLabel(8,"no. of out centrality events");
205
206 fHistNEvents->GetXaxis()->SetNdivisions(1,kFALSE);
207 fHistNEvents->Sumw2();
208 fHistNEvents->SetMinimum(0);
209 fOutput->Add(fHistNEvents);
210
211 fHistTrackStatus = new TH1F("hTrackStatus", "",8,-0.5,7.5);
212 fHistTrackStatus->GetXaxis()->SetBinLabel(1,"Not OK");
213 fHistTrackStatus->GetXaxis()->SetBinLabel(2,"Track OK");
214 fHistTrackStatus->GetXaxis()->SetBinLabel(3,"Kaon, Not OK");
215 fHistTrackStatus->GetXaxis()->SetBinLabel(4,"Kaon OK");
216 fHistTrackStatus->GetXaxis()->SetBinLabel(5,"Pion, Not OK");
217 fHistTrackStatus->GetXaxis()->SetBinLabel(6,"Pion OK");
218 fHistTrackStatus->GetXaxis()->SetBinLabel(7,"Kaon||Pion, Not OK");
219 fHistTrackStatus->GetXaxis()->SetBinLabel(8,"Kaon||Pion OK");
2b6a376c 220 fHistTrackStatus->GetXaxis()->SetNdivisions(1,kFALSE);
811e9cbd 221 fHistTrackStatus->Sumw2();
222 fHistTrackStatus->SetMinimum(0);
223 fOutput->Add(fHistTrackStatus);
224
ebc460de 225 if(fReadMC){
226
227 fHistCheckOrigin=new TH1F("hCheckOrigin","",7,-1.5,5.5);
228 fHistCheckOrigin->Sumw2();
229 fHistCheckOrigin->SetMinimum(0);
230 fOutput->Add(fHistCheckOrigin);
231
232 fHistCheckOriginSel=new TH1F("hCheckOriginSel","",7,-1.5,5.5);
233 fHistCheckOriginSel->Sumw2();
234 fHistCheckOriginSel->SetMinimum(0);
235 fOutput->Add(fHistCheckOriginSel);
236
237 fHistCheckDecChan=new TH1F("hCheckDecChan","",7,-2.5,4.5);
238 fHistCheckDecChan->Sumw2();
239 fHistCheckDecChan->SetMinimum(0);
240 fOutput->Add(fHistCheckDecChan);
241
242 fHistCheckDecChanAcc=new TH1F("hCheckDecChanAcc","",7,-2.5,4.5);
243 fHistCheckDecChanAcc->Sumw2();
244 fHistCheckDecChanAcc->SetMinimum(0);
245 fOutput->Add(fHistCheckDecChanAcc);
246
247 fPtVsYGen= new TH2F("hPtVsYGen","",20,0.,10.,20,-1.,1.);
248 fPtVsYGen->Sumw2();
249 fPtVsYGen->SetMinimum(0);
250 fOutput->Add(fPtVsYGen);
251
252 fPtVsYGenLimAcc= new TH2F("hPtVsYGenLimAcc","",20,0.,10.,20,-1.,1.);
253 fPtVsYGenLimAcc->Sumw2();
254 fPtVsYGenLimAcc->SetMinimum(0);
255 fOutput->Add(fPtVsYGenLimAcc);
256
257 fPtVsYGenAcc= new TH2F("hPtVsYGenAcc","",20,0.,10.,20,-1.,1.);
258 fPtVsYGenAcc->Sumw2();
259 fPtVsYGenAcc->SetMinimum(0);
260 fOutput->Add(fPtVsYGenAcc);
261
262 fPtVsYReco= new TH2F("hPtVsYReco","",20,0.,10.,20,-1.,1.);
263 fPtVsYReco->Sumw2();
264 fPtVsYReco->SetMinimum(0);
265 fOutput->Add(fPtVsYReco);
266 }
267
268
811e9cbd 269 Int_t nMassBins=fMaxMass*1000.-fMinMass*1000.;
270 Double_t maxm=fMinMass+nMassBins*0.001;
271 fMassVsPtVsY=new TH3F("hMassVsPtVsY","",nMassBins,fMinMass,maxm,20,0.,10.,20,-1.,1.);
272 fMassVsPtVsY->Sumw2();
273 fMassVsPtVsY->SetMinimum(0);
274 fOutput->Add(fMassVsPtVsY);
275
276 fMassVsPtVsYRot=new TH3F("hMassVsPtVsYRot","",nMassBins,fMinMass,maxm,20,0.,10.,20,-1.,1.);
277 fMassVsPtVsYRot->Sumw2();
278 fMassVsPtVsYRot->SetMinimum(0);
279 fOutput->Add(fMassVsPtVsYRot);
280
281 if(fMeson==kDzero){
282 fMassVsPtVsYLSpp=new TH3F("hMassVsPtVsYLSpp","",nMassBins,fMinMass,maxm,20,0.,10.,20,-1.,1.);
283 fMassVsPtVsYLSpp->Sumw2();
284 fMassVsPtVsYLSpp->SetMinimum(0);
285 fOutput->Add(fMassVsPtVsYLSpp);
286 fMassVsPtVsYLSmm=new TH3F("hMassVsPtVsYLSmm","",nMassBins,fMinMass,maxm,20,0.,10.,20,-1.,1.);
287 fMassVsPtVsYLSmm->Sumw2();
288 fMassVsPtVsYLSmm->SetMinimum(0);
289 fOutput->Add(fMassVsPtVsYLSmm);
290 }
291
292 fMassVsPtVsYSig=new TH3F("hMassVsPtVsYSig","",nMassBins,fMinMass,maxm,20,0.,10.,20,-1.,1.);
293 fMassVsPtVsYSig->Sumw2();
294 fMassVsPtVsYSig->SetMinimum(0);
295 fOutput->Add(fMassVsPtVsYSig);
296
297 fMassVsPtVsYRefl=new TH3F("hMassVsPtVsYRefl","",nMassBins,fMinMass,maxm,20,0.,10.,20,-1.,1.);
298 fMassVsPtVsYRefl->Sumw2();
299 fMassVsPtVsYRefl->SetMinimum(0);
300 fOutput->Add(fMassVsPtVsYRefl);
301
302 fMassVsPtVsYBkg=new TH3F("hMassVsPtVsYBkg","",nMassBins,fMinMass,maxm,20,0.,10.,20,-1.,1.);
303 fMassVsPtVsYBkg->Sumw2();
304 fMassVsPtVsYBkg->SetMinimum(0);
305 fOutput->Add(fMassVsPtVsYBkg);
306
307 fNSelected=new TH1F("hNSelected","",100,-0.5,99.5);
308 fNSelected->Sumw2();
309 fNSelected->SetMinimum(0);
310 fOutput->Add(fNSelected);
311
312 fNormRotated=new TH1F("hNormRotated","",11,-0.5,10.5);
313 fNormRotated->Sumw2();
314 fNormRotated->SetMinimum(0);
315 fOutput->Add(fNormRotated);
316
317 fDeltaMass=new TH1F("hDeltaMass","",100,-0.4,0.4);
318 fDeltaMass->Sumw2();
319 fDeltaMass->SetMinimum(0);
320 fOutput->Add(fDeltaMass);
321
322 Int_t binSparseDMassRot[5]={nMassBins,100,24,40,20};
323 Double_t edgeLowSparseDMassRot[5]={fMinMass,-0.4,0.,-4.,0};
324 Double_t edgeHighSparseDMassRot[5]={maxm,0.4,12.,4.,3.14};
325 fDeltaMassFullAnalysis=new THnSparseF("fDeltaMassFullAnalysis","fDeltaMassFullAnalysis;inv mass (GeV/c);#Delta inv mass (GeV/c) ; p_{T}^{D} (GeV/c); #Delta p_{T} (GeV/c); daughter angle (2prongs) (rad);",5,binSparseDMassRot,edgeLowSparseDMassRot,edgeHighSparseDMassRot);
326 fOutput->Add(fDeltaMassFullAnalysis);
327
328 //Counter for Normalization
329 fCounter = new AliNormalizationCounter("NormalizationCounter");
330 fCounter->Init();
331
332 PostData(1,fOutput);
333 PostData(2,fCounter);
334}
335
336//________________________________________________________________________
337void AliAnalysisTaskCombinHF::UserExec(Option_t */*option*/){
338 //Build the 3-track combinatorics (+-+ and -+-) for D+->Kpipi decays
339
340 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
341 if(!aod && AODEvent() && IsStandardAOD()) {
342 // In case there is an AOD handler writing a standard AOD, use the AOD
343 // event in memory rather than the input (ESD) event.
344 aod = dynamic_cast<AliAODEvent*> (AODEvent());
345 }
346 if(!aod){
347 printf("AliAnalysisTaskCombinHF::UserExec: AOD not found!\n");
348 return;
349 }
350
351 // fix for temporary bug in ESDfilter
352 // the AODs with null vertex pointer didn't pass the PhysSel
353 if(!aod->GetPrimaryVertex() || TMath::Abs(aod->GetMagneticField())<0.001) return;
354 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
355 AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
356 AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
357 fPidHF->SetPidResponse(pidResp);
358
359
360 fHistNEvents->Fill(0); // count event
361 // Post the data already here
362 PostData(1,fOutput);
363
364 fCounter->StoreEvent(aod,fAnalysisCuts,fReadMC);
365
366 Bool_t isEvSel=fAnalysisCuts->IsEventSelected(aod);
367 if(fAnalysisCuts->IsEventRejectedDueToTrigger())fHistNEvents->Fill(2);
368 if(fAnalysisCuts->IsEventRejectedDueToNotRecoVertex())fHistNEvents->Fill(3);
369 if(fAnalysisCuts->IsEventRejectedDueToVertexContributors())fHistNEvents->Fill(4);
370 if(fAnalysisCuts->IsEventRejectedDueToZVertexOutsideFiducialRegion())fHistNEvents->Fill(5);
371 // if(fAnalysisCuts->IsEventRejectedDueToPileup())fHistNEvents->Fill(6);
372 if(fAnalysisCuts->IsEventRejectedDueToCentrality()) fHistNEvents->Fill(7);
373
374 if(!isEvSel)return;
375
376 fHistNEvents->Fill(1);
377
378 TClonesArray *arrayMC=0;
379 AliAODMCHeader *mcHeader=0;
380 if(fReadMC){
381 arrayMC = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
382 if(!arrayMC) {
383 printf("AliAnalysisTaskCombinHF::UserExec: MC particles branch not found!\n");
384 return;
385 }
386
387 // load MC header
388 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
389 if(!mcHeader) {
390 printf("AliAnalysisTaskCombinHF::UserExec: MC header branch not found!\n");
391 return;
392 }
ebc460de 393 FillGenHistos(arrayMC);
811e9cbd 394 }
395
396
397 Int_t ntracks=aod->GetNTracks();
398
399 // select and flag tracks
400 UChar_t* status = new UChar_t[ntracks];
401 for(Int_t iTr=0; iTr<ntracks; iTr++){
402 status[iTr]=0;
403 AliAODTrack* track=aod->GetTrack(iTr);
404 if(IsTrackSelected(track)) status[iTr]+=1;
405 if(IsKaon(track)) status[iTr]+=2;
406 if(IsPion(track)) status[iTr]+=4;
407 fHistTrackStatus->Fill(status[iTr]);
408 }
409
410 // build the combinatorics
411 Int_t nSelected=0;
412 Int_t nFiltered=0;
04731395 413 Double_t dummypos[3]={0.,0.,0.};
414 AliAODVertex* v2=new AliAODVertex(dummypos,999.,-1,2);
415 AliAODVertex* v3=new AliAODVertex(dummypos,999.,-1,3);
811e9cbd 416 // dummy values of track impact parameter, needed to build an AliAODRecoDecay object
417 Double_t d02[2]={0.,0.};
418 Double_t d03[3]={0.,0.,0.};
419 AliAODRecoDecay* tmpRD2 = new AliAODRecoDecay(0x0,2,0,d02);
420 AliAODRecoDecay* tmpRD3 = new AliAODRecoDecay(0x0,3,1,d03);
421 UInt_t pdg0[2]={321,211};
422 UInt_t pdgp[3]={321,211,211};
423 // UInt_t pdgs[3]={321,321,211};
424 Double_t tmpp[3];
425 Double_t px[3],py[3],pz[3];
426 Int_t dgLabels[3];
427
428 for(Int_t iTr1=0; iTr1<ntracks; iTr1++){
429 AliAODTrack* trK=aod->GetTrack(iTr1);
430 if((status[iTr1] & 1)==0) continue;
431 if((status[iTr1] & 2)==0) continue;
432 Int_t chargeK=trK->Charge();
433 trK->GetPxPyPz(tmpp);
434 px[0] = tmpp[0];
435 py[0] = tmpp[1];
436 pz[0] = tmpp[2];
437 dgLabels[0]=trK->GetLabel();
438 for(Int_t iTr2=0; iTr2<ntracks; iTr2++){
439 if((status[iTr2] & 1)==0) continue;
440 if((status[iTr2] & 4)==0) continue;
441 if(iTr1==iTr2) continue;
442 AliAODTrack* trPi1=aod->GetTrack(iTr2);
443 Int_t chargePi1=trPi1->Charge();
444 trPi1->GetPxPyPz(tmpp);
445 px[1] = tmpp[0];
446 py[1] = tmpp[1];
447 pz[1] = tmpp[2];
448 dgLabels[1]=trPi1->GetLabel();
449 if(chargePi1==chargeK){
450 if(fMeson==kDzero) FillLSHistos(421,2,tmpRD2,px,py,pz,pdg0,chargePi1);
451 continue;
452 }
453 if(fMeson==kDzero){
454 nFiltered++;
455 v2->AddDaughter(trK);
456 v2->AddDaughter(trPi1);
457 tmpRD2->SetSecondaryVtx(v2);
458 Bool_t ok=FillHistos(421,2,tmpRD2,px,py,pz,pdg0,arrayMC,dgLabels);
459 v2->RemoveDaughters();
460 if(ok) nSelected++;
461 }else{
462 for(Int_t iTr3=iTr2+1; iTr3<ntracks; iTr3++){
463 if((status[iTr3] & 1)==0) continue;
464 if((status[iTr3] & 4)==0) continue;
465 if(iTr1==iTr3) continue;
466 AliAODTrack* trPi2=aod->GetTrack(iTr3);
467 Int_t chargePi2=trPi2->Charge();
468 if(chargePi2==chargeK) continue;
469 nFiltered++;
470 trPi2->GetPxPyPz(tmpp);
471 px[2] = tmpp[0];
472 py[2] = tmpp[1];
473 pz[2] = tmpp[2];
474 dgLabels[2]=trPi2->GetLabel();
475 v3->AddDaughter(trK);
476 v3->AddDaughter(trPi1);
477 v3->AddDaughter(trPi2);
478 tmpRD3->SetSecondaryVtx(v3);
479 Bool_t ok=FillHistos(411,3,tmpRD3,px,py,pz,pdgp,arrayMC,dgLabels);
480 v3->RemoveDaughters();
481 if(ok) nSelected++;
482 }
483 }
484 }
485 }
486
487 delete [] status;
488 delete v2;
489 delete v3;
490 delete tmpRD2;
491 delete tmpRD3;
492
493 fNSelected->Fill(nSelected);
494
495 fCounter->StoreCandidates(aod,nFiltered,kTRUE);
496 fCounter->StoreCandidates(aod,nSelected,kFALSE);
497
498 PostData(1,fOutput);
499 PostData(2,fCounter);
500
501 return;
502}
503
504//________________________________________________________________________
505void AliAnalysisTaskCombinHF::FillLSHistos(Int_t pdgD,Int_t nProngs, AliAODRecoDecay* tmpRD, Double_t* px, Double_t* py, Double_t* pz, UInt_t *pdgdau, Int_t charge){
506 // Fill histos for LS candidates
507
508 tmpRD->SetPxPyPzProngs(nProngs,px,py,pz);
509 Double_t pt = tmpRD->Pt();
510 Double_t minv2 = tmpRD->InvMass2(nProngs,pdgdau);
511 if(minv2>fMinMass*fMinMass && minv2<fMaxMass*fMaxMass){
512 Double_t rapid = tmpRD->Y(pdgD);
513 if(fAnalysisCuts->IsInFiducialAcceptance(pt,rapid)){
514 if(charge>0) fMassVsPtVsYLSpp->Fill(TMath::Sqrt(minv2),pt,rapid);
515 else fMassVsPtVsYLSmm->Fill(TMath::Sqrt(minv2),pt,rapid);
516 }
517 }
518 return;
519}
520
ebc460de 521//________________________________________________________________________
522void AliAnalysisTaskCombinHF::FillGenHistos(TClonesArray* arrayMC){
523 // Fill histos with generated quantities
524 Int_t totPart=arrayMC->GetEntriesFast();
525 Int_t thePDG=411;
526 Int_t nProng=3;
527 if(fMeson==kDzero){
528 thePDG=421;
529 nProng=2;
530 }
ebc460de 531 for(Int_t ip=0; ip<totPart; ip++){
532 AliAODMCParticle *part = (AliAODMCParticle*)arrayMC->At(ip);
533 if(TMath::Abs(part->GetPdgCode())==thePDG){
534 Int_t orig=AliVertexingHFUtils::CheckOrigin(arrayMC,part,fGoUpToQuark);
535 fHistCheckOrigin->Fill(orig);
536 if(fPromptFeeddown==kFeeddown && orig!=5) continue;
537 else if(fPromptFeeddown==kPrompt && orig!=4) continue;
538 else if(fPromptFeeddown==kBoth && orig<4) continue;
539 fHistCheckOriginSel->Fill(orig);
540 Int_t deca=0;
541 Bool_t isGoodDecay=kFALSE;
f1f4f53a 542 Int_t labDau[4]={-1,-1,-1,-1};
ebc460de 543 if(fMeson==kDzero){
544 deca=AliVertexingHFUtils::CheckD0Decay(arrayMC,part,labDau);
545 if(part->GetNDaughters()!=2) continue;
546 if(deca==1) isGoodDecay=kTRUE;
547 }else if(fMeson==kDplus){
548 deca=AliVertexingHFUtils::CheckDplusDecay(arrayMC,part,labDau);
549 if(deca>0) isGoodDecay=kTRUE;
550 }
ebc460de 551 fHistCheckDecChan->Fill(deca);
f1f4f53a 552 if(labDau[0]==-1){
553 // printf(Form("Meson %d Label of daughters not filled correctly -- %d\n",fMeson,isGoodDecay));
554 continue; //protection against unfilled array of labels
555 }
556 Bool_t isInAcc=CheckAcceptance(arrayMC,nProng,labDau);
ebc460de 557 if(isInAcc) fHistCheckDecChanAcc->Fill(deca);
558 if(isGoodDecay){
559 Double_t ptgen=part->Pt();
560 Double_t ygen=part->Y();
561 if(fAnalysisCuts->IsInFiducialAcceptance(ptgen,ygen)){
562 fPtVsYGen->Fill(ptgen,ygen);
563 if(TMath::Abs(ygen)<0.5) fPtVsYGenLimAcc->Fill(ptgen,ygen);
564 if(isInAcc) fPtVsYGenAcc->Fill(ptgen,ygen);
565 }
566 }
567 }
568 }
569}
570
811e9cbd 571//________________________________________________________________________
572Bool_t AliAnalysisTaskCombinHF::FillHistos(Int_t pdgD,Int_t nProngs, AliAODRecoDecay* tmpRD, Double_t* px, Double_t* py, Double_t* pz, UInt_t *pdgdau, TClonesArray *arrayMC, Int_t* dgLabels){
573 // Fill histos for candidates with proper charge sign
574
575 Bool_t accept=kFALSE;
576
577 tmpRD->SetPxPyPzProngs(nProngs,px,py,pz);
578 Double_t pt = tmpRD->Pt();
579 Double_t minv2 = tmpRD->InvMass2(nProngs,pdgdau);
580 Double_t mass=TMath::Sqrt(minv2);
581
582 if(minv2>fMinMass*fMinMass && minv2<fMaxMass*fMaxMass){
583 Double_t rapid = tmpRD->Y(pdgD);
584 if(fAnalysisCuts->IsInFiducialAcceptance(pt,rapid)){
585 fMassVsPtVsY->Fill(mass,pt,rapid);
586 accept=kTRUE;
587 if(fReadMC){
588 Int_t signPdg[3]={0,0,0};
589 for(Int_t iii=0; iii<nProngs; iii++) signPdg[iii]=pdgdau[iii];
590 Int_t labD = tmpRD->MatchToMC(pdgD,arrayMC,nProngs,signPdg);
591 if(labD>=0){
ebc460de 592 AliAODMCParticle* part = dynamic_cast<AliAODMCParticle*>(arrayMC->At(TMath::Abs(dgLabels[0])));
811e9cbd 593 if(part){
594 Int_t pdgCode = TMath::Abs( part->GetPdgCode() );
ebc460de 595 if(pdgCode==321){
ebc460de 596 AliAODMCParticle* dmes = dynamic_cast<AliAODMCParticle*>(arrayMC->At(labD));
fcd3588f 597 if(dmes){
598 Int_t orig=AliVertexingHFUtils::CheckOrigin(arrayMC,dmes,fGoUpToQuark);
599 if((fPromptFeeddown==kFeeddown && orig==5)|| (fPromptFeeddown==kPrompt && orig==4) || (fPromptFeeddown==kBoth && orig>=4)) {
600 fPtVsYReco->Fill(dmes->Pt(),dmes->Y());
601 }
602 }
f1f4f53a 603 fMassVsPtVsYSig->Fill(mass,pt,rapid);
ebc460de 604 }else{
605 fMassVsPtVsYRefl->Fill(mass,pt,rapid);
606 }
811e9cbd 607 }
608 }else{
609 fMassVsPtVsYBkg->Fill(mass,pt,rapid);
610 }
611 }
612 }
613 }
614
615 Int_t nRotated=0;
616 Double_t massRot=0;// calculated later only if candidate is acceptable
617 Double_t angleProngXY;
618 if(TMath::Abs(pdgD)==421)angleProngXY=TMath::ACos((px[0]*px[1]+py[0]*py[1])/TMath::Sqrt((px[0]*px[0]+py[0]*py[0])*(px[1]*px[1]+py[1]*py[1])));
619 else {
620 angleProngXY=TMath::ACos(((px[0]+px[1])*px[2]+(py[0]+py[1])*py[2])/TMath::Sqrt(((px[0]+px[1])*(px[0]+px[1])+(py[0]+py[1])*(py[0]+py[1]))*(px[2]*px[2]+py[2]*py[2])));
621 }
622 Double_t ptOrig=pt;
623
624
625 Double_t rotStep=(fMaxAngleForRot-fMinAngleForRot)/(fNRotations-1); // -1 is to ensure that the last rotation is done with angle=fMaxAngleForRot
626 Double_t rotStep3=(fMaxAngleForRot3-fMinAngleForRot3)/(fNRotations-1); // -1 is to ensure that the last rotation is done with angle=fMaxAngleForRot
627
628 for(Int_t irot=0; irot<fNRotations; irot++){
629 Double_t phirot=fMinAngleForRot+rotStep*irot;
630 Double_t tmpx=px[0];
631 Double_t tmpy=py[0];
632 Double_t tmpx2=px[2];
633 Double_t tmpy2=py[2];
634 px[0]=tmpx*TMath::Cos(phirot)-tmpy*TMath::Sin(phirot);
635 py[0]=tmpx*TMath::Sin(phirot)+tmpy*TMath::Cos(phirot);
636 if(pdgD==411 || pdgD==431){
637 Double_t phirot2=fMinAngleForRot3+rotStep3*irot;
638 px[2]=tmpx*TMath::Cos(phirot2)-tmpy*TMath::Sin(phirot2);
639 py[2]=tmpx*TMath::Sin(phirot2)+tmpy*TMath::Cos(phirot2);
640 }
641 tmpRD->SetPxPyPzProngs(nProngs,px,py,pz);
642 pt = tmpRD->Pt();
643 minv2 = tmpRD->InvMass2(nProngs,pdgdau);
644 if(minv2>fMinMass*fMinMass && minv2<fMaxMass*fMaxMass){
645 Double_t rapid = tmpRD->Y(pdgD);
646 if(fAnalysisCuts->IsInFiducialAcceptance(pt,rapid)){
647 massRot=TMath::Sqrt(minv2);
648 fMassVsPtVsYRot->Fill(massRot,pt,rapid);
649 nRotated++;
650 fDeltaMass->Fill(massRot-mass);
651 if(fFullAnalysis){
652 Double_t pointRot[5]={mass,massRot-mass,ptOrig,pt-ptOrig,angleProngXY};
653 fDeltaMassFullAnalysis->Fill(pointRot);
654 }
655 }
656 }
657 px[0]=tmpx;
658 py[0]=tmpy;
659 if(pdgD==411 || pdgD==431){
660 px[2]=tmpx2;
661 py[2]=tmpy2;
662 }
663 }
664 fNormRotated->Fill(nRotated);
665
666 return accept;
667
668}
669//________________________________________________________________________
670Bool_t AliAnalysisTaskCombinHF::IsTrackSelected(AliAODTrack* track){
671 // track selection cuts
672
673 if(track->Charge()==0) return kFALSE;
674 if(!(track->TestFilterMask(fFilterMask))) return kFALSE;
675 if(!SelectAODTrack(track,fTrackCutsAll)) return kFALSE;
676 return kTRUE;
677}
678//________________________________________________________________________
679Bool_t AliAnalysisTaskCombinHF::IsKaon(AliAODTrack* track){
680 // kaon selection cuts
681
682 if(!fPidHF) return kTRUE;
683 Int_t isKaon=fPidHF->MakeRawPid(track,AliPID::kKaon);
684 if(isKaon>=0 && SelectAODTrack(track,fTrackCutsKaon)) return kTRUE;
685 return kFALSE;
686}
687//________________________________________________________________________
688Bool_t AliAnalysisTaskCombinHF::IsPion(AliAODTrack* track){
689 // pion selection cuts
690
691 if(!fPidHF) return kTRUE;
692 Int_t isPion=fPidHF->MakeRawPid(track,AliPID::kPion);
693 if(isPion>=0&& SelectAODTrack(track,fTrackCutsPion)) return kTRUE;
694 return kFALSE;
695}
696//________________________________________________________________________
697Bool_t AliAnalysisTaskCombinHF::SelectAODTrack(AliAODTrack *track, AliESDtrackCuts *cuts){
698 // AOD track selection
699
700 if(!cuts) return kTRUE;
701
702 AliESDtrack esdTrack(track);
703 // set the TPC cluster info
704 esdTrack.SetTPCClusterMap(track->GetTPCClusterMap());
705 esdTrack.SetTPCSharedMap(track->GetTPCSharedMap());
706 esdTrack.SetTPCPointsF(track->GetTPCNclsF());
707 if(!cuts->IsSelected(&esdTrack)) return kFALSE;
708
709 return kTRUE;
710}
711
ebc460de 712//_________________________________________________________________
713Bool_t AliAnalysisTaskCombinHF::CheckAcceptance(TClonesArray* arrayMC,Int_t nProng, Int_t *labDau){
714 // check if the decay products are in the good eta and pt range
715 for (Int_t iProng = 0; iProng<nProng; iProng++){
716 AliAODMCParticle* mcPartDaughter=dynamic_cast<AliAODMCParticle*>(arrayMC->At(labDau[iProng]));
717 if(!mcPartDaughter) return kFALSE;
718 Double_t eta = mcPartDaughter->Eta();
719 Double_t pt = mcPartDaughter->Pt();
720 if (TMath::Abs(eta) > fEtaAccCut || pt < fPtAccCut) return kFALSE;
721 }
722 return kTRUE;
723}
724
811e9cbd 725//_________________________________________________________________
726void AliAnalysisTaskCombinHF::Terminate(Option_t */*option*/)
727{
728 // Terminate analysis
729 //
730 if(fDebug > 1) printf("AliAnalysisTaskCombinHF: Terminate() \n");
731 fOutput = dynamic_cast<TList*> (GetOutputData(1));
732 if (!fOutput) {
733 printf("ERROR: fOutput not available\n");
734 return;
735 }
736 fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("hNEvents"));
737 if(fHistNEvents){
738 printf("Number of analyzed events = %d\n",(Int_t)fHistNEvents->GetBinContent(2));
739 }else{
740 printf("ERROR: fHistNEvents not available\n");
741 return;
742 }
743 return;
744}
745