1 /**************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 /////////////////////////////////////////////////////////////
18 // AliAnalysisTaskSE for the extraction of signal(e.g D+) of heavy flavor
19 // decay candidates with the MC truth.
20 // Authors: Renu Bala, bala@to.infn.it
21 // F. Prino, prino@to.infn.it
22 // G. Ortona, ortona@to.infn.it
23 /////////////////////////////////////////////////////////////
25 #include <TClonesArray.h>
32 #include <TDatabasePDG.h>
34 #include "AliAnalysisManager.h"
35 #include "AliAODHandler.h"
36 #include "AliAODEvent.h"
37 #include "AliAODVertex.h"
38 #include "AliAODTrack.h"
39 #include "AliAODMCHeader.h"
40 #include "AliAODMCParticle.h"
41 #include "AliAODRecoDecayHF3Prong.h"
42 #include "AliAnalysisVertexingHF.h"
43 #include "AliAnalysisTaskSE.h"
44 #include "AliAnalysisTaskSELambdac.h"
45 #include "AliKFParticle.h"
46 #include "AliAODPidHF.h"
47 #include "AliRDHFCutsLctopKpi.h"
48 #include "AliRDHFCuts.h"
50 ClassImp(AliAnalysisTaskSELambdac)
53 //________________________________________________________________________
54 AliAnalysisTaskSELambdac::AliAnalysisTaskSELambdac():
85 fTPCSignal3SigmaReK(0),
86 fTPCSignal3SigmaSedK(0),
87 fTPCSignal3SigmaRep(0),
88 fTPCSignal3SigmaSedp(0),
90 fTPCSignal2SigmaReK(0),
91 fTPCSignal2SigmaSedK(0),
92 fTPCSignal2SigmaRep(0),
93 fTPCSignal2SigmaSedp(0),
109 // Default constructor
112 //________________________________________________________________________
113 AliAnalysisTaskSELambdac::AliAnalysisTaskSELambdac(const char *name,Bool_t fillNtuple,AliRDHFCutsLctopKpi *lccutsana,AliRDHFCutsLctopKpi *lccutsprod):
114 AliAnalysisTaskSE(name),
119 fWellIDProt3sigSe(0),
120 fWellIDProt2sigSe(0),
121 fWellIDProt3sigRe(0),
122 fWellIDProt2sigRe(0),
124 fWellIDKaon3sigSe(0),
125 fWellIDKaon3sigRe(0),
127 fWellIDKaon2sigSe(0),
128 fWellIDKaon2sigRe(0),
144 fTPCSignal3SigmaReK(0),
145 fTPCSignal3SigmaSedK(0),
146 fTPCSignal3SigmaRep(0),
147 fTPCSignal3SigmaSedp(0),
149 fTPCSignal2SigmaReK(0),
150 fTPCSignal2SigmaSedK(0),
151 fTPCSignal2SigmaRep(0),
152 fTPCSignal2SigmaSedp(0),
155 fLowmasslimit(2.086),
157 fRDCutsAnalysis(lccutsana),
158 fRDCutsProduction(lccutsprod),
160 fFillNtuple(fillNtuple),
168 //Double_t ptlim[5]={0.,2.,3.,5,9999999.};
169 //SetPtBinLimit(5, ptlim);
170 SetPtBinLimit(fRDCutsAnalysis->GetNPtBins()+1,fRDCutsAnalysis->GetPtBinLimits());
171 // Default constructor
172 // Output slot #1 writes into a TList container
173 DefineOutput(1,TList::Class()); //My private output
174 DefineOutput(2,TList::Class());
177 // Output slot #2 writes into a TNtuple container
178 DefineOutput(3,TNtuple::Class()); //My private output
182 //________________________________________________________________________
183 AliAnalysisTaskSELambdac::~AliAnalysisTaskSELambdac()
196 delete fRDCutsAnalysis;
199 if(fRDCutsProduction){
200 delete fRDCutsProduction;
201 fRDCutsProduction = 0;
209 //_________________________________________________________________
210 void AliAnalysisTaskSELambdac::SetMassLimits(Float_t range){
211 fUpmasslimit = 2.286+range;
212 fLowmasslimit = 2.286-range;
214 //_________________________________________________________________
215 void AliAnalysisTaskSELambdac::SetMassLimits(Float_t lowlimit, Float_t uplimit){
218 fUpmasslimit = lowlimit;
219 fLowmasslimit = uplimit;
224 //________________________________________________________________________
225 void AliAnalysisTaskSELambdac::SetPtBinLimit(Int_t n, Float_t* lim){
226 // define pt bins for analysis
228 printf("Max. number of Pt bins = %d\n",kMaxPtBins);
230 fArrayBinLimits[0]=0.;
231 fArrayBinLimits[1]=2.;
232 fArrayBinLimits[2]=3.;
233 fArrayBinLimits[3]=5.;
234 for(Int_t i=4; i<kMaxPtBins+1; i++) fArrayBinLimits[i]=99999999.;
237 fArrayBinLimits[0]=lim[0];
238 for(Int_t i=1; i<fNPtBins+1; i++)
239 if(lim[i]>fArrayBinLimits[i-1]){
240 fArrayBinLimits[i]=lim[i];
243 fArrayBinLimits[i]=fArrayBinLimits[i-1];
245 for(Int_t i=fNPtBins; i<kMaxPtBins+1; i++) fArrayBinLimits[i]=99999999.;
248 printf("Number of Pt bins = %d\n",fNPtBins);
249 for(Int_t i=0; i<fNPtBins+1; i++) printf(" Bin%d = %8.2f-%8.2f\n",i,fArrayBinLimits[i],fArrayBinLimits[i+1]);
252 //_________________________________________________________________
253 Double_t AliAnalysisTaskSELambdac::GetPtBinLimit(Int_t ibin){
254 if(ibin>fNPtBins)return -1;
255 return fArrayBinLimits[ibin];
258 //_________________________________________________________________
259 void AliAnalysisTaskSELambdac::Init()
263 if(fDebug > 1) printf("AnalysisTaskSELambdac::Init() \n");
265 fListCuts=new TList();
266 AliRDHFCutsLctopKpi *production = new AliRDHFCutsLctopKpi();
267 production=fRDCutsProduction;
269 AliRDHFCutsLctopKpi *analysis = new AliRDHFCutsLctopKpi();
270 analysis=fRDCutsAnalysis;
272 fListCuts->Add(analysis);
273 fListCuts->Add(production);
274 PostData(2,fListCuts);
278 //________________________________________________________________________
279 void AliAnalysisTaskSELambdac::UserCreateOutputObjects()
281 // Create the output container
283 if(fDebug > 1) printf("AnalysisTaskSELambdac::UserCreateOutputObjects() \n");
285 // Several histograms are more conveniently managed in a TList
286 fOutput = new TList();
288 fOutput->SetName("OutputHistos");
293 for(Int_t i=0;i<fNPtBins;i++){
295 index=GetHistoIndex(i);
296 indexLS=GetLSHistoIndex(i);
298 hisname.Form("hMassPt%d",i);
299 fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
300 fMassHist[index]->Sumw2();
301 hisname.Form("hCosPAllPt%d",i);
302 fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
303 fCosPHist[index]->Sumw2();
304 hisname.Form("hDLenAllPt%d",i);
305 fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
306 fDLenHist[index]->Sumw2();
307 hisname.Form("hSumd02AllPt%d",i);
308 fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
309 fSumd02Hist[index]->Sumw2();
310 hisname.Form("hSigVertAllPt%d",i);
311 fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
312 fSigVertHist[index]->Sumw2();
313 hisname.Form("hPtMaxAllPt%d",i);
314 fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
315 fPtMaxHist[index]->Sumw2();
317 hisname.Form("hDCAAllPt%d",i);
318 fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
319 fDCAHist[index]->Sumw2();
323 hisname.Form("hMassPt%dTC",i);
324 fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
325 fMassHistTC[index]->Sumw2();
331 hisname.Form("hCosPAllPt%dLS",i);
332 fCosPHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
333 fCosPHistLS[index]->Sumw2();
334 hisname.Form("hDLenAllPt%dLS",i);
335 fDLenHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
336 fDLenHistLS[index]->Sumw2();
337 hisname.Form("hSumd02AllPt%dLS",i);
338 fSumd02HistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
339 fSumd02HistLS[index]->Sumw2();
340 hisname.Form("hSigVertAllPt%dLS",i);
341 fSigVertHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
342 fSigVertHistLS[index]->Sumw2();
343 hisname.Form("hPtMaxAllPt%dLS",i);
344 fPtMaxHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
345 fPtMaxHistLS[index]->Sumw2();
347 hisname.Form("hDCAAllPt%dLS",i);
348 fDCAHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
349 fDCAHistLS[index]->Sumw2();
351 hisname.Form("hLSPt%dLC",i);
352 fMassHistLS[indexLS] = new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
353 fMassHistLS[indexLS]->Sumw2();
355 hisname.Form("hLSPt%dTC",i);
356 fMassHistLSTC[indexLS] = new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
357 fMassHistLSTC[indexLS]->Sumw2();
361 index=GetSignalHistoIndex(i);
363 hisname.Form("hSigPt%d",i);
364 fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
365 fMassHist[index]->Sumw2();
366 hisname.Form("hCosPSigPt%d",i);
367 fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
368 fCosPHist[index]->Sumw2();
369 hisname.Form("hDLenSigPt%d",i);
370 fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
371 fDLenHist[index]->Sumw2();
372 hisname.Form("hSumd02SigPt%d",i);
373 fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
374 fSumd02Hist[index]->Sumw2();
375 hisname.Form("hSigVertSigPt%d",i);
376 fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
377 fSigVertHist[index]->Sumw2();
378 hisname.Form("hPtMaxSigPt%d",i);
379 fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
380 fPtMaxHist[index]->Sumw2();
382 hisname.Form("hDCASigPt%d",i);
383 fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
384 fDCAHist[index]->Sumw2();
387 hisname.Form("hSigPt%dTC",i);
388 fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
389 fMassHistTC[index]->Sumw2();
391 hisname.Form("hLSPt%dLCnw",i);
392 fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
393 fMassHistLS[indexLS]->Sumw2();
394 hisname.Form("hLSPt%dTCnw",i);
395 fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
396 fMassHistLSTC[indexLS]->Sumw2();
400 hisname.Form("hCosPSigPt%dLS",i);
401 fCosPHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
402 fCosPHistLS[index]->Sumw2();
403 hisname.Form("hDLenSigPt%dLS",i);
404 fDLenHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
405 fDLenHistLS[index]->Sumw2();
406 hisname.Form("hSumd02SigPt%dLS",i);
407 fSumd02HistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
408 fSumd02HistLS[index]->Sumw2();
409 hisname.Form("hSigVertSigPt%dLS",i);
410 fSigVertHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
411 fSigVertHistLS[index]->Sumw2();
412 hisname.Form("hPtMaxSigPt%dLS",i);
413 fPtMaxHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
414 fPtMaxHistLS[index]->Sumw2();
416 hisname.Form("hDCASigPt%dLS",i);
417 fDCAHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
418 fDCAHistLS[index]->Sumw2();
422 index=GetBackgroundHistoIndex(i);
424 hisname.Form("hBkgPt%d",i);
425 fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
426 fMassHist[index]->Sumw2();
427 hisname.Form("hCosPBkgPt%d",i);
428 fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
429 fCosPHist[index]->Sumw2();
430 hisname.Form("hDLenBkgPt%d",i);
431 fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
432 fDLenHist[index]->Sumw2();
433 hisname.Form("hSumd02BkgPt%d",i);
434 fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
435 fSumd02Hist[index]->Sumw2();
436 hisname.Form("hSigVertBkgPt%d",i);
437 fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
438 fSigVertHist[index]->Sumw2();
439 hisname.Form("hPtMaxBkgPt%d",i);
440 fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
441 fPtMaxHist[index]->Sumw2();
443 hisname.Form("hDCABkgPt%d",i);
444 fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
445 fDCAHist[index]->Sumw2();
448 hisname.Form("hBkgPt%dTC",i);
449 fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
450 fMassHistTC[index]->Sumw2();
452 hisname.Form("hLSPt%dLCntrip",i);
453 fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
454 fMassHistLS[indexLS]->Sumw2();
455 hisname.Form("hLSPt%dTCntrip",i);
456 fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
457 fMassHistLSTC[indexLS]->Sumw2();
460 hisname.Form("hCosPBkgPt%dLS",i);
461 fCosPHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
462 fCosPHistLS[index]->Sumw2();
463 hisname.Form("hDLenBkgPt%dLS",i);
464 fDLenHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
465 fDLenHistLS[index]->Sumw2();
466 hisname.Form("hSumd02BkgPt%dLS",i);
467 fSumd02HistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
468 fSumd02HistLS[index]->Sumw2();
469 hisname.Form("hSigVertBkgPt%dLS",i);
470 fSigVertHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
471 fSigVertHistLS[index]->Sumw2();
472 hisname.Form("hPtMaxBkgPt%dLS",i);
473 fPtMaxHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
474 fPtMaxHistLS[index]->Sumw2();
475 hisname.Form("hDCABkgPt%dLS",i);
476 fDCAHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
477 fDCAHistLS[index]->Sumw2();
481 hisname.Form("hLSPt%dLCntripsinglecut",i);
482 fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
483 fMassHistLS[indexLS]->Sumw2();
484 hisname.Form("hLSPt%dTCntripsinglecut",i);
485 fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
486 fMassHistLSTC[indexLS]->Sumw2();
489 hisname.Form("hLSPt%dLCspc",i);
490 fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
491 fMassHistLS[indexLS]->Sumw2();
492 hisname.Form("hLSPt%dTCspc",i);
493 fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
494 fMassHistLSTC[indexLS]->Sumw2();
497 fRealProt=new TH1F("fRealProt","fRealProt",100,0, 10);
498 fRealKaon=new TH1F("fRealKaon","fRealKaon",100,0, 10);
499 fWellIDProt3sig=new TH1F("fWellIDProt3sig","fWellIDProt3sig",100,0, 10);
500 fWellIDProt2sig=new TH1F("fWellIDProt2sig","fWellIDProt2sig",100,0, 10);
501 fWellIDProt3sigSe=new TH1F("fWellIDProt3sigSe","fWellIDProt3sigSe",100,0, 10);
502 fWellIDProt2sigSe=new TH1F("fWellIDProt2sigSe","fWellIDProt2sigSe",100,0, 10);
503 fWellIDProt3sigRe=new TH1F("fWellIDProt3sigRe","fWellIDProt3sigRe",100,0, 10);
504 fWellIDProt2sigRe=new TH1F("fWellIDProt2sigRe","fWellIDProt2sigRe",100,0, 10);
505 fWellIDKaon3sig=new TH1F("fWellIDKaon3sig","fWellIDKaon3sig",100,0, 10);
506 fWellIDKaon2sig=new TH1F("fWellIDKaon2sig","fWellIDKaon2sig",100,0, 10);
507 fWellIDKaon3sigSe=new TH1F("fWellIDKaon3sigSe","fWellIDKaon3sigSe",100,0, 10);
508 fWellIDKaon2sigSe=new TH1F("fWellIDKaon2sigSe","fWellIDKaon2sigSe",100,0, 10);
509 fWellIDKaon3sigRe=new TH1F("fWellIDKaon3sigRe","fWellIDKaon3sigRe",100,0, 10);
510 fWellIDKaon2sigRe=new TH1F("fWellIDKaon2sigRe","fWellIDKaon2sigRe",100,0, 10);
511 fTPCSignal3Sigma=new TH2F("fTPCSignal3Sigma","fTPCSignal3Sigma",100,0, 10, 100, 0, 100);
512 fTPCSignal3SigmaReK=new TH2F("fTPCSignal3SigmaReK","fTPCSignal3SigmaReK",100,0, 10, 100, 0, 100);
513 fTPCSignal3SigmaSedK=new TH2F("fTPCSignal3SigmaSedK","fTPCSignal3SigmaSedK",100,0, 10, 100, 0, 100);
514 fTPCSignal3SigmaRep=new TH2F("fTPCSignal3SigmaRep","fTPCSignal3SigmaRep",100,0, 10, 100, 0, 100);
515 fTPCSignal3SigmaSedp=new TH2F("fTPCSignal3SigmaSedp","fTPCSignal3SigmaSedp",100,0, 10, 100, 0, 100);
516 fTPCSignal2Sigma=new TH2F("fTPCSignal2Sigma","fTPCSignal2Sigma",100,0, 10, 100, 0, 100);
517 fTPCSignal2SigmaReK=new TH2F("fTPCSignal2SigmaReK","fTPCSignal2SigmaReK",100,0, 10, 100, 0, 100);
518 fTPCSignal2SigmaSedK=new TH2F("fTPCSignal2SigmaSedK","fTPCSignal2SigmaSedK",100,0, 10, 100, 0, 100);
519 fTPCSignal2SigmaRep=new TH2F("fTPCSignal2SigmaRep","fTPCSignal2SigmaRep",100,0, 10, 100, 0, 100);
520 fTPCSignal2SigmaSedp=new TH2F("fTPCSignal2SigmaSedp","fTPCSignal2SigmaSedp",100,0, 10, 100, 0, 100);
521 fFakeProt3sig=new TH1F("fFakeProt3sig","fFakeProt3sig",100,0, 10);
522 fFakeProt2sig=new TH1F("fFakeProt2sig","fFakeProt2sig",100,0, 10);
523 fFakeProt3sigSe=new TH1F("fFakeProt3sigSe","fFakeProt3sigSe",100,0, 10);
524 fFakeProt2sigSe=new TH1F("fFakeProt2sigSe","fFakeProt2sigSe",100,0, 10);
525 fFakeProt3sigRe=new TH1F("fFakeProt3sigRe","fFakeProt3sigRe",100,0, 10);
526 fFakeProt2sigRe=new TH1F("fFakeProt2sigRe","fFakeProt2sigRe",100,0, 10);
527 fFakeKaon3sig=new TH1F("fFakeKaon3sig","fFakeKaon3sig",100,0, 10);
528 fFakeKaon2sig=new TH1F("fFakeKaon2sig","fFakeKaon2sig",100,0, 10);
529 fFakeKaon3sigSe=new TH1F("fFakeKaon3sigSe","fFakeKaon3sigSe",100,0, 10);
530 fFakeKaon2sigSe=new TH1F("fFakeKaon2sigSe","fFakeKaon2sigSe",100,0, 10);
531 fFakeKaon3sigRe=new TH1F("fFakeKaon3sigRe","fFakeKaon3sigRe",100,0, 10);
532 fFakeKaon2sigRe=new TH1F("fFakeKaon2sigRe","fFakeKaon2sigRe",100,0, 10);
534 for(Int_t i=0; i<3*fNPtBins; i++){
535 fOutput->Add(fMassHist[i]);
536 fOutput->Add(fMassHistTC[i]);
537 fOutput->Add(fCosPHist[i]);
538 fOutput->Add(fDLenHist[i]);
539 fOutput->Add(fSumd02Hist[i]);
540 fOutput->Add(fSigVertHist[i]);
541 fOutput->Add(fPtMaxHist[i]);
542 fOutput->Add(fDCAHist[i]);
545 fOutput->Add(fTPCSignal3Sigma);
546 fOutput->Add(fTPCSignal3SigmaReK);
547 fOutput->Add(fTPCSignal3SigmaSedK);
548 fOutput->Add(fTPCSignal3SigmaRep);
549 fOutput->Add(fTPCSignal3SigmaSedp);
550 fOutput->Add(fTPCSignal2Sigma);
551 fOutput->Add(fTPCSignal2SigmaReK);
552 fOutput->Add(fTPCSignal2SigmaSedK);
553 fOutput->Add(fTPCSignal2SigmaRep);
554 fOutput->Add(fTPCSignal2SigmaSedp);
555 fOutput->Add(fWellIDProt3sig);
556 fOutput->Add(fWellIDProt2sig);
557 fOutput->Add(fWellIDProt3sigSe);
558 fOutput->Add(fWellIDProt2sigSe);
559 fOutput->Add(fWellIDProt3sigRe);
560 fOutput->Add(fWellIDProt2sigRe);
561 fOutput->Add(fWellIDKaon3sig);
562 fOutput->Add(fWellIDKaon2sig);
563 fOutput->Add(fWellIDKaon3sigSe);
564 fOutput->Add(fWellIDKaon2sigSe);
565 fOutput->Add(fWellIDKaon3sigRe);
566 fOutput->Add(fWellIDKaon2sigRe);
568 fOutput->Add(fRealProt);
569 fOutput->Add(fRealKaon);
570 fOutput->Add(fFakeProt3sig);
571 fOutput->Add(fFakeProt2sig);
572 fOutput->Add(fFakeProt3sigSe);
573 fOutput->Add(fFakeProt2sigSe);
574 fOutput->Add(fFakeProt3sigRe);
575 fOutput->Add(fFakeProt2sigRe);
576 fOutput->Add(fFakeKaon3sig);
577 fOutput->Add(fFakeKaon2sig);
578 fOutput->Add(fFakeKaon3sigSe);
579 fOutput->Add(fFakeKaon2sigSe);
580 fOutput->Add(fFakeKaon3sigRe);
581 fOutput->Add(fFakeKaon2sigRe);
582 fHistNEvents = new TH1F("fHistNEvents", "Number of processed events; ; Events",3,-1.5,1.5);
583 fHistNEvents->Sumw2();
584 fHistNEvents->SetMinimum(0);
585 fOutput->Add(fHistNEvents);
590 //OpenFile(3); // 2 is the slot number of the ntuple
592 fNtupleLambdac = new TNtuple("fNtupleLambdac","D +","pdg:Px:Py:Pz:PtTrue:VxTrue:VyTrue:VzTrue:Ptpi:PtK:Ptpi2:PtRec:PointingAngle:DecLeng:VxRec:VyRec:VzRec:InvMass:sigvert:d0Pi:d0K:d0Pi2:dca:d0square");
599 //________________________________________________________________________
600 void AliAnalysisTaskSELambdac::UserExec(Option_t */*option*/)
602 // Execute analysis for current event:
603 // heavy flavor candidates association to MC truth
605 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
607 fHistNEvents->Fill(0); // count event
608 // Post the data already here
611 TClonesArray *array3Prong = 0;
612 TClonesArray *arrayLikeSign =0;
613 if(!aod && AODEvent() && IsStandardAOD()) {
614 // In case there is an AOD handler writing a standard AOD, use the AOD
615 // event in memory rather than the input (ESD) event.
616 aod = dynamic_cast<AliAODEvent*> (AODEvent());
617 // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
618 // have to taken from the AOD event hold by the AliAODExtension
619 AliAODHandler* aodHandler = (AliAODHandler*)
620 ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
621 if(aodHandler->GetExtensions()) {
622 AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
623 AliAODEvent *aodFromExt = ext->GetAOD();
624 array3Prong=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
625 arrayLikeSign=(TClonesArray*)aodFromExt->GetList()->FindObject("LikeSign3Prong");
628 array3Prong=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
629 arrayLikeSign=(TClonesArray*)aod->GetList()->FindObject("LikeSign3Prong");
633 printf("AliAnalysisTaskSELambdac::UserExec: Charm3Prong branch not found!\n");
637 printf("AliAnalysisTaskSELambdac::UserExec: LikeSign3Prong branch not found!\n");
642 TClonesArray *arrayMC=0;
643 AliAODMCHeader *mcHeader=0;
645 // AOD primary vertex
646 AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
651 arrayMC = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
653 printf("AliAnalysisTaskSELambdac::UserExec: MC particles branch not found!\n");
658 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
660 printf("AliAnalysisTaskSELambdac::UserExec: MC header branch not found!\n");
665 Int_t n3Prong = array3Prong->GetEntriesFast();
670 for (Int_t i3Prong = 0; i3Prong < n3Prong; i3Prong++) {
671 AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)array3Prong->UncheckedAt(i3Prong);
674 Bool_t unsetvtx=kFALSE;
675 if(!d->GetOwnPrimaryVtx()){
676 d->SetOwnPrimaryVtx(vtx1);
680 if(fRDCutsProduction->IsSelected(d,AliRDHFCuts::kCandidate)) {
682 Double_t ptCand = d->Pt();
684 for(Int_t ibin=0;ibin<fNPtBins&&iPtBin<0&&ptCand>fArrayBinLimits[0]&&ptCand<fArrayBinLimits[fNPtBins];ibin++){
685 if(ptCand<fArrayBinLimits[ibin+1])iPtBin=ibin;
688 Bool_t passTightCuts=fRDCutsAnalysis->IsSelected(d,AliRDHFCuts::kCandidate);
699 labDp = MatchToMCLambdac(d,arrayMC);
702 AliAODMCParticle *partDp = (AliAODMCParticle*)arrayMC->At(labDp);
703 AliAODMCParticle *dg0 = (AliAODMCParticle*)arrayMC->At(partDp->GetDaughter(0));
704 deltaPx=partDp->Px()-d->Px();
705 deltaPy=partDp->Py()-d->Py();
706 deltaPz=partDp->Pz()-d->Pz();
711 pdgCode=TMath::Abs(partDp->GetPdgCode());
717 Double_t invMasspKpi=-1.;
718 Double_t invMasspiKp=-1.;
719 Int_t pdgs[3]={0,0,0};
720 Double_t field=aod->GetMagneticField();
721 if(fReadMC && fMCPid){
723 if(IspKpiMC(d,arrayMC,pdgs)) {
724 invMasspKpi=d->InvMassLcpKpi();
726 pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
727 if(!VertexingKF(d,pdgs,field)) invMasspKpi=-1.;
730 if(IspiKpMC(d,arrayMC)) {
731 invMasspiKp=d->InvMassLcpiKp();
733 pdgs[0]=211;pdgs[1]=321;pdgs[2]=2212;
734 if(!VertexingKF(d,pdgs,field)) invMasspiKp=-1.;
740 invMasspKpi=d->InvMassLcpKpi();
742 pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
743 if(!VertexingKF(d,pdgs,field)) invMasspKpi=-1.;
746 if(fReadMC) IspKpiMC(d,arrayMC,pdgs);
747 if(IspiKpReal(d,pdgs)) {
748 invMasspiKp=d->InvMassLcpiKp();
750 pdgs[0]=211;pdgs[1]=321;pdgs[2]=2212;
751 if(!VertexingKF(d,pdgs,field)) invMasspiKp=-1.;
756 if(IspKpiResonant(d,field)) {
757 invMasspKpi=d->InvMassLcpKpi();
759 pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
760 if(!VertexingKF(d,pdgs,field)) invMasspKpi=-1.;
763 if(IspiKpResonant(d,field)) {
764 invMasspiKp=d->InvMassLcpiKp();
766 pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
767 if(!VertexingKF(d,pdgs,field)) invMasspiKp=-1.;
771 if(!fResPid && !fRealPid && !fMCPid){
772 invMasspiKp=d->InvMassLcpiKp();
774 pdgs[0]=211;pdgs[1]=321;pdgs[2]=2212;
775 if(!VertexingKF(d,pdgs,field)) invMasspiKp=-1.;
777 invMasspKpi=d->InvMassLcpKpi();
779 pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
780 if(!VertexingKF(d,pdgs,field)) invMasspKpi=-1.;
784 if(invMasspiKp<0. && invMasspKpi<0.) continue;
797 tmp[8]=d->PtProng(0);
798 tmp[9]=d->PtProng(1);
799 tmp[10]=d->PtProng(2);
801 tmp[12]=d->CosPointingAngle();
802 tmp[13]=d->DecayLength();
806 if(invMasspiKp>0.) tmp[17]=invMasspiKp;
807 if(invMasspKpi>0.) tmp[17]=invMasspKpi;
808 tmp[18]=d->GetSigmaVert();
809 tmp[19]=d->Getd0Prong(0);
810 tmp[20]=d->Getd0Prong(1);
811 tmp[21]=d->Getd0Prong(2);
813 tmp[23]=d->Prodd0d0();
814 fNtupleLambdac->Fill(tmp);
815 PostData(3,fNtupleLambdac);
817 Double_t dlen=d->DecayLength();
818 Double_t cosp=d->CosPointingAngle();
819 Double_t sumD02=d->Getd0Prong(0)*d->Getd0Prong(0)+d->Getd0Prong(1)*d->Getd0Prong(1)+d->Getd0Prong(2)*d->Getd0Prong(2);
820 Double_t dca=d->GetDCA();
822 for(Int_t i=0;i<3;i++){
823 if(d->PtProng(i)>ptmax)ptmax=d->PtProng(i);
827 index=GetHistoIndex(iPtBin);
828 if(invMasspiKp>0. && invMasspKpi>0.){
829 if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp,0.5);
830 if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi,0.5);
832 if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp);
833 if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi);
835 fCosPHist[index]->Fill(cosp);
836 fDLenHist[index]->Fill(dlen);
837 fSumd02Hist[index]->Fill(sumD02);
838 fPtMaxHist[index]->Fill(ptmax);
839 fDCAHist[index]->Fill(dca);
842 if(invMasspiKp>0. && invMasspKpi>0.){
843 if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp,0.5);
844 if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi,0.5);
846 if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp);
847 if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi);
853 index=GetSignalHistoIndex(iPtBin);
854 if(invMasspiKp>0. && invMasspKpi>0.){
855 if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp,0.5);
856 if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi,0.5);
858 if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp);
859 if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi);
861 fCosPHist[index]->Fill(cosp);
862 fDLenHist[index]->Fill(dlen);
863 fSumd02Hist[index]->Fill(sumD02);
864 fPtMaxHist[index]->Fill(ptmax);
865 fDCAHist[index]->Fill(dca);
867 if(invMasspiKp>0. && invMasspKpi>0.){
868 if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp,0.5);
869 if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi,0.5);
871 if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp);
872 if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi);
877 index=GetBackgroundHistoIndex(iPtBin);
878 if(invMasspiKp>0. && invMasspKpi>0.){
879 fMassHist[index]->Fill(invMasspiKp,0.5);
880 fMassHist[index]->Fill(invMasspKpi,0.5);
882 if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp);
883 if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi);
885 fCosPHist[index]->Fill(cosp);
886 fDLenHist[index]->Fill(dlen);
887 fSumd02Hist[index]->Fill(sumD02);
888 fPtMaxHist[index]->Fill(ptmax);
889 fDCAHist[index]->Fill(dca);
891 if(invMasspiKp>0. && invMasspKpi>0.){
892 fMassHistTC[index]->Fill(invMasspiKp,0.5);
893 fMassHistTC[index]->Fill(invMasspKpi,0.5);
895 if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp);
896 if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi);
904 if(labDp<0)fHistOSbkg->Fill(d->InvMassDplus());
905 fHistOS->Fill(d->InvMassDplus());
909 if(unsetvtx) d->UnsetOwnPrimaryVtx();
918 //________________________________________________________________________
919 void AliAnalysisTaskSELambdac::Terminate(Option_t */*option*/)
921 // Terminate analysis
923 if(fDebug > 1) printf("AnalysisTaskSELambdac: Terminate() \n");
925 fOutput = dynamic_cast<TList*> (GetOutputData(1));
927 printf("ERROR: fOutput not available\n");
930 fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));
937 for(Int_t i=0;i<fNPtBins;i++){
938 index=GetHistoIndex(i);
939 hisname.Form("hMassPt%d",i);
940 fMassHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
941 hisname.Form("hCosPAllPt%d",i);
942 fCosPHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
943 hisname.Form("hDLenAllPt%d",i);
944 fDLenHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
945 hisname.Form("hSumd02AllPt%d",i);
946 fSumd02Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
947 hisname.Form("hSigVertAllPt%d",i);
948 fSigVertHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
949 hisname.Form("hPtMaxAllPt%d",i);
950 fPtMaxHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
951 hisname.Form("hDCAAllPt%d",i);
952 fDCAHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
953 hisname.Form("hMassPt%dTC",i);
954 fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
956 index=GetSignalHistoIndex(i);
957 hisname.Form("hSigPt%d",i);
958 fMassHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
959 hisname.Form("hCosPSigPt%d",i);
960 fCosPHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
961 hisname.Form("hDLenSigPt%d",i);
962 fDLenHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
963 hisname.Form("hSumd02SigPt%d",i);
964 fSumd02Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
965 hisname.Form("hSigVertSigPt%d",i);
966 fSigVertHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
967 hisname.Form("hPtMaxSigPt%d",i);
968 fPtMaxHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
969 hisname.Form("hDCASigPt%d",i);
970 fDCAHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
972 hisname.Form("hSigPt%dTC",i);
973 fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
975 index=GetBackgroundHistoIndex(i);
976 hisname.Form("hBkgPt%d",i);
977 fMassHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
978 hisname.Form("hCosPBkgPt%d",i);
979 fCosPHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
980 hisname.Form("hDLenBkgPt%d",i);
981 fDLenHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
982 hisname.Form("hSumd02BkgPt%d",i);
983 fSumd02Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
984 hisname.Form("hSigVertBkgPt%d",i);
985 fSigVertHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
986 hisname.Form("hPtMaxBkgPt%d",i);
987 fPtMaxHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
988 hisname.Form("hDCABkgPt%d",i);
989 fDCAHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
990 hisname.Form("hBkgPt%dTC",i);
991 fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
994 fWellIDProt3sig=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDProt3sig));
995 fWellIDProt2sig=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDProt2sig));
996 fWellIDProt3sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDProt3sigSe));
997 fWellIDProt2sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDProt2sigSe));
998 fWellIDProt2sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDProt2sigRe));
999 fWellIDProt3sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDProt3sigRe));
1000 fWellIDKaon3sig=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDKaon3sig));
1001 fWellIDKaon2sig=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDKaon2sig));
1002 fWellIDKaon3sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDKaon3sigSe));
1003 fWellIDKaon2sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDKaon2sigSe));
1004 fWellIDKaon2sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDKaon2sigRe));
1005 fWellIDKaon3sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDKaon3sigRe));
1006 fTPCSignal3Sigma=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3Sigma));
1007 fTPCSignal3SigmaReK=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3SigmaReK));
1008 fTPCSignal3SigmaSedK=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3SigmaSedK));
1009 fTPCSignal3SigmaRep=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3SigmaRep));
1010 fTPCSignal3SigmaSedp=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3SigmaSedp));
1011 fTPCSignal2Sigma=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2Sigma));
1012 fTPCSignal2SigmaReK=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2SigmaReK));
1013 fTPCSignal2SigmaSedK=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2SigmaSedK));
1014 fTPCSignal2SigmaSedp=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2SigmaSedp));
1015 fTPCSignal2SigmaRep=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2SigmaRep));
1016 fFakeProt3sig=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeProt3sig));
1017 fFakeProt2sig=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeProt2sig));
1018 fFakeProt3sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeProt3sigSe));
1019 fFakeProt2sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeProt2sigSe));
1020 fFakeProt2sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeProt2sigRe));
1021 fFakeProt3sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeProt3sigRe));
1022 fFakeKaon3sig=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeKaon3sig));
1023 fFakeKaon2sig=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeKaon2sig));
1024 fFakeKaon3sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeKaon3sigSe));
1025 fFakeKaon2sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeKaon2sigSe));
1026 fFakeKaon2sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeKaon2sigRe));
1027 fFakeKaon3sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeKaon3sigRe));
1029 fRealProt=dynamic_cast<TH1F*>(fOutput->FindObject(fRealProt));
1030 fRealKaon=dynamic_cast<TH1F*>(fOutput->FindObject(fRealKaon));
1033 fNtupleLambdac = dynamic_cast<TNtuple*>(GetOutputData(3));
1036 TCanvas *c1=new TCanvas("c1","D+ invariant mass distribution",500,500);
1038 //TH1F *hMassPt=(TH1F*)fOutput->FindObject("hMassPt3TC");
1044 //________________________________________________________________________
1045 Int_t AliAnalysisTaskSELambdac::MatchToMCLambdac(AliAODRecoDecayHF3Prong *d,TClonesArray *arrayMC){
1047 Int_t lambdacLab[3]={0,0,0};
1048 Int_t pdgs[3]={0,0,0};
1049 for(Int_t i=0;i<3;i++){
1050 AliAODTrack *daugh=(AliAODTrack*)d->GetDaughter(i);
1051 Int_t lab=daugh->GetLabel();
1053 AliAODMCParticle *part= (AliAODMCParticle*)arrayMC->At(lab);
1055 pdgs[i]=part->GetPdgCode();
1056 Int_t partPdgcode = TMath::Abs(part->GetPdgCode());
1057 if(partPdgcode==211 || partPdgcode==321 || partPdgcode==2212){
1058 Int_t motherLabel=part->GetMother();
1059 if(motherLabel<0) return 0;
1060 AliAODMCParticle *motherPart = (AliAODMCParticle*)arrayMC->At(motherLabel);
1061 if(!motherPart) continue;
1062 Int_t motherPdg = TMath::Abs(motherPart->GetPdgCode());
1063 if(motherPdg==4122) {
1064 if(GetLambdacDaugh(motherPart,arrayMC)){lambdacLab[i]=motherLabel;continue;}
1066 if(motherPdg==313 || motherPdg==2224 || motherPdg==3124){
1067 Int_t GmotherLabel=motherPart->GetMother();
1068 if(GmotherLabel<0) return 0;
1069 AliAODMCParticle *GmotherPart = (AliAODMCParticle*)arrayMC->At(GmotherLabel);
1070 if(!GmotherPart) continue;
1071 Int_t GmotherPdg = TMath::Abs(GmotherPart->GetPdgCode());
1072 if(GmotherPdg==4122) {
1073 if(GetLambdacDaugh(GmotherPart,arrayMC)) {lambdacLab[i]=GmotherLabel;continue;}
1079 if(lambdacLab[0]==lambdacLab[1] && lambdacLab[1]==lambdacLab[2]) {return lambdacLab[0];}
1083 //------------------------
1084 Bool_t AliAnalysisTaskSELambdac::GetLambdacDaugh(AliAODMCParticle *part,TClonesArray *arrayMC){
1086 Int_t numberOfLambdac=0;
1087 if(TMath::Abs(part->GetPdgCode())!=4122) return kFALSE;
1089 daugh_tmp[0]=part->GetDaughter(0);
1090 daugh_tmp[1]=part->GetDaughter(1);
1091 Int_t nDaugh = (Int_t)part->GetNDaughters();
1092 if(nDaugh<2) return kFALSE;
1093 if(nDaugh>3) return kFALSE;
1094 AliAODMCParticle* pdaugh1 = (AliAODMCParticle*)arrayMC->At(part->GetDaughter(0));
1095 if(!pdaugh1) {return kFALSE;}
1096 Int_t number1 = TMath::Abs(pdaugh1->GetPdgCode());
1097 AliAODMCParticle* pdaugh2 = (AliAODMCParticle*)arrayMC->At(part->GetDaughter(1));
1098 if(!pdaugh2) {return kFALSE;}
1099 Int_t number2 = TMath::Abs(pdaugh2->GetPdgCode());
1102 Int_t thirdDaugh=part->GetDaughter(1)-1;
1103 AliAODMCParticle* pdaugh3 = (AliAODMCParticle*)arrayMC->At(thirdDaugh);
1104 Int_t number3 = TMath::Abs(pdaugh3->GetPdgCode());
1105 if((number1==321 && number2==211 && number3==2212) || (number1==211 && number2==321 && number3==2212) || (number1==211 && number2==2212 && number3==321) || (number1==321 && number2==2212 && number3==211) || (number1==2212 && number2==321 && number3==211) || (number1==2212 && number2==211 && number3==321)) {
1112 //Lambda in Lambda pi
1113 // if((number1==211 && number2==3122)|| (number1==3122 && number2==211))return kFALSE;
1121 if((number1==2212 && number2==313)){
1122 nfiglieK=pdaugh2->GetNDaughters();
1123 if(nfiglieK!=2) return kFALSE;
1124 AliAODMCParticle* pdaughK1 = (AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(0));
1125 AliAODMCParticle* pdaughK2 = (AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(1));
1126 if(!pdaughK1) return kFALSE;
1127 if(!pdaughK2) return kFALSE;
1128 if((TMath::Abs(pdaughK1->GetPdgCode())==211 && TMath::Abs(pdaughK2->GetPdgCode())==321) || (TMath::Abs(pdaughK1->GetPdgCode())==321 && TMath::Abs(pdaughK2->GetPdgCode())==211)) numberOfLambdac++;
1131 if((number1==313 && number2==2212)){
1132 nfiglieK=pdaugh1->GetNDaughters();
1133 if(nfiglieK!=2) return kFALSE;
1134 AliAODMCParticle* pdaughK1 = (AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(0));
1135 AliAODMCParticle* pdaughK2 = (AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(1));
1136 if(!pdaughK1) return kFALSE;
1137 if(!pdaughK2) return kFALSE;
1138 if((TMath::Abs(pdaughK1->GetPdgCode())==211 && TMath::Abs(pdaughK2->GetPdgCode())==321) || (TMath::Abs(pdaughK1->GetPdgCode())==321 && TMath::Abs(pdaughK2->GetPdgCode())==211)) numberOfLambdac++;
1141 //Lambda -> Delta++ k
1142 Int_t nfiglieDelta=0;
1143 if(number1==321 && number2==2224){
1144 nfiglieDelta=pdaugh2->GetNDaughters();
1145 if(nfiglieDelta!=2) return kFALSE;
1146 AliAODMCParticle *pdaughD1=(AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(0));
1147 AliAODMCParticle *pdaughD2=(AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(1));
1148 if(!pdaughD1) return kFALSE;
1149 if(!pdaughD2) return kFALSE;
1150 if((TMath::Abs(pdaughD1->GetPdgCode())==211 && TMath::Abs(pdaughD2->GetPdgCode())==2212) || (TMath::Abs(pdaughD1->GetPdgCode())==2212 && TMath::Abs(pdaughD2->GetPdgCode())==211)) numberOfLambdac++;
1152 if(number1==2224 && number2==321){
1153 nfiglieDelta=pdaugh1->GetNDaughters();
1154 if(nfiglieDelta!=2) return kFALSE;
1155 AliAODMCParticle* pdaughD1 = (AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(0));
1156 AliAODMCParticle* pdaughD2 = (AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(1));
1157 if(!pdaughD1) return kFALSE;
1158 if(!pdaughD2) return kFALSE;
1159 if((TMath::Abs(pdaughD1->GetPdgCode())==211 && TMath::Abs(pdaughD2->GetPdgCode())==2212) || (TMath::Abs(pdaughD1->GetPdgCode())==2212 && TMath::Abs(pdaughD2->GetPdgCode())==211)) numberOfLambdac++;
1163 //Lambdac -> Lambda(1520) pi
1165 if(number1==3124 && number2==211){
1166 nfiglieLa=pdaugh1->GetNDaughters();
1167 if(nfiglieLa!=2) return kFALSE;
1168 AliAODMCParticle *pdaughL1=(AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(0));
1169 AliAODMCParticle *pdaughL2=(AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(1));
1170 if(!pdaughL1) return kFALSE;
1171 if(!pdaughL2) return kFALSE;
1172 if((TMath::Abs(pdaughL1->GetPdgCode())==321 && TMath::Abs(pdaughL2->GetPdgCode())==2212) || (TMath::Abs(pdaughL1->GetPdgCode())==2212 && TMath::Abs(pdaughL2->GetPdgCode())==321)) numberOfLambdac++;
1174 if(number1==211 && number2==3124){
1175 nfiglieLa=pdaugh2->GetNDaughters();
1176 if(nfiglieLa!=2) return kFALSE;
1177 AliAODMCParticle *pdaughL1=(AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(0));
1178 AliAODMCParticle *pdaughL2=(AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(1));
1179 if(!pdaughL1) return kFALSE;
1180 if(!pdaughL2) return kFALSE;
1181 if((TMath::Abs(pdaughL1->GetPdgCode())==321 && TMath::Abs(pdaughL2->GetPdgCode())==2212) || (TMath::Abs(pdaughL1->GetPdgCode())==2212 && TMath::Abs(pdaughL2->GetPdgCode())==321)) numberOfLambdac++;
1186 if(numberOfLambdac>0) {return kTRUE;}
1189 //-----------------------------
1190 Bool_t AliAnalysisTaskSELambdac::IspKpiMC(AliAODRecoDecayHF3Prong *d,TClonesArray *arrayMC,Int_t *pdgs){
1192 Int_t lab[3]={0,0,0};//,pdgs[3]={0,0,0};
1193 for(Int_t i=0;i<3;i++){
1194 AliAODTrack *daugh=(AliAODTrack*)d->GetDaughter(i);
1195 lab[i]=daugh->GetLabel();
1196 if(lab[i]<0) return kFALSE;
1197 AliAODMCParticle *part= (AliAODMCParticle*)arrayMC->At(lab[i]);
1198 if(!part) return kFALSE;
1199 pdgs[i]=TMath::Abs(part->GetPdgCode());
1202 if(pdgs[0]==2212 && pdgs[1]==321 && pdgs[2]==211) return kTRUE;
1206 //-----------------------------
1207 Bool_t AliAnalysisTaskSELambdac::IspiKpMC(AliAODRecoDecayHF3Prong *d,TClonesArray *arrayMC){
1209 Int_t lab[3]={0,0,0},pdgs[3]={0,0,0};
1210 for(Int_t i=0;i<3;i++){
1211 AliAODTrack *daugh=(AliAODTrack*)d->GetDaughter(i);
1212 lab[i]=daugh->GetLabel();
1213 if(lab[i]<0) return kFALSE;
1214 AliAODMCParticle *part= (AliAODMCParticle*)arrayMC->At(lab[i]);
1215 if(!part) return kFALSE;
1216 pdgs[i]=TMath::Abs(part->GetPdgCode());
1219 if(pdgs[2]==2212 && pdgs[1]==321 && pdgs[0]==211) {return kTRUE;}
1223 //--------------------------------------
1224 Bool_t AliAnalysisTaskSELambdac::IspiKpReal(AliAODRecoDecayHF3Prong *d,Int_t *pdgs){
1226 AliAODPidHF* pidObj=new AliAODPidHF();
1227 Bool_t type[2]={kFALSE,kFALSE};
1228 for(Int_t i=0;i<3;i++){
1229 // Bool_t pid[3]={kFALSE,kFALSE,kFALSE};
1230 AliAODTrack *track=(AliAODTrack*)d->GetDaughter(i);
1231 //test asymmetric TPC PID
1232 Double_t plim[2]={0.6,2.};
1233 Double_t sigmas[5]={2.,1.,2.,3.,0.};
1234 pidObj->SetPLimit(plim);
1235 pidObj->SetSigma(sigmas);
1236 if(i==1) type[i]=pidObj->MatchTPCTOF(track,2,3,kFALSE);
1237 if(i==2) type[i]=pidObj->MatchTPCTOF(track,2,4,kFALSE);
1238 //pidinfo=pidObj->MatchTPCTOF(track,3,3,kFALSE);
1242 if(type[0] && type[1]) {delete pidObj;return kTRUE;}
1246 //--------------------------------------
1247 Bool_t AliAnalysisTaskSELambdac::IspKpiReal(AliAODRecoDecayHF3Prong *d){
1249 AliAODPidHF* pidObj=new AliAODPidHF();
1250 Bool_t type[2]={kFALSE,kFALSE};
1251 for(Int_t i=0;i<2;i++){
1252 //Bool_t pid[3]={kFALSE,kFALSE,kFALSE};
1253 AliAODTrack *track=(AliAODTrack*)d->GetDaughter(i);
1254 Double_t plim[2]={0.6,2.};
1255 Double_t sigmas[5]={2.,1.,2.,3.,0.};
1256 pidObj->SetPLimit(plim);
1257 pidObj->SetSigma(sigmas);
1258 if(i==1) type[i]=pidObj->MatchTPCTOF(track,2,3,kFALSE);
1259 if(i==0) type[i]=pidObj->MatchTPCTOF(track,2,4,kFALSE);
1261 if(type[0] && type[1]) {delete pidObj;return kTRUE;}
1265 //------------------------------------
1266 Bool_t AliAnalysisTaskSELambdac::VertexingKF(AliAODRecoDecayHF3Prong *d,Int_t *pdgs,Double_t field){
1268 Int_t iprongs[3]={0,1,2};
1269 Double_t mass[2]={0.,0.};
1270 //topological constr
1271 AliKFParticle *Lambdac=d->ApplyVertexingKF(iprongs,3,pdgs,kTRUE,field,mass);
1272 if(!Lambdac) return kFALSE;
1273 Double_t probTot=TMath::Prob(Lambdac->GetChi2(),Lambdac->GetNDF());
1274 if(probTot<fCutsKF[0]) return kFALSE;
1276 if(d->PtProng(0)<1.5 && d->Getd0Prong(0)>0.02) return kFALSE;
1277 if(d->PtProng(1)<1.5 && d->Getd0Prong(1)>0.02) return kFALSE;
1278 if(d->PtProng(2)<1.5 && d->Getd0Prong(2)>0.02) return kFALSE;
1279 //mass constr for K*
1282 mass[0]=0.8961;mass[1]=0.03;
1283 if(TMath::Abs(pdgs[0])==211){
1284 ipRes[0]=0;ipRes[1]=1;
1285 pdgres[0]=pdgs[0];pdgres[1]=321;
1287 if(TMath::Abs(pdgs[2])==211){
1288 ipRes[0]=2;ipRes[1]=1;
1289 pdgres[0]=pdgs[2];pdgres[1]=321;
1291 AliKFParticle *KappaStar=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1293 Double_t probKstar=TMath::Prob(KappaStar->GetChi2(),KappaStar->GetNDF());
1294 if(probKstar>fCutsKF[1]) {
1295 AliAODTrack *esdProng1=(AliAODTrack*)d->GetDaughter(ipRes[0]);
1296 AliAODTrack *esdProng2=(AliAODTrack*)d->GetDaughter(ipRes[1]);
1297 AliKFParticle prong1(*esdProng1,pdgres[0]);
1298 AliKFParticle prong2(*esdProng2,pdgres[1]);
1299 if(KappaStar->GetPt()<fCutsKF[2] && prong1.GetAngle(prong2)>fCutsKF[3]) return kFALSE;
1301 //mass constr for Lambda
1302 mass[0]=1.520;mass[1]=0.005;
1303 if(TMath::Abs(pdgs[0])==2212){
1304 ipRes[0]=0;ipRes[1]=1;
1305 pdgres[0]=pdgs[0];pdgres[1]=pdgs[1];
1307 if(TMath::Abs(pdgs[2])==2212){
1308 ipRes[0]=2;ipRes[1]=1;
1309 pdgres[0]=pdgs[2];pdgres[1]=pdgs[1];
1311 AliKFParticle *Lambda1520=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1312 Double_t probLa=TMath::Prob(Lambda1520->GetChi2(),Lambda1520->GetNDF());
1313 if(probLa>fCutsKF[4]) {
1314 AliAODTrack *esdProng1=(AliAODTrack*)d->GetDaughter(ipRes[0]);
1315 AliAODTrack *esdProng2=(AliAODTrack*)d->GetDaughter(ipRes[1]);
1316 AliKFParticle prong1(*esdProng1,pdgres[0]);
1317 AliKFParticle prong2(*esdProng2,pdgres[1]);
1318 if(Lambda1520->GetPt()<fCutsKF[5] && prong1.GetAngle(prong2)>fCutsKF[6]) return kFALSE;
1320 //mass constr for Delta
1321 mass[0]=1.2;mass[1]=0.15;
1322 ipRes[0]=0;ipRes[1]=2;
1323 pdgres[0]=pdgs[0];pdgres[2]=pdgs[2];
1324 AliKFParticle *Delta=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1325 Double_t probDelta=TMath::Prob(Delta->GetChi2(),Delta->GetNDF());
1326 if(probDelta>fCutsKF[7]) {
1327 AliAODTrack *esdProng1=(AliAODTrack*)d->GetDaughter(ipRes[0]);
1328 AliAODTrack *esdProng2=(AliAODTrack*)d->GetDaughter(ipRes[1]);
1329 AliKFParticle prong1(*esdProng1,pdgres[0]);
1330 AliKFParticle prong2(*esdProng2,pdgres[1]);
1331 if(Delta->GetPt()<fCutsKF[8] && prong1.GetAngle(prong2)>fCutsKF[9]) return kFALSE;
1335 //-------------------------------------
1336 Bool_t AliAnalysisTaskSELambdac::IspiKpResonant(AliAODRecoDecayHF3Prong *d,Double_t field){
1339 Double_t mass[2]={1.520,0.005};
1340 Int_t ipRes[2]={1,2};
1341 Int_t pdgres[2]={321,2212};
1342 AliKFParticle *Lambda1520=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1343 Double_t probLa=TMath::Prob(Lambda1520->GetChi2(),Lambda1520->GetNDF());
1344 if(probLa>0.1) return kTRUE;
1346 // mass[0]=0.8961;mass[1]=0.03;
1347 // ipRes[0]=0;ipRes[1]=1;
1348 // pdgres[0]=211;pdgres[1]=321;
1349 // AliKFParticle *Kstar=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1350 // Double_t probKa=TMath::Prob(Kstar->GetChi2(),Kstar->GetNDF());
1351 // if(probKa>0.1) return kTRUE;
1356 //-------------------------------------
1357 Bool_t AliAnalysisTaskSELambdac::IspKpiResonant(AliAODRecoDecayHF3Prong *d,Double_t field){
1360 Double_t mass[2]={1.520,0.005};
1361 Int_t ipRes[2]={0,1};
1362 Int_t pdgres[2]={2212,321};
1363 AliKFParticle *Lambda1520=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1364 Double_t probLa=TMath::Prob(Lambda1520->GetChi2(),Lambda1520->GetNDF());
1365 if(probLa>0.1) return kTRUE;
1367 // mass[0]=0.8961;mass[1]=0.03;
1368 // ipRes[0]=1;ipRes[1]=2;
1369 // pdgres[1]=211;pdgres[0]=321;
1370 // AliKFParticle *Kstar=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1371 // Double_t probKa=TMath::Prob(Kstar->GetChi2(),Kstar->GetNDF());
1372 // if(probKa>0.1) return kTRUE;