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