From: arossi Date: Thu, 1 Aug 2013 09:26:03 +0000 (+0000) Subject: 2D D meson efficiency correction + updates for running on MC samples X-Git-Url: http://git.uio.no/git/?a=commitdiff_plain;h=241994abf14dcf6404ecd9df655929c4023925af;p=u%2Fmrichter%2FAliRoot.git 2D D meson efficiency correction + updates for running on MC samples --- diff --git a/PWGHF/correlationHF/AliAnalysisTaskSED0Correlations.cxx b/PWGHF/correlationHF/AliAnalysisTaskSED0Correlations.cxx index 54a1304972a..dfacaabeb1a 100644 --- a/PWGHF/correlationHF/AliAnalysisTaskSED0Correlations.cxx +++ b/PWGHF/correlationHF/AliAnalysisTaskSED0Correlations.cxx @@ -52,6 +52,7 @@ #include "AliAnalysisTaskSE.h" #include "AliAnalysisTaskSED0Correlations.h" #include "AliNormalizationCounter.h" +#include "AliVertexingHFUtils.h" using std::cout; using std::endl; @@ -66,7 +67,6 @@ AliAnalysisTaskSE(), fBinLimsCorr(), fPtThreshLow(), fPtThreshUp(), - fD0Eff(), fEvents(0), fAlreadyFilled(kFALSE), fOutputMass(0), @@ -90,7 +90,10 @@ AliAnalysisTaskSE(), fSys(0), fEtaForCorrel(0), fIsRejectSDDClusters(0), - fFillGlobal(kTRUE) + fFillGlobal(kTRUE), + fMultEv(0.), + fSoftPiCut(kTRUE), + fMEAxisThresh(kFALSE) { // Default constructor @@ -103,7 +106,6 @@ AliAnalysisTaskSED0Correlations::AliAnalysisTaskSED0Correlations(const char *nam fBinLimsCorr(), fPtThreshLow(), fPtThreshUp(), - fD0Eff(), fEvents(0), fAlreadyFilled(kFALSE), fOutputMass(0), @@ -127,7 +129,10 @@ AliAnalysisTaskSED0Correlations::AliAnalysisTaskSED0Correlations(const char *nam fSys(0), fEtaForCorrel(0), fIsRejectSDDClusters(0), - fFillGlobal(kTRUE) + fFillGlobal(kTRUE), + fMultEv(0.), + fSoftPiCut(kTRUE), + fMEAxisThresh(kFALSE) { // Default constructor @@ -158,7 +163,6 @@ AliAnalysisTaskSED0Correlations::AliAnalysisTaskSED0Correlations(const AliAnalys fBinLimsCorr(source.fBinLimsCorr), fPtThreshLow(source.fPtThreshLow), fPtThreshUp(source.fPtThreshUp), - fD0Eff(source.fD0Eff), fEvents(source.fEvents), fAlreadyFilled(source.fAlreadyFilled), fOutputMass(source.fOutputMass), @@ -182,7 +186,10 @@ AliAnalysisTaskSED0Correlations::AliAnalysisTaskSED0Correlations(const AliAnalys fSys(source.fSys), fEtaForCorrel(source.fEtaForCorrel), fIsRejectSDDClusters(source.fIsRejectSDDClusters), - fFillGlobal(source.fFillGlobal) + fFillGlobal(source.fFillGlobal), + fMultEv(source.fMultEv), + fSoftPiCut(source.fSoftPiCut), + fMEAxisThresh(source.fMEAxisThresh) { // Copy constructor } @@ -239,7 +246,6 @@ AliAnalysisTaskSED0Correlations& AliAnalysisTaskSED0Correlations::operator=(cons fBinLimsCorr = orig.fBinLimsCorr; fPtThreshLow = orig.fPtThreshLow; fPtThreshUp = orig.fPtThreshUp; - fD0Eff = orig.fD0Eff; fEvents = orig.fEvents; fAlreadyFilled = orig.fAlreadyFilled; fOutputMass = orig.fOutputMass; @@ -264,6 +270,9 @@ AliAnalysisTaskSED0Correlations& AliAnalysisTaskSED0Correlations::operator=(cons fEtaForCorrel = orig.fEtaForCorrel; fIsRejectSDDClusters = orig.fIsRejectSDDClusters; fFillGlobal = orig.fFillGlobal; + fMultEv = orig.fMultEv; + fSoftPiCut = orig.fSoftPiCut; + fMEAxisThresh = orig.fMEAxisThresh; return *this; //returns pointer of the class } @@ -300,9 +309,9 @@ void AliAnalysisTaskSED0Correlations::UserCreateOutputObjects() fCorrelatorTr = new AliHFCorrelator("CorrelatorTr",fCutsTracks,fSys); fCorrelatorKc = new AliHFCorrelator("CorrelatorKc",fCutsTracks,fSys); fCorrelatorK0 = new AliHFCorrelator("CorrelatorK0",fCutsTracks,fSys); - fCorrelatorTr->SetDeltaPhiInterval(-TMath::Pi()/2+TMath::Pi()/32.,3*TMath::Pi()/2+TMath::Pi()/32.);// set the Delta Phi Interval you want - fCorrelatorKc->SetDeltaPhiInterval(-TMath::Pi()/2+TMath::Pi()/32.,3*TMath::Pi()/2+TMath::Pi()/32.); - fCorrelatorK0->SetDeltaPhiInterval(-TMath::Pi()/2+TMath::Pi()/32.,3*TMath::Pi()/2+TMath::Pi()/32.); + fCorrelatorTr->SetDeltaPhiInterval(-TMath::Pi()/2,3*TMath::Pi()/2);// set the Delta Phi Interval you want (in this case -0.5Pi to 1.5 Pi) + fCorrelatorKc->SetDeltaPhiInterval(-TMath::Pi()/2,3*TMath::Pi()/2); + fCorrelatorK0->SetDeltaPhiInterval(-TMath::Pi()/2,3*TMath::Pi()/2); fCorrelatorTr->SetEventMixing(fMixing);// sets the analysis on a single event (kFALSE) or mixed events (kTRUE) fCorrelatorKc->SetEventMixing(fMixing); fCorrelatorK0->SetEventMixing(fMixing); @@ -341,39 +350,40 @@ void AliAnalysisTaskSED0Correlations::UserCreateOutputObjects() TString nameMass=" ",nameSgn=" ", nameBkg=" ", nameRfl=" ",nameMassWg=" ",nameSgnWg=" ", nameBkgWg=" ", nameRflWg=" "; +//for origin c case (or for data) for(Int_t i=0;iGetNPtBins();i++){ - nameMass="histMass_"; + nameMass="histMass_"; if(fReadMC) nameMass+="c_"; nameMass+=i; - nameMassWg="histMass_WeigD0Eff_"; + nameMassWg="histMass_WeigD0Eff_"; if(fReadMC) nameMassWg+="c_"; nameMassWg+=i; - nameSgn="histSgn_"; + nameSgn="histSgn_"; if(fReadMC) nameSgn+="c_"; nameSgn+=i; - nameSgnWg="histSgn_WeigD0Eff_"; + nameSgnWg="histSgn_WeigD0Eff_"; if(fReadMC) nameSgnWg+="c_"; nameSgnWg+=i; - nameBkg="histBkg_"; + nameBkg="histBkg_"; if(fReadMC) nameBkg+="c_"; nameBkg+=i; - nameBkgWg="histBkg_WeigD0Eff_"; + nameBkgWg="histBkg_WeigD0Eff_"; if(fReadMC) nameBkgWg+="c_"; nameBkgWg+=i; - nameRfl="histRfl_"; + nameRfl="histRfl_"; if(fReadMC) nameRfl+="c_"; nameRfl+=i; - nameRflWg="histRfl_WeigD0Eff_"; + nameRflWg="histRfl_WeigD0Eff_"; if(fReadMC) nameRflWg+="c_"; nameRflWg+=i; //histograms of invariant mass distributions //MC signal if(fReadMC){ - TH1F* tmpSt = new TH1F(nameSgn.Data(), "D^{0} invariant mass - MC; M [GeV]; Entries",120,1.5648,2.1648); - TH1F* tmpStWg = new TH1F(nameSgnWg.Data(), "D^{0} invariant mass - MC; M [GeV] - weight 1/D0eff; Entries",120,1.5648,2.1648); + TH1F* tmpSt = new TH1F(nameSgn.Data(), "D^{0} invariant mass c - MC; M [GeV]; Entries",120,1.5648,2.1648); + TH1F* tmpStWg = new TH1F(nameSgnWg.Data(), "D^{0} invariant mass c - MC; M [GeV] - weight 1/D0eff; Entries",120,1.5648,2.1648); tmpSt->Sumw2(); tmpStWg->Sumw2(); //Reflection: histo filled with D0Mass which pass the cut (also) as D0bar and with D0bar which pass (also) the cut as D0 - TH1F* tmpRt = new TH1F(nameRfl.Data(), "Reflected signal invariant mass - MC; M [GeV]; Entries",120,1.5648,2.1648); - TH1F* tmpRtWg = new TH1F(nameRflWg.Data(), "Reflected signal invariant mass - MC - weight 1/D0eff; M [GeV]; Entries",120,1.5648,2.1648); - TH1F* tmpBt = new TH1F(nameBkg.Data(), "Background invariant mass - MC; M [GeV]; Entries",120,1.5648,2.1648); - TH1F* tmpBtWg = new TH1F(nameBkgWg.Data(), "Background invariant mass - MC - weight 1/D0eff; M [GeV]; Entries",120,1.5648,2.1648); + TH1F* tmpRt = new TH1F(nameRfl.Data(), "Reflected signal invariant mass c - MC; M [GeV]; Entries",120,1.5648,2.1648); + TH1F* tmpRtWg = new TH1F(nameRflWg.Data(), "Reflected signal invariant mass c - MC - weight 1/D0eff; M [GeV]; Entries",120,1.5648,2.1648); + TH1F* tmpBt = new TH1F(nameBkg.Data(), "Background invariant mass c - MC; M [GeV]; Entries",120,1.5648,2.1648); + TH1F* tmpBtWg = new TH1F(nameBkgWg.Data(), "Background invariant mass c - MC - weight 1/D0eff; M [GeV]; Entries",120,1.5648,2.1648); tmpBt->Sumw2(); tmpBtWg->Sumw2(); tmpRt->Sumw2(); @@ -387,15 +397,48 @@ void AliAnalysisTaskSED0Correlations::UserCreateOutputObjects() } //mass - TH1F* tmpMt = new TH1F(nameMass.Data(),"D^{0} invariant mass; M [GeV]; Entries",120,1.5648,2.1648); + TH1F* tmpMt = new TH1F(nameMass.Data(),"D^{0} invariant mass c; M [GeV]; Entries",120,1.5648,2.1648); tmpMt->Sumw2(); fOutputMass->Add(tmpMt); //mass weighted by 1/D0eff - TH1F* tmpMtwg = new TH1F(nameMassWg.Data(),"D^{0} invariant mass - weight 1/D0eff; M [GeV]; Entries",120,1.5648,2.1648); + TH1F* tmpMtwg = new TH1F(nameMassWg.Data(),"D^{0} invariant mass c - weight 1/D0eff; M [GeV]; Entries",120,1.5648,2.1648); tmpMtwg->Sumw2(); fOutputMass->Add(tmpMtwg); } +//for origin b case (no Bkg and Mass histos, here for weights you should use c+b efficiencies, while on data (on MC they're useless)) + for(Int_t i=0;iGetNPtBins();i++){ + + nameSgn="histSgn_b_"; + nameSgn+=i; + nameSgnWg="histSgn_WeigD0Eff_b_"; + nameSgnWg+=i; + nameRfl="histRfl_b_"; + nameRfl+=i; + nameRflWg="histRfl_WeigD0Eff_b_"; + nameRflWg+=i; + + //histograms of invariant mass distributions + + //MC signal + if(fReadMC){ + TH1F* tmpSt = new TH1F(nameSgn.Data(), "D^{0} invariant mass b - MC; M [GeV]; Entries",120,1.5648,2.1648); + TH1F* tmpStWg = new TH1F(nameSgnWg.Data(), "D^{0} invariant mass b - MC; M [GeV] - weight 1/D0eff; Entries",120,1.5648,2.1648); + tmpSt->Sumw2(); + tmpStWg->Sumw2(); + + //Reflection: histo filled with D0Mass which pass the cut (also) as D0bar and with D0bar which pass (also) the cut as D0 + TH1F* tmpRt = new TH1F(nameRfl.Data(), "Reflected signal invariant mass b - MC; M [GeV]; Entries",120,1.5648,2.1648); + TH1F* tmpRtWg = new TH1F(nameRflWg.Data(), "Reflected signal invariant mass b - MC - weight 1/D0eff; M [GeV]; Entries",120,1.5648,2.1648); + tmpRt->Sumw2(); + tmpRtWg->Sumw2(); + fOutputMass->Add(tmpSt); + fOutputMass->Add(tmpStWg); + fOutputMass->Add(tmpRt); + fOutputMass->Add(tmpRtWg); + } + } + const char* nameoutput=GetOutputSlot(2)->GetContainer()->GetName(); fNentries=new TH1F(nameoutput, "Integral(1,2) = number of AODs *** Integral(2,3) = number of candidates selected with cuts *** Integral(3,4) = number of D0 selected with cuts *** Integral(4,5) = events with good vertex *** Integral(5,6) = pt out of bounds", 20,-0.5,19.5); @@ -461,6 +504,8 @@ void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/) TClonesArray *inputArray=0; + fMultEv = 0.; //reset event multiplicity + if(!aod && AODEvent() && IsStandardAOD()) { // In case there is an AOD handler writing a standard AOD, use the AOD // event in memory rather than the input (ESD) event. @@ -546,7 +591,7 @@ void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/) } if(NMCevents && !isMCeventgood){ - std::cout << "The MC event " << eventType << " not interesting for this analysis: skipping" << std::endl; + if(fDebug>2)std::cout << "The MC event " << eventType << " not interesting for this analysis: skipping" << std::endl; return; } fNentries->Fill(19); //event with particular production type @@ -634,12 +679,15 @@ void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/) Int_t nSelectedloose=0,nSelectedtight=0; + //Fill Event Multiplicity (needed only in Reco) + fMultEv = (Double_t)(AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aod,-1.,1.)); + //RecoD0 case ************************************************ if(fRecoD0) { for (Int_t iD0toKpi = 0; iD0toKpi < nInD0toKpi; iD0toKpi++) { AliAODRecoDecayHF2Prong *d = (AliAODRecoDecayHF2Prong*)inputArray->UncheckedAt(iD0toKpi); - + if(d->Pt()<2.) continue; //to save time and merging memory... if(d->GetSelectionMap()) if(!d->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts)){ @@ -671,12 +719,18 @@ void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/) fCorrelatorK0->SetD0Properties(d,fIsSelectedCandidate); if(!fReadMC) { - if (TMath::Abs(d->Eta())Eta())FindObject(Form("hMultiplEvt_Bin%d",ptbin)))->Fill(fMultEv); + } } else { //correlations on MC -> association of selected D0 to MCinfo with MCtruth if (TMath::Abs(d->Eta())MatchToMC(421,mcArray,2,pdgDgD0toKpi); //return MC particle label if the array corresponds to a D0, -1 if not - if (labD0>-1) CalculateCorrelations(d,labD0,mcArray); + if (labD0>-1) { + CalculateCorrelations(d,labD0,mcArray); + if(!fMixing) ((TH1F*)fOutputStudy->FindObject(Form("hMultiplEvt_Bin%d",ptbin)))->Fill(fMultEv); //Fill multiplicity histo + } } } @@ -721,15 +775,20 @@ void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/) /* Int_t mother = mcPart->GetMother(); AliAODMCParticle* mcMoth = dynamic_cast(mcArray->At(mother)); if(!mcMoth) continue; - if(TMath::Abs(mcMoth->GetPdgCode())==413) continue; - */ + if(TMath::Abs(mcMoth->GetPdgCode())==413) continue;*/ + if (mcPart->GetPdgCode()==421) fIsSelectedCandidate = 1; else fIsSelectedCandidate = 2; - TString fillthis="histSgn_"; fillthis+=ptbin; + TString fillthis="histSgn_"; + if(CheckD0Origin(mcArray,mcPart)==4) fillthis+="c_"; + else if(CheckD0Origin(mcArray,mcPart)==5) fillthis+="b_"; + else continue; + fillthis+=ptbin; ((TH1F*)(fOutputMass->FindObject(fillthis)))->Fill(1.864); CalculateCorrelationsMCKine(mcPart,mcArray); + if(!fMixing) ((TH1F*)fOutputStudy->FindObject(Form("hMultiplEvt_Bin%d",ptbin)))->Fill(fMultEv); //Fill multiplicity histo } } } @@ -782,33 +841,53 @@ void AliAnalysisTaskSED0Correlations::FillMassHists(AliAODRecoDecayHF2Prong *par if ((fIsSelectedCandidate==1 || fIsSelectedCandidate==3) && fFillOnlyD0D0bar<2) { //D0 if(fReadMC){ - if(labD0>=0) { - AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0); + if(labD0>=0 && CheckD0Origin(arrMC,(AliAODMCParticle*)arrMC->At(labD0))==4) { + AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0); Int_t pdgD0 = partD0->GetPdgCode(); if (pdgD0==421){ //D0 - fillthis="histSgn_"; + fillthis="histSgn_c_"; fillthis+=ptbin; ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0); - fillthis="histSgn_WeigD0Eff_"; + fillthis="histSgn_WeigD0Eff_c_"; fillthis+=ptbin; - Double_t effD0 = fD0Eff.at(ptbin); + Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv); ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0); - } else{ //it was a D0bar - fillthis="histRfl_"; + } else{ //it was a D0bar + fillthis="histRfl_c_"; fillthis+=ptbin; - ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0); - fillthis="histRfl_WeigD0Eff_"; + ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0); + fillthis="histRfl_WeigD0Eff_c_"; fillthis+=ptbin; - Double_t effD0 = fD0Eff.at(ptbin); + Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv); ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0); - } + } + } else if(labD0>=0 && CheckD0Origin(arrMC,(AliAODMCParticle*)arrMC->At(labD0))==5) { + AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0); + Int_t pdgD0 = partD0->GetPdgCode(); + if (pdgD0==421){ //D0 + fillthis="histSgn_b_"; + fillthis+=ptbin; + ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0); + fillthis="histSgn_WeigD0Eff_b_"; + fillthis+=ptbin; + Double_t effD0 = fCutsTracks->GetTrigWeightB(part->Pt(),fMultEv); + ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0); + } else{ //it was a D0bar + fillthis="histRfl_b_"; + fillthis+=ptbin; + ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0); + fillthis="histRfl_WeigD0Eff_b_"; + fillthis+=ptbin; + Double_t effD0 = fCutsTracks->GetTrigWeightB(part->Pt(),fMultEv); + ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0); + } } else {//background - fillthis="histBkg_"; + fillthis="histBkg_c_"; fillthis+=ptbin; ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0); - fillthis="histBkg_WeigD0Eff_"; + fillthis="histBkg_WeigD0Eff_c_"; fillthis+=ptbin; - Double_t effD0 = fD0Eff.at(ptbin); + Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv); ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0); } }else{ @@ -817,7 +896,7 @@ void AliAnalysisTaskSED0Correlations::FillMassHists(AliAODRecoDecayHF2Prong *par ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0); fillthis="histMass_WeigD0Eff_"; fillthis+=ptbin; - Double_t effD0 = fD0Eff.at(ptbin); + Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv); ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0,1./effD0); } @@ -825,45 +904,65 @@ void AliAnalysisTaskSED0Correlations::FillMassHists(AliAODRecoDecayHF2Prong *par if (fIsSelectedCandidate>1 && (fFillOnlyD0D0bar==0 || fFillOnlyD0D0bar==2)) { //D0bar if(fReadMC){ - if(labD0>=0) { - AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0); + if(labD0>=0 && CheckD0Origin(arrMC,(AliAODMCParticle*)arrMC->At(labD0))==4) { + AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0); Int_t pdgD0 = partD0->GetPdgCode(); - - if (pdgD0==-421){ //D0bar - fillthis="histSgn_"; + if (pdgD0==-421){ //D0 + fillthis="histSgn_c_"; fillthis+=ptbin; ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar); - fillthis="histSgn_WeigD0Eff_"; + fillthis="histSgn_WeigD0Eff_c_"; fillthis+=ptbin; - Double_t effD0 = fD0Eff.at(ptbin); + Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv); ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0); - } else{ - fillthis="histRfl_"; + } else{ //it was a D0bar + fillthis="histRfl_c_"; fillthis+=ptbin; ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar); - fillthis="histRfl_WeigD0Eff_"; + fillthis="histRfl_WeigD0Eff_c_"; fillthis+=ptbin; - Double_t effD0 = fD0Eff.at(ptbin); + Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv); ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0); - } - } else {//background or LS - fillthis="histBkg_"; + } + } else if(labD0>=0 && CheckD0Origin(arrMC,(AliAODMCParticle*)arrMC->At(labD0))==5) { + AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0); + Int_t pdgD0 = partD0->GetPdgCode(); + if (pdgD0==-421){ //D0 + fillthis="histSgn_b_"; + fillthis+=ptbin; + ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar); + fillthis="histSgn_WeigD0Eff_b_"; + fillthis+=ptbin; + Double_t effD0 = fCutsTracks->GetTrigWeightB(part->Pt(),fMultEv); + ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0); + } else{ //it was a D0bar + fillthis="histRfl_b_"; + fillthis+=ptbin; + ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar); + fillthis="histRfl_WeigD0Eff_b_"; + fillthis+=ptbin; + Double_t effD0 = fCutsTracks->GetTrigWeightB(part->Pt(),fMultEv); + ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0); + } + } else {//background + fillthis="histBkg_c_"; fillthis+=ptbin; ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar); - fillthis="histBkg_WeigD0Eff_"; + fillthis="histBkg_WeigD0Eff_c_"; fillthis+=ptbin; - Double_t effD0 = fD0Eff.at(ptbin); + Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv); ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0); } }else{ fillthis="histMass_"; fillthis+=ptbin; - ((TH1F*)listout->FindObject(fillthis))->Fill(invmassD0bar); + ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar); fillthis="histMass_WeigD0Eff_"; fillthis+=ptbin; - Double_t effD0 = fD0Eff.at(ptbin); + Double_t effD0 = fCutsTracks->GetTrigWeight(part->Pt(),fMultEv); ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar,1./effD0); } + } return; @@ -965,13 +1064,13 @@ void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() { //These for limits in THnSparse (one per bin, same limits). //Vars: DeltaPhi, InvMass, PtTrack, Displacement, DeltaEta Int_t nBinsPhi[5] = {32,150,6,3,16}; - Double_t binMinPhi[5] = {-TMath::Pi()/2.+TMath::Pi()/32.,1.6,0.,0.,-1.6}; //is the minimum for all the bins - Double_t binMaxPhi[5] = {3*TMath::Pi()/2.+TMath::Pi()/32.,2.2,3.0,3.,1.6}; //is the maximum for all the bins + Double_t binMinPhi[5] = {-TMath::Pi()/2.,1.6,0.,0.,-1.6}; //is the minimum for all the bins + Double_t binMaxPhi[5] = {3.*TMath::Pi()/2.,2.2,3.0,3.,1.6}; //is the maximum for all the bins //Vars: DeltaPhi, InvMass, DeltaEta - Int_t nBinsMix[3] = {32,150,16}; - Double_t binMinMix[3] = {-TMath::Pi()/2.+TMath::Pi()/32.,1.6,-1.6}; //is the minimum for all the bins - Double_t binMaxMix[3] = {3*TMath::Pi()/2.+TMath::Pi()/32.,2.2,1.6}; //is the maximum for all the bins + Int_t nBinsMix[4] = {32,150,16,6}; + Double_t binMinMix[4] = {-TMath::Pi()/2.,1.6,-1.6,0.}; //is the minimum for all the bins + Double_t binMaxMix[4] = {3.*TMath::Pi()/2.,2.2,1.6,3.}; //is the maximum for all the bins for(Int_t i=0;iSumw2(); fOutputCorr->Add(hCorrLead); @@ -1121,35 +1220,35 @@ void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() { namePlot="hPhi_Lead_From_c_Bin"; namePlot+=i; - THnSparseF *hCorrLead_c = new THnSparseF(namePlot.Data(), "Leading particle correlations - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hCorrLead_c = new THnSparseF(namePlot.Data(), "Leading particle correlations - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); hCorrLead_c->Sumw2(); fOutputCorr->Add(hCorrLead_c); namePlot="hPhi_Lead_From_b_Bin"; namePlot+=i; - THnSparseF *hCorrLead_b = new THnSparseF(namePlot.Data(), "Leading particle correlations - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hCorrLead_b = new THnSparseF(namePlot.Data(), "Leading particle correlations - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); hCorrLead_b->Sumw2(); fOutputCorr->Add(hCorrLead_b); namePlot="hPhi_Lead_HF_From_c_Bin"; namePlot+=i; - THnSparseF *hCorrLead_HF_c = new THnSparseF(namePlot.Data(), "Leading particle correlations HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hCorrLead_HF_c = new THnSparseF(namePlot.Data(), "Leading particle correlations HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); hCorrLead_HF_c->Sumw2(); fOutputCorr->Add(hCorrLead_HF_c); namePlot="hPhi_Lead_HF_From_b_Bin"; namePlot+=i; - THnSparseF *hCorrLead_HF_b = new THnSparseF(namePlot.Data(), "Leading particle correlations HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hCorrLead_HF_b = new THnSparseF(namePlot.Data(), "Leading particle correlations HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); hCorrLead_HF_b->Sumw2(); fOutputCorr->Add(hCorrLead_HF_b); namePlot="hPhi_Lead_NonHF_Bin"; namePlot+=i; - THnSparseF *hCorrLead_Non = new THnSparseF(namePlot.Data(), "Leading particle correlations - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hCorrLead_Non = new THnSparseF(namePlot.Data(), "Leading particle correlations - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); hCorrLead_Non->Sumw2(); fOutputCorr->Add(hCorrLead_Non); } @@ -1158,38 +1257,38 @@ void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() { namePlot="hPhi_Weig_Bin"; namePlot+=i; - THnSparseF *hCorrWeig = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted); #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hCorrWeig = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted); #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); fOutputCorr->Add(hCorrWeig); if (fReadMC) { namePlot="hPhi_Weig_From_c_Bin"; namePlot+=i; - THnSparseF *hCorrWeig_c = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hCorrWeig_c = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); fOutputCorr->Add(hCorrWeig_c); namePlot="hPhi_Weig_From_b_Bin"; namePlot+=i; - THnSparseF *hCorrWeig_b = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hCorrWeig_b = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); fOutputCorr->Add(hCorrWeig_b); namePlot="hPhi_Weig_HF_From_c_Bin"; namePlot+=i; - THnSparseF *hCorrWeig_HF_c = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hCorrWeig_HF_c = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); fOutputCorr->Add(hCorrWeig_HF_c); namePlot="hPhi_Weig_HF_From_b_Bin"; namePlot+=i; - THnSparseF *hCorrWeig_HF_b = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hCorrWeig_HF_b = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); fOutputCorr->Add(hCorrWeig_HF_b); namePlot="hPhi_Weig_NonHF_Bin"; namePlot+=i; - THnSparseF *hCorrWeig_Non = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hCorrWeig_Non = new THnSparseF(namePlot.Data(), "Charged particle correlations (pT weighted) - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); fOutputCorr->Add(hCorrWeig_Non); } @@ -1216,6 +1315,13 @@ void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() { hDstarPions->GetXaxis()->SetBinLabel(2,"Rejected"); hDstarPions->SetMinimum(0); fOutputStudy->Add(hDstarPions); + + //Events multiplicity + Double_t yAxisMult[13] = {0, 4, 8, 12, 16, 20, 28, 36, 44, 100}; + namePlot = "hMultiplEvt_Bin"; namePlot+=i; + TH1F *hMultEv = new TH1F(namePlot.Data(), "Event multiplicity",9,yAxisMult); + hMultEv->SetMinimum(0); + fOutputStudy->Add(hMultEv); } @@ -1224,21 +1330,21 @@ void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() { namePlot="hPhi_K0_Bin"; namePlot+=i;namePlot+="_EvMix"; - THnSparseF *hPhiK_EvMix = new THnSparseF(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hPhiK_EvMix = new THnSparseF(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); hPhiK_EvMix->Sumw2(); fOutputCorr->Add(hPhiK_EvMix); namePlot="hPhi_Kcharg_Bin"; namePlot+=i;namePlot+="_EvMix"; - THnSparseF *hPhiH_EvMix = new THnSparseF(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hPhiH_EvMix = new THnSparseF(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); hPhiH_EvMix->Sumw2(); fOutputCorr->Add(hPhiH_EvMix); namePlot="hPhi_Charg_Bin"; namePlot+=i;namePlot+="_EvMix"; - THnSparseF *hPhiC_EvMix = new THnSparseF(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + THnSparseF *hPhiC_EvMix = new THnSparseF(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsMix,binMinMix,binMaxMix); hPhiC_EvMix->Sumw2(); fOutputCorr->Add(hPhiC_EvMix); } @@ -1452,6 +1558,14 @@ void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() { hOrigin_D0->GetXaxis()->SetBinLabel(1,"From c"); hOrigin_D0->GetXaxis()->SetBinLabel(2,"From b"); fOutputStudy->Add(hOrigin_D0); + + //primary tracks (Kine & Reco) + namePlot="hPhysPrim_Bin"; namePlot+=i; + TH1F *hPhysPrim = new TH1F(namePlot.Data(), "Origin of hadrons",2,0.,2.); + hPhysPrim->SetMinimum(0); + hPhysPrim->GetXaxis()->SetBinLabel(1,"OK"); + hPhysPrim->GetXaxis()->SetBinLabel(2,"NO"); + fOutputStudy->Add(hPhysPrim); } } } @@ -1524,7 +1638,7 @@ void AliAnalysisTaskSED0Correlations::CalculateCorrelations(AliAODRecoDecayHF2Pr AliReducedParticle* track = fCorrelatorTr->GetAssociatedParticle(); if(!fMixing) { - if(!track->CheckSoftPi()) { //removal of soft pions + if(fSoftPiCut && !track->CheckSoftPi()) { //removal of soft pions if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0); if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0bar); continue; @@ -1537,11 +1651,24 @@ void AliAnalysisTaskSED0Correlations::CalculateCorrelations(AliAODRecoDecayHF2Pr } if(track->Pt() < fPtThreshLow.at(ptbin) || track->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K + if(fReadMC) { + AliAODMCParticle* trkKine = (AliAODMCParticle*)mcArray->At(track->GetLabel()); + if (!trkKine) continue; + if (!trkKine->IsPhysicalPrimary()) { + ((TH1F*)fOutputStudy->FindObject(Form("hPhysPrim_Bin%d",ptbin)))->Fill(1.); + continue; //reject the Reco track if correspondent Kine track is not primary + } else ((TH1F*)fOutputStudy->FindObject(Form("hPhysPrim_Bin%d",ptbin)))->Fill(0.); + } + Double_t effTr = track->GetWeight(); //extract track efficiency - Double_t effD0 = fD0Eff.at(ptbin); + Double_t effD0 = 1.; + if(fReadMC) { + if(origD0==4) effD0 = fCutsTracks->GetTrigWeight(d->Pt(),fMultEv); + if(origD0==5) effD0 = fCutsTracks->GetTrigWeightB(d->Pt(),fMultEv); + } else effD0 = fCutsTracks->GetTrigWeight(d->Pt(),fMultEv); Double_t eff = effTr*effD0; - - FillSparsePlots(mcArray,mInv,origD0,PDGD0,track,ptbin,kTrack,1./eff); //fills for charged tracks, weight = 1./eff + + FillSparsePlots(mcArray,mInv,origD0,PDGD0,track,ptbin,kTrack,1./eff); //fills for charged tracks if(!fMixing) N_Charg++; @@ -1552,7 +1679,10 @@ void AliAnalysisTaskSED0Correlations::CalculateCorrelations(AliAODRecoDecayHF2Pr lead[0] = fCorrelatorTr->GetDeltaPhi(); lead[1] = fCorrelatorTr->GetDeltaEta(); lead[2] = fReadMC ? CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(track->GetLabel())) : 0; - lead[3] = 1./track->GetWeight(); //weight is 1./efficiency + if(fReadMC) { + if(origD0==4) lead[3] = 1./(track->GetWeight()*fCutsTracks->GetTrigWeight(d->Pt(),fMultEv)); //weight is 1./efficiency + if(origD0==5) lead[3] = 1./(track->GetWeight()*fCutsTracks->GetTrigWeightB(d->Pt(),fMultEv)); //weight is 1./efficiency + } else lead[3] = 1./(track->GetWeight()*fCutsTracks->GetTrigWeight(d->Pt(),fMultEv)); highPt = track->Pt(); } @@ -1583,7 +1713,7 @@ void AliAnalysisTaskSED0Correlations::CalculateCorrelations(AliAODRecoDecayHF2Pr AliReducedParticle* kCharg = fCorrelatorKc->GetAssociatedParticle(); if(!fMixing) { - if(!kCharg->CheckSoftPi()) { //removal of soft pions + if(fSoftPiCut && !kCharg->CheckSoftPi()) { //removal of soft pions if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0); if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0bar); continue; @@ -1635,20 +1765,20 @@ void AliAnalysisTaskSED0Correlations::CalculateCorrelations(AliAODRecoDecayHF2Pr } // end of charged kaons loop } //end of event loop for fCorrelatorK0 - Double_t fillSpLeadD0[3] = {lead[0],mD0,lead[1]}; - Double_t fillSpLeadD0bar[3] = {lead[0],mD0bar,lead[1]}; + Double_t fillSpLeadD0[4] = {lead[0],mD0,lead[1],0.4}; //dummy value for threshold of leading! + Double_t fillSpLeadD0bar[4] = {lead[0],mD0bar,lead[1],0.4}; //leading track correlations fill if(!fMixing) { if(fReadMC) { - if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==421) { //D0 + if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==421 && (fIsSelectedCandidate==1||fIsSelectedCandidate==3)) { //D0 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadD0,lead[3]); //c and b D0 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0,lead[3]); //c or b D0 if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0,lead[3]); if(origD0==5&&(int)lead[2]>=4&&(int)lead[2]<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0,lead[3]); if((int)lead[2]==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_NonHF_Bin%d",ptbin)))->Fill(fillSpLeadD0,lead[3]); //non HF } - if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==-421) { //D0bar + if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==-421 && fIsSelectedCandidate>1) { //D0bar ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadD0bar,lead[3]); ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0bar,lead[3]); //c or b D0 if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadD0bar,lead[3]); @@ -1744,8 +1874,14 @@ void AliAnalysisTaskSED0Correlations::CalculateCorrelationsMCKine(AliAODMCPartic AliAODMCParticle *trkMC = (AliAODMCParticle*)mcArray->At(track->GetLabel()); if(!trkMC) continue; - if (!trkMC->IsPhysicalPrimary()) continue; //reject material budget, or other fake tracks + if (!trkMC->IsPhysicalPrimary()) { //reject material budget, or other fake tracks + ((TH1F*)fOutputStudy->FindObject(Form("hPhysPrim_Bin%d",ptbin)))->Fill(1.); + continue; + } else ((TH1F*)fOutputStudy->FindObject(Form("hPhysPrim_Bin%d",ptbin)))->Fill(0.); + if (IsDDaughter(d,trkMC,mcArray)) continue; + if (fSoftPiCut && IsSoftPion_MCKine(d,trkMC,mcArray)) continue; //remove soft pions (if requestes, e.g. for templates) + FillSparsePlots(mcArray,mInv,origD0,PDGD0,track,ptbin,kTrack); //fills for charged tracks //retrieving leading info... @@ -1832,18 +1968,18 @@ void AliAnalysisTaskSED0Correlations::CalculateCorrelationsMCKine(AliAODMCPartic } // end of charged kaons loop } //end of event loop for fCorrelatorK0 - Double_t fillSpLeadMC[3] = {lead[0],mD0,lead[1]}; //mD0 = mD0bar = 1.864 + Double_t fillSpLeadMC[4] = {lead[0],mD0,lead[1],0.4}; //mD0 = mD0bar = 1.864 //leading track correlations fill if(!fMixing) { - if(d->GetPdgCode()==421) { //D0 + if(d->GetPdgCode()==421 && (fIsSelectedCandidate==1||fIsSelectedCandidate==3)) { //D0 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadMC); //c and b D0 ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC); //c or b D0 if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC); if(origD0==5&&(int)lead[2]>=4&&(int)lead[2]<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC); if((int)lead[2]==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_NonHF_Bin%d",ptbin)))->Fill(fillSpLeadMC); //non HF } - if(d->GetPdgCode()==-421) { //D0bar + if(d->GetPdgCode()==-421 && fIsSelectedCandidate>1) { //D0bar ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(fillSpLeadMC); ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC); //c or b D0 if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpLeadMC); @@ -1920,8 +2056,8 @@ void AliAnalysisTaskSED0Correlations::FillSparsePlots(TClonesArray* mcArray, Dou //variables for filling histos Double_t fillSpPhiD0[5] = {deltaphi,mD0,ptTrack,d0Track,deltaeta}; Double_t fillSpPhiD0bar[5] = {deltaphi,mD0bar,ptTrack,d0Track,deltaeta}; - Double_t fillSpWeigD0[3] = {deltaphi,mD0,deltaeta}; - Double_t fillSpWeigD0bar[3] = {deltaphi,mD0bar,deltaeta}; + Double_t fillSpWeigD0[4] = {deltaphi,mD0,deltaeta,ptTrack}; + Double_t fillSpWeigD0bar[4] = {deltaphi,mD0bar,deltaeta,ptTrack}; if(fReadMC == 0) { //sparse fill for data (tracks, K+-, K0) + weighted @@ -1944,7 +2080,7 @@ void AliAnalysisTaskSED0Correlations::FillSparsePlots(TClonesArray* mcArray, Dou if(origD0==4) {orig = "_From_c";} else {orig = "_From_b";} //sparse fill for data (tracks, K+-, K0) + weighted - if(PdgD0==421) { //D0 (from MCTruth) + if(PdgD0==421 && (fIsSelectedCandidate==1||fIsSelectedCandidate==3)) { //D0 (from MCTruth) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg); ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0,wg); if(origD0==4&&origTr>=1&&origTr<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0,wg); @@ -1956,7 +2092,7 @@ void AliAnalysisTaskSED0Correlations::FillSparsePlots(TClonesArray* mcArray, Dou if(origD0==5&&origTr>=4&&origTr<=8) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpWeigD0,pTorig*wg); if(origTr==0) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Weig_NonHF_Bin%d",ptbin)))->Fill(fillSpWeigD0,pTorig*wg); } - if(PdgD0==-421) { //D0bar + if(PdgD0==-421 && fIsSelectedCandidate>1) { //D0bar ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg); ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg); if(origD0==4&&origTr>=1&&origTr<=3) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg); @@ -1987,10 +2123,16 @@ void AliAnalysisTaskSED0Correlations::FillSparsePlots(TClonesArray* mcArray, Dou Double_t EtaLim_Sparse = ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d_EvMix",ptbin)))->GetAxis(2)->GetXmax(); if(deltaeta > EtaLim_Sparse) deltaeta = EtaLim_Sparse-0.01; if(deltaeta < -EtaLim_Sparse) deltaeta = -EtaLim_Sparse+0.01; + Double_t ptLim_Sparse = ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d_EvMix",ptbin)))->GetAxis(3)->GetXmax(); //all plots have same axes... + if(ptTrack > ptLim_Sparse) ptTrack = ptLim_Sparse-0.01; //variables for filling histos - Double_t fillSpPhiD0[3] = {deltaphi,mD0,deltaeta}; - Double_t fillSpPhiD0bar[3] = {deltaphi,mD0bar,deltaeta}; + Double_t fillSpPhiD0[4] = {deltaphi,mD0,deltaeta,0.4}; //dummy for ME threshold! unless explicitly set by flag... + Double_t fillSpPhiD0bar[4] = {deltaphi,mD0bar,deltaeta,0.4}; + if(fMEAxisThresh) { + fillSpPhiD0[3] = ptTrack; + fillSpPhiD0bar[3] = ptTrack; + } if(fReadMC == 0) { //sparse fill for data (tracks, K+-, K0) @@ -1999,8 +2141,8 @@ void AliAnalysisTaskSED0Correlations::FillSparsePlots(TClonesArray* mcArray, Dou } if(fReadMC == 1) { //sparse fill for data (tracks, K+-, K0) - if(PdgD0==421) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg); - if(PdgD0==-421) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg); + if(PdgD0==421 && (fIsSelectedCandidate==1||fIsSelectedCandidate==3)) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0,wg); + if(PdgD0==-421 && fIsSelectedCandidate>1) ((THnSparseF*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0bar,wg); }//end MC case } //end of ME fill @@ -2152,6 +2294,37 @@ Bool_t AliAnalysisTaskSED0Correlations::SelectV0(AliAODv0* v0, AliAODVertex *vtx return kTRUE; } +//--------------------------------------------------------------------------- +Bool_t AliAnalysisTaskSED0Correlations::IsSoftPion_MCKine(AliAODMCParticle* d, AliAODMCParticle* track, TClonesArray* arrayMC) const +{ + // + // Removes soft pions in Kine + + //Daughter removal in MCKine case + Bool_t isSoftPi = kFALSE; + Int_t labelD0 = d->GetLabel(); + + Int_t mother = track->GetMother(); + if(!mother) return isSoftPi; //safety check + + AliAODMCParticle* mcMoth = dynamic_cast(arrayMC->At(mother)); //it's the mother of the track! + if(TMath::Abs(mcMoth->GetPdgCode())==413 && mcMoth->GetNDaughters()==2) { //mother is D* with 2 daughs + Int_t labdau1 = mcMoth->GetDaughter(0); + Int_t labdau2 = mcMoth->GetDaughter(1); + AliAODMCParticle* dau1 = dynamic_cast(arrayMC->At(labdau1)); + AliAODMCParticle* dau2 = dynamic_cast(arrayMC->At(labdau2)); + if(!dau1 || !dau2) return isSoftPi; //safety check + if(dau1->GetLabel()==labelD0 || dau2->GetLabel()==labelD0) { //one of the daughs is the D0 trigger + if((TMath::Abs(dau1->GetPdgCode())==421 && TMath::Abs(dau2->GetPdgCode())==211)||(TMath::Abs(dau1->GetPdgCode())==211 && TMath::Abs(dau2->GetPdgCode())==421)) { + isSoftPi = kTRUE; //ok, soft pion was found + return isSoftPi; + } + } + } + + return isSoftPi; +} + //________________________________________________________________________ void AliAnalysisTaskSED0Correlations::PrintBinsAndLimits() { @@ -2170,10 +2343,6 @@ void AliAnalysisTaskSED0Correlations::PrintBinsAndLimits() { for (int i=0; i fBinLimsCorr; // limits of pt bins per correlations std::vector fPtThreshLow; // pT treshold of hadrons - low std::vector fPtThreshUp; // pT treshold of hadrons - up - std::vector fD0Eff; // D0 efficiencies Int_t fEvents; // EventCounter Bool_t fAlreadyFilled; // D0 in an event already analyzed (for track distribution plots) @@ -126,8 +130,11 @@ class AliAnalysisTaskSED0Correlations : public AliAnalysisTaskSE Double_t fEtaForCorrel; // cut for D0 eta to enable correlation with associated particles Bool_t fIsRejectSDDClusters; // flag to reject events with SDD clusters Bool_t fFillGlobal; // flag to fill global plots (in loops on tracks and V0 for each event) + Double_t fMultEv; // event multiplicity (for trigger eff) + Bool_t fSoftPiCut; // flag to activate soft pion cut on Data + Bool_t fMEAxisThresh; // flag to fill threshold axis in ME plots - ClassDef(AliAnalysisTaskSED0Correlations,3); // AliAnalysisTaskSE for D0->Kpi + ClassDef(AliAnalysisTaskSED0Correlations,4); // AliAnalysisTaskSE for D0->Kpi }; #endif diff --git a/PWGHF/correlationHF/macros/AddTaskD0Correlations.C b/PWGHF/correlationHF/macros/AddTaskD0Correlations.C index fd225111791..f739ce306a3 100644 --- a/PWGHF/correlationHF/macros/AddTaskD0Correlations.C +++ b/PWGHF/correlationHF/macros/AddTaskD0Correlations.C @@ -1,6 +1,4 @@ -AliAnalysisTaskSED0Correlations *AddTaskD0Correlations(Bool_t readMC=kFALSE, Bool_t mixing=kFALSE, Bool_t recoTrMC=kFALSE, Bool_t recoD0MC = kFALSE, Double_t etacorr=1.5, Int_t -system=0/*0=pp,1=PbPb*/, Int_t flagD0D0bar=0, Float_t minC=0, Float_t maxC=0, TString finDirname="Output", TString cutsfilename="D0toKpiCuts.root", TString cutsfilename2="AssocPartCuts_fBit0_woITS.root", TString -cutsD0name="D0toKpiCuts", TString cutsTrkname="AssociatedTrkCuts", TString effName = "3D_eff_wo_ITScls2_f0_p8eta.root", Bool_t flagAOD049=kFALSE, Int_t standardbins=1, Bool_t stdcuts=kFALSE, TString effD0name="D0Eff_From_c_wLimAcc.root") +AliAnalysisTaskSED0Correlations *AddTaskD0Correlations(Bool_t readMC=kFALSE, Bool_t mixing=kFALSE, Bool_t recoTrMC=kFALSE, Bool_t recoD0MC = kFALSE, Bool_t flagsoftpicut = kTRUE, Bool_t MEthresh = kFALSE, TString cutsfilename="D0toKpiCuts.root", TString cutsfilename2="AssocPartCuts_fBit0_woITS_EffMap2D.root", TString effD0namec="D0Eff_From_c_2D.root", TString effD0nameb="D0Eff_From_b_2D.root", TString effName = "3D_eff_wo_ITScls2_f0_p8eta_FineBins.root", TString cutsD0name="D0toKpiCuts", TString cutsTrkname="AssociatedTrkCuts", Double_t etacorr=1.5, Int_t system=0/*0=pp,1=PbPb*/, Int_t flagD0D0bar=0, Float_t minC=0, Float_t maxC=0, TString finDirname="Output", Bool_t flagAOD049=kFALSE, Int_t standardbins=1, Bool_t stdcuts=kFALSE) { // // AddTask for the AliAnalysisTaskSE for D0 candidates @@ -79,14 +77,12 @@ cutsD0name="D0toKpiCuts", TString cutsTrkname="AssociatedTrkCuts", TString effNa // printf(" d0d0 [cm^2] < %f\n",fD0CorrCuts[7]); // printf(" cosThetaPoint > %f\n",fD0CorrCuts[8]); - TFile* filecuts; - filecuts=TFile::Open(cutsfilename.Data()); + TFile* filecuts=new TFile(cutsfilename.Data()); if(!filecuts->IsOpen()){ cout<<"Input file not found for D0 cuts: using std cut object"<IsOpen()){ cout<<"Input file not found for tracks cuts!"<IsOpen()){ cout<<"Input file not found for efficiency! Exiting..."<Get("c"); TH3D *h3D = (TH3D*)c->FindObject("heff_rebin"); + //Load D0 efficiency map + TFile* fileeffD0c=new TFile(effD0namec.Data()); + if(!fileeffD0c->IsOpen()){ + cout<<"Input file not found for efficiency! Exiting..."<Get("c1"); + TH2D *hEffD0c = (TH2D*)cc->FindObject("h_Eff"); + + //Load D0 efficiency map from b + if(readMC) { + TFile* fileeffD0b=new TFile(effD0nameb.Data()); + if(!fileeffD0b->IsOpen()){ + cout<<"Input file not found for efficiency! Exiting..."<Get("c1"); + TH2D *hEffD0b = (TH2D*)cb->FindObject("h_Eff"); + } + //Cuts for correlated tracks/K0 AliHFAssociatedTrackCuts* corrCuts=new AliHFAssociatedTrackCuts(); corrCuts = (AliHFAssociatedTrackCuts*)filecuts2->Get(cutsTrkname.Data()); @@ -140,7 +155,9 @@ cutsD0name="D0toKpiCuts", TString cutsTrkname="AssociatedTrkCuts", TString effNa } corrCuts->SetTrackCutsNames(); corrCuts->SetvZeroCutsNames(); - corrCuts->SetEfficiencyWeightMap(h3D); + if(recoD0MC) corrCuts->SetEfficiencyWeightMap(h3D); //data and MC Reco + if(recoD0MC) corrCuts->SetTriggerEffWeightMap(hEffD0c); //data and MC Reco + if(recoD0MC && readMC) corrCuts->SetTriggerEffWeightMapB(hEffD0b); //MC Reco corrCuts->PrintAll(); TString centr=""; @@ -166,6 +183,8 @@ cutsD0name="D0toKpiCuts", TString cutsTrkname="AssociatedTrkCuts", TString effNa massD0Task->SetFillOnlyD0D0bar(flagD0D0bar); massD0Task->SetSystem(system); //0=pp, 1=PbPb massD0Task->SetEtaForCorrel(etacorr); + massD0Task->SetSoftPiFlag(flagsoftpicut); + massD0Task->SetMEAxisThresh(MEthresh); //********************* //correlation settings @@ -185,22 +204,6 @@ cutsD0name="D0toKpiCuts", TString cutsTrkname="AssociatedTrkCuts", TString effNa massD0Task->SetPtTreshLow(pttreshlow); massD0Task->SetPtTreshUp(pttreshup); -//******************** -//D0 efficiency -//******************** - - TFile* fileeffD0; - fileeffD0=TFile::Open(effD0name.Data()); - if(!fileeffD0->IsOpen()){ - cout<<"Input file not found for efficiency! Exiting..."<Get("c1"); - TH1D *hEff = (TH1D*)c1->FindObject("h_Eff"); - Double_t effD0values[15] = {1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.}; //only 4 to 11 bins are filled with efficiencies for now! - for(int i=3; i<11; i++) {effD0values[i]=hEff->GetBinContent(i+1);} - massD0Task->SetD0Efficiency(effD0values); - // massD0Task->SetRejectSDDClusters(kTRUE); // massD0Task->SetWriteVariableTree(kTRUE);