1 #include "AliAnalysisTaskSE.h"
2 #include "AliAnalysisManager.h"
3 #include "AliAnalysisDataContainer.h"
4 #include "AliESDEvent.h"
5 #include "AliESDtrackCuts.h"
7 #include "AliCentrality.h"
8 #include "AliMCEventHandler.h"
9 #include "AliMCEvent.h"
10 #include "AliHeader.h"
11 #include "AliGenCocktailEventHeader.h"
12 #include "AliGenEventHeader.h"
13 #include "AliGenerator.h"
15 #include "AliMultiplicity.h"
16 #include <TParticle.h>
24 #include "AliESDInputHandlerRP.h"
65 fHistoTrackletsEta1(0),
80 fHistOriginFeeddown(0),
86 fHistEtaPhiPtGenEle(0),
87 fHistEtaPhiPtGenPi(0),
89 fHistEtaPhiPtGenPro(0),
90 fHistEtaPhiPtRecEle(0),
91 fHistEtaPhiPtRecPi(0),
93 fHistEtaPhiPtRecPro(0),
94 fSearchUpToQuark(kFALSE),
100 for(Int_t i=0; i<5; i++){
108 for(Int_t i=0; i<2; i++){
115 DefineInput(0, TChain::Class());
116 DefineOutput(1, TList::Class());
123 if (AliAnalysisManager::GetAnalysisManager()->IsProofMode())
return;
138 fOutput->SetName(
"OutputHistos");
140 fHistoNEvents =
new TH1F(
"hNEvents",
"Number of processed events",3,-0.5,2.5);
151 fHistoTracks =
new TH1F(
"hTracks",
"",100,-0.5,maxMult*2-0.5);
204 Double_t maxncn=nBinscb-0.5;
205 fHistoNcharmed =
new TH2F(
"hncharmed",
"",100,-0.5,maxMult-0.5,nBinscb,-0.5,maxncn);
208 fHistoNbVsNc =
new TH2F(
"hnbvsnc",
"",nBinscb,-0.5,maxncn,nBinscb,-0.5,maxncn);
212 fHistYPtPrompt[0] =
new TH2F(
"hyptD0prompt",
"D0 - Prompt",40,0.,40.,20,-2.,2.);
213 fHistYPtPrompt[1] =
new TH2F(
"hyptDplusprompt",
"Dplus - Prompt",40,0.,40.,20,-2.,2.);
214 fHistYPtPrompt[2] =
new TH2F(
"hyptDstarprompt",
"Dstar - Prompt",40,0.,40.,20,-2.,2.);
215 fHistYPtPrompt[3] =
new TH2F(
"hyptDsprompt",
"Ds - Prompt",40,0.,40.,20,-2.,2.);
216 fHistYPtPrompt[4] =
new TH2F(
"hyptLcprompt",
"Lc - Prompt",40,0.,40.,20,-2.,2.);
218 fHistBYPtAllDecay[0] =
new TH2F(
"hyptB0AllDecay",
"B0 - All",40,0.,40.,40,-2.,2.);
219 fHistBYPtAllDecay[1] =
new TH2F(
"hyptBplusAllDecay",
"Bplus - All",40,0.,40.,40,-2.,2.);
220 fHistBYPtAllDecay[2] =
new TH2F(
"hyptBstarAllDecay",
"Bstar - All",40,0.,40.,40,-2.,2.);
221 fHistBYPtAllDecay[3] =
new TH2F(
"hyptBsAllDecay",
"Bs - All",40,0.,40.,40,-2.,2.);
222 fHistBYPtAllDecay[4] =
new TH2F(
"hyptLbAllDecay",
"LB - All",40,0.,40.,40,-2.,2.);
224 fHistYPtAllDecay[0] =
new TH2F(
"hyptD0AllDecay",
"D0 - All",40,0.,40.,40,-2.,2.);
225 fHistYPtAllDecay[1] =
new TH2F(
"hyptDplusAllDecay",
"Dplus - All",40,0.,40.,40,-2.,2.);
226 fHistYPtAllDecay[2] =
new TH2F(
"hyptDstarAllDecay",
"Dstar - All",40,0.,40.,40,-2.,2.);
227 fHistYPtAllDecay[3] =
new TH2F(
"hyptDsAllDecay",
"Ds - All",40,0.,40.,40,-2.,2.);
228 fHistYPtAllDecay[4] =
new TH2F(
"hyptLcAllDecay",
"Lc - All",40,0.,40.,40,-2.,2.);
243 fHistYPtFeeddown[0] =
new TH2F(
"hyptD0feeddown",
"D0 - Feeddown",40,0.,40.,20,-2.,2.);
244 fHistYPtFeeddown[1] =
new TH2F(
"hyptDplusfeeddown",
"Dplus - Feeddown",40,0.,40.,20,-2.,2.);
245 fHistYPtFeeddown[2] =
new TH2F(
"hyptDstarfeedown",
"Dstar - Feeddown",40,0.,40.,20,-2.,2.);
246 fHistYPtFeeddown[3] =
new TH2F(
"hyptDsfeedown",
"Ds - Feeddown",40,0.,40.,20,-2.,2.);
247 fHistYPtFeeddown[4] =
new TH2F(
"hyptLcfeedown",
"Lc - Feeddown",40,0.,40.,20,-2.,2.);
249 for(Int_t ih=0; ih<5; ih++){
278 for(Int_t ih=0; ih<2; ih++){
340 fHistLcDecayChan->GetXaxis()->SetBinLabel(9,
"#pi#Lambda#rightarrowp#pi#pi");
347 Double_t binseta[11]={-1.0,-0.8,-0.6,-0.4,-0.2,0.,0.2,0.4,0.6,0.8,1.0};
348 const Int_t nBinsPhi=40;
349 Double_t binsphi[nBinsPhi+1];
350 for(Int_t ib=0; ib<=nBinsPhi; ib++) binsphi[ib]=ib*TMath::Pi()/20.;
351 const Int_t nBinsPt=24;
352 Double_t binspt[nBinsPt+1]={0.,0.10,0.15,0.2,0.25,
358 fHistEtaPhiPtGenEle=
new TH3F(
"hEtaPhiPtGenEle",
"",10,binseta,nBinsPhi,binsphi,nBinsPt,binspt);
360 fHistEtaPhiPtGenPi=
new TH3F(
"hEtaPhiPtGenPi",
"",10,binseta,nBinsPhi,binsphi,nBinsPt,binspt);
362 fHistEtaPhiPtGenK=
new TH3F(
"hEtaPhiPtGenK",
"",10,binseta,nBinsPhi,binsphi,nBinsPt,binspt);
364 fHistEtaPhiPtGenPro=
new TH3F(
"hEtaPhiPtGenPro",
"",10,binseta,nBinsPhi,binsphi,nBinsPt,binspt);
368 fHistEtaPhiPtRecEle=
new TH3F(
"hEtaPhiPtRecEle",
"",10,binseta,nBinsPhi,binsphi,nBinsPt,binspt);
370 fHistEtaPhiPtRecPi=
new TH3F(
"hEtaPhiPtRecPi",
"",10,binseta,nBinsPhi,binsphi,nBinsPt,binspt);
372 fHistEtaPhiPtRecK=
new TH3F(
"hEtaPhiPtRecK",
"",10,binseta,nBinsPhi,binsphi,nBinsPt,binspt);
374 fHistEtaPhiPtRecPro=
new TH3F(
"hEtaPhiPtRecPro",
"",10,binseta,nBinsPhi,binsphi,nBinsPt,binspt);
387 AliESDEvent *esd = (AliESDEvent*) (InputEvent());
391 printf(
"AliAnalysisTaskSDDRP::Exec(): bad ESD\n");
399 if(esd->GetRunNumber()<=139517) year=2010;
400 if(year==2010)
fESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kFALSE);
401 else fESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(kFALSE);
405 fESDtrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
406 AliESDtrackCuts::kAny);
409 Int_t nTracks=esd->GetNumberOfTracks();
412 for(Int_t it=0; it<nTracks; it++){
413 AliESDtrack* tr=esd->GetTrack(it);
414 UInt_t status=tr->GetStatus();
415 if(!(status&AliESDtrack::kITSrefit))
continue;
416 if(!(status&AliESDtrack::kTPCin))
continue;
421 const AliMultiplicity* mult=esd->GetMultiplicity();
422 Int_t nTracklets=mult->GetNumberOfTracklets();
424 for(Int_t it=0; it<nTracklets; it++){
425 Double_t eta=TMath::Abs(mult->GetEta(it));
426 if(eta<1) nTracklets1++;
431 const AliESDVertex *spdv=esd->GetVertex();
432 if(spdv && spdv->IsFromVertexer3D()){
437 if(spdv && spdv->IsFromVertexerZ()){
442 const AliESDVertex *trkv=esd->GetPrimaryVertex();
443 if(trkv && trkv->GetNContributors()>1){
451 AliMCEventHandler* eventHandler =
dynamic_cast<AliMCEventHandler*
> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
453 Printf(
"ERROR: Could not retrieve MC event handler");
456 AliMCEvent* mcEvent = eventHandler->MCEvent();
458 Printf(
"ERROR: Could not retrieve MC event");
461 stack = mcEvent->Stack();
463 Printf(
"ERROR: stack not available");
466 const AliVVertex* mcVert=mcEvent->GetPrimaryVertex();
468 Printf(
"ERROR: generated vertex not available");
471 if(TMath::Abs(mcVert->GetZ())>10)
return;
475 TString genname=mcEvent->GenEventHeader()->ClassName();
479 if(genname.Contains(
"CocktailEventHeader")){
480 AliGenCocktailEventHeader *cockhead=(AliGenCocktailEventHeader*)mcEvent->GenEventHeader();
481 lgen=cockhead->GetHeaders();
482 for(Int_t ig=0; ig<lgen->GetEntries(); ig++){
483 AliGenerator* gen=(AliGenerator*)lgen->At(ig);
484 TString
title=gen->GetName();
485 if(title.Contains(
"bchadr")) typeHF=1;
486 else if(title.Contains(
"chadr")) typeHF=0;
487 else if(title.Contains(
"bele")) typeHF=3;
488 else if(title.Contains(
"cele")) typeHF=2;
490 nColl=lgen->GetEntries();
493 TString genTitle=mcEvent->GenEventHeader()->GetTitle();
494 if(genTitle.Contains(
"bchadr")) typeHF=1;
495 else if(genTitle.Contains(
"chadr")) typeHF=0;
496 else if(genTitle.Contains(
"bele")) typeHF=3;
497 else if(genTitle.Contains(
"cele")) typeHF=2;
500 Int_t nParticles=stack->GetNtrack();
501 Double_t dNchdy = 0.;
505 for (Int_t i=0;i<nParticles;i++){
506 TParticle* part = (TParticle*)stack->Particle(i);
507 Int_t absPdg=TMath::Abs(part->GetPdgCode());
508 Int_t
pdg=part->GetPdgCode();
511 if(stack->IsPhysicalPrimary(i)){
512 Double_t eta=part->Eta();
519 if(TMath::Abs(eta)<0.5){
523 if(TMath::Abs(eta)<0.9){
528 if (part->Energy() != TMath::Abs(part->Pz())){
529 rapid=0.5*TMath::Log((part->Energy()+part->Pz())/(part->Energy()-part->Pz()));
542 else if(absPdg==411){
547 if(iTypeKKpi>0) iType=3;
551 else if(absPdg==413){
556 else if(absPdg==431){
559 if(iType==1 || iType==2) iPart=3;
561 else if(absPdg==4122){
566 if(iTypeV0==1) iType=5;
567 if(iTypeV0==2) iType=6;
570 if(iType>=0) iPart=4;
592 if(iSpecies<0)
continue;
605 Double_t distx=part->Vx()-mcVert->GetX();
606 Double_t disty=part->Vy()-mcVert->GetY();
607 Double_t distz=part->Vz()-mcVert->GetZ();
608 Double_t distToVert=TMath::Sqrt(distx*distx+disty*disty+distz*distz);
620 if(iPart<0)
continue;
621 if(iType<0)
continue;
623 if(iPart==0 && iType>0 && iType<=2){
625 }
else if(iPart==1 && iType>0 && iType<=3){
627 }
else if(iPart==3 && iType>0 && iType<=2){
631 if(iFromB==4 && iPart>=0 && iPart<5)
fHistYPtPrompt[iPart]->Fill(part->Pt(),rapid);
632 else if(iFromB==5 && iPart>=0 && iPart<5)
fHistYPtFeeddown[iPart]->Fill(part->Pt(),rapid);
635 for(Int_t i=0; i<nTracks; i++){
636 AliESDtrack* track=esd->GetTrack(i);
638 Int_t label=TMath::Abs(track->GetLabel());
639 if(stack->IsPhysicalPrimary(label)){
640 TParticle* part = (TParticle*)stack->Particle(label);
641 Int_t absPdg=TMath::Abs(part->GetPdgCode());
642 Double_t eta=part->Eta();
662 fOutput =
dynamic_cast<TList*
> (GetOutputData(1));
664 printf(
"ERROR: fOutput not available\n");
676 if(labLc<0)
return -1;
677 TParticle* dp = (TParticle*)stack->Particle(labLc);
678 Int_t pdgdp=dp->GetPdgCode();
679 Int_t nDau=dp->GetNDaughters();
685 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
686 if(iDau<0)
return -1;
687 TParticle* dau=(TParticle*)stack->Particle(iDau);
688 Int_t pdgdau=dau->GetPdgCode();
689 if(TMath::Abs(pdgdau)==321){
690 if(pdgdp>0 && pdgdau>0)
return -1;
691 if(pdgdp<0 && pdgdau<0)
return -1;
693 }
else if(TMath::Abs(pdgdau)==211){
694 if(pdgdp<0 && pdgdau>0)
return -1;
695 if(pdgdp>0 && pdgdau<0)
return -1;
697 }
else if(TMath::Abs(pdgdau)==2212){
698 if(pdgdp<0 && pdgdau>0)
return -1;
699 if(pdgdp>0 && pdgdau<0)
return -1;
705 if(nPions!=1)
return -1;
706 if(nKaons!=1)
return -1;
707 if(nProtons!=1)
return -1;
708 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
709 if(iDau<0)
return -1;
718 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
719 if(iDau<0)
return -1;
720 TParticle* dau=(TParticle*)stack->Particle(iDau);
721 Int_t pdgdau=dau->GetPdgCode();
722 if(TMath::Abs(pdgdau)==313 || TMath::Abs(pdgdau)==3124 || TMath::Abs(pdgdau)==2224
723 || TMath::Abs(pdgdau)==3122 || TMath::Abs(pdgdau)==311){
724 Int_t nDauRes=dau->GetNDaughters();
725 if(nDauRes!=2)
return -1;
726 for(Int_t resDau=dau->GetFirstDaughter(); resDau<=dau->GetLastDaughter(); resDau++){
727 if(resDau<0)
return -1;
728 TParticle* resdau=(TParticle*)stack->Particle(resDau);
729 Int_t pdgresdau=resdau->GetPdgCode();
730 if(TMath::Abs(pdgresdau)==321){
731 if(pdgdp>0 && pdgresdau>0)
return -1;
732 if(pdgdp<0 && pdgresdau<0)
return -1;
735 if(TMath::Abs(pdgresdau)==211){
736 if(pdgdp<0 && pdgresdau>0)
return -1;
737 if(pdgdp>0 && pdgresdau<0)
return -1;
740 if(TMath::Abs(pdgresdau)==2212){
741 if(pdgdp<0 && pdgresdau>0)
return -1;
742 if(pdgdp>0 && pdgresdau<0)
return -1;
746 }
else if(TMath::Abs(pdgdau)==321){
747 if(pdgdp>0 && pdgdau>0)
return -1;
748 if(pdgdp<0 && pdgdau<0)
return -1;
750 }
else if(TMath::Abs(pdgdau)==211){
751 if(pdgdp<0 && pdgdau>0)
return -1;
752 if(pdgdp>0 && pdgdau<0)
return -1;
754 }
else if(TMath::Abs(pdgdau)==2212){
755 if(pdgdp<0 && pdgdau>0)
return -1;
756 if(pdgdp>0 && pdgdau<0)
return -1;
762 if(nPions!=1)
return -1;
763 if(nKaons!=1)
return -1;
764 if(nProtons!=1)
return -1;
765 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
766 if(iDau<0)
return -1;
TH3F * fHistEtaPhiPtRecK
! histo of generated kaons
AliESDtrackCuts * fESDtrackCuts
0=pp, 1=PbPb, 2=pPb
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
static Int_t CheckDplusDecay(AliStack *stack, Int_t label, Int_t *arrayDauLab)
TH1F * fHistBSpecies
! histo of B hadron species
TH1F * fHistoSPDZVtxY
! histo with distr. of y coord. of SPDZ vertex
TH2F * fHistYPtPrompt[5]
! histo of y vs. pt from prompt D0, D+, D*, Ds, Lc
Int_t fSystem
c/b separation using quarks
TH3F * fHistEtaPhiPtGenEle
! histo of generated electrons
TH1F * fHistoSPDZVtxX
! histo with distr. of x coord. of SPDZ vertex
TH2F * fHistYPtDplusbyDecChannel[3]
! histo of y vs. pt for D+->Kpipi and D+->K0*pi
TH1F * fHistMotherID
! histo of mother ID
TH3F * fHistEtaPhiPtGenK
! histo of generated kaons
TH2F * fHistBYPtAllDecay[5]
! histo of y vs. pt from prompt B0, B+, B*, Bs, Lb
static Int_t CheckDstarDecay(AliStack *stack, Int_t label, Int_t *arrayDauLab)
TH1F * fHistOriginFeeddown
! histo of D production point (feeddown)
virtual ~AliAnalysisTaskCheckHFMCProd()
static Int_t CheckDsDecay(AliStack *stack, Int_t label, Int_t *arrayDauLab)
static Int_t CheckLcpKpiDecay(AliStack *stack, Int_t label, Int_t *arrayDauLab)
TH1F * fHistoTracklets
! histo with number of SPD tracklets
virtual void UserCreateOutputObjects()
static Int_t CheckOrigin(TClonesArray *arrayMC, AliAODMCParticle *mcPart, Bool_t searchUpToQuark=kTRUE)
virtual void UserExec(Option_t *option)
TH2F * fHistYPtDsbyDecChannel[2]
! histo of y vs. pt for Ds->phipi and Ds->K0*K
TH1F * fHistoSPD3DVtxX
! histo with distr. of x coord. of SPD3D vertex
TH2F * fHistoNbVsNc
! histo of n. b quarks vs. n c. quarks
TH1F * fHistoTRKVtxX
! histo with distr. of x coord. of TRK vertex
TH2F * fHistNcollHFtype
! histo of B hadron species
TH2F * fHistYPtFeeddown[5]
! histo of y vs. pt from feeddown D0, D+, D*, Ds, Lc
TH3F * fHistEtaPhiPtRecPro
! histo of generated protons
TH1F * fHistoSelTracks
! histo with number of TPC+ITS refit ESD tracks
static Int_t CheckLcV0bachelorDecay(AliStack *stack, Int_t label, Int_t *arrayDauLab)
TH1F * fHistoPhysPrim
! histo of n. of physical primaries in |eta|<0.5
TH1F * fHistoNEvents
! histo with N of events
TH3F * fHistEtaPhiPtGenPi
! histo of generated pions
TH3F * fHistEtaPhiPtGenPro
! histo of generated protons
TH1F * fHistOriginPrompt
! histo of D production point (prompt)
static Int_t CheckD0Decay(AliStack *stack, Int_t label, Int_t *arrayDauLab)
TH2F * fHistYPtPromptAllDecay[5]
! histo of y vs. pt from prompt D0, D+, D*, Ds, Lc, no selection on decay channel ...
TH2F * fHistoNcharmed
! histo of D mesons vs. dN/dy
static Int_t CheckDplusKKpiDecay(AliStack *stack, Int_t label, Int_t *arrayDauLab)
TH1F * fHistoTRKVtxZ
! histo with distr. of z coord. of TRK vertex
TH1F * fHistDSpecies
! histo of D hadron species
AliAnalysisTaskCheckHFMCProd()
TH2F * fHistYPtAllDecay[5]
! histo of y vs. pt from prompt D0, D+, D*, Ds, Lc, no selection on decay channel ...
TH1F * fHistoEtaPhysPrim
! histo of eta distribution of physical primaries
TH1F * fHistoSPD3DVtxY
! histo with distr. of y coord. of SPD3D vertex
TH1F * fHistoTrackletsEta1
! histo with number of SPD tracklets in |eta|<1
TH2F * fHistYPtD0byDecChannel[2]
! histo of y vs. pt for D0->Kpi and D0->Kpipipi
TH1F * fHistoTracks
! histo with number of ESD tracks
TH1F * fHistoPtPhysPrim
! histo of pt distribution of phys primaries
TH1F * fHistoSPDZVtxZ
! histo with distr. of z coord. of SPDZ vertex
TH1F * fHistoTRKVtxY
! histo with distr. of y coord. of TRK vertex
TH3F * fHistEtaPhiPtRecEle
! histo of generated electrons
TH2F * fHistYPtFeeddownAllDecay[5]
! histo of y vs. pt from prompt D0, D+, D*, Ds, Lc, no selection on decay channel ...
TList * fOutput
! list of output histos
TH1F * fHistLcDecayChan
! histo of Lc decay modes
TH1F * fHistoSPD3DVtxZ
! histo with distr. of z coord. of SPD3D vertex
TH3F * fHistEtaPhiPtRecPi
! histo of generated pions
Int_t CheckLcDecay(Int_t labLc, AliStack *stack) const
Bool_t fReadMC
track selection
virtual void Terminate(Option_t *option)