From bce70c9644fcea320d91195a3ed8a71ff09ab118 Mon Sep 17 00:00:00 2001 From: prino Date: Tue, 10 Jul 2012 08:49:05 +0000 Subject: [PATCH] New class for steering D-hadron correlation analysis (Sandro) --- PWGHF/CMakelibPWGHFcorrelationHF.pkg | 1 + .../AliAnalysisTaskDStarCorrelations.cxx | 546 +++----- .../AliAnalysisTaskDStarCorrelations.h | 22 +- .../AliAnalysisTaskSED0Correlations.cxx | 1181 ++++++++--------- .../AliAnalysisTaskSED0Correlations.h | 21 +- .../AliHFAssociatedTrackCuts.cxx | 168 ++- .../correlationHF/AliHFAssociatedTrackCuts.h | 48 +- PWGHF/correlationHF/AliHFCorrelator.cxx | 425 ++++++ PWGHF/correlationHF/AliHFCorrelator.h | 160 +++ PWGHF/correlationHF/AliReducedParticle.cxx | 24 +- PWGHF/correlationHF/AliReducedParticle.h | 23 +- .../macros/AddTaskD0Correlations.C | 5 +- .../macros/AddTaskDStarCorrelations.C | 240 ++-- .../macros/makeTFileAssociatedTrackCuts.C | 197 +-- 14 files changed, 1815 insertions(+), 1246 deletions(-) create mode 100644 PWGHF/correlationHF/AliHFCorrelator.cxx create mode 100644 PWGHF/correlationHF/AliHFCorrelator.h diff --git a/PWGHF/CMakelibPWGHFcorrelationHF.pkg b/PWGHF/CMakelibPWGHFcorrelationHF.pkg index 765c1605382..47d4c1c1d49 100644 --- a/PWGHF/CMakelibPWGHFcorrelationHF.pkg +++ b/PWGHF/CMakelibPWGHFcorrelationHF.pkg @@ -34,6 +34,7 @@ set ( CLASS_HDRS AliAnalysisTaskDxHFEParticleSelection.h AliAnalysisTaskDxHFECorrelation.h AliHFAssociatedTrackCuts.h + AliHFCorrelator.h AliReducedParticle.h AliAnalysisTaskDStarCorrelations.h AliAnalysisTaskSED0Correlations.h diff --git a/PWGHF/correlationHF/AliAnalysisTaskDStarCorrelations.cxx b/PWGHF/correlationHF/AliAnalysisTaskDStarCorrelations.cxx index e9bc77ec2f6..a9886a251ec 100644 --- a/PWGHF/correlationHF/AliAnalysisTaskDStarCorrelations.cxx +++ b/PWGHF/correlationHF/AliAnalysisTaskDStarCorrelations.cxx @@ -27,7 +27,7 @@ /* $Id$ */ -#include +//#include #include #include #include @@ -40,7 +40,7 @@ #include "AliAODRecoCascadeHF.h" #include "AliAODRecoDecayHF2Prong.h" #include "AliAODPidHF.h" -#include "AliEventPoolManager.h" +//#include "AliEventPoolManager.h" #include "AliVParticle.h" #include "AliAnalysisManager.h" #include "AliAODInputHandler.h" @@ -49,11 +49,10 @@ #include "AliAODMCParticle.h" #include "AliNormalizationCounter.h" #include "AliReducedParticle.h" +#include "AliHFCorrelator.h" -using std::cout; -using std::endl; ClassImp(AliAnalysisTaskDStarCorrelations) @@ -62,14 +61,17 @@ ClassImp(AliAnalysisTaskDStarCorrelations) AliAnalysisTaskDStarCorrelations::AliAnalysisTaskDStarCorrelations() : AliAnalysisTaskSE(), fhandler(0x0), -fPoolMgr(0x0), +//fPoolMgr(0x0), fmcArray(0x0), fCounter(0x0), +fCorrelator(0x0), fselect(0), fmontecarlo(kFALSE), fmixing(kFALSE), +fSystem(kFALSE), fEvents(0), fDebug(0), +fDisplacement(0), fOutput(0x0), fCuts(0), @@ -83,14 +85,17 @@ AliAnalysisTaskDStarCorrelations::AliAnalysisTaskDStarCorrelations(const Char_t* AliAnalysisTaskSE(name), fhandler(0x0), -fPoolMgr(0x0), +//fPoolMgr(0x0), fmcArray(0x0), fCounter(0x0), +fCorrelator(0x0), fselect(0), fmontecarlo(kFALSE), fmixing(kFALSE), +fSystem(kFALSE), fEvents(0), fDebug(0), +fDisplacement(0), fOutput(0x0), fCuts(0), @@ -115,9 +120,10 @@ AliAnalysisTaskDStarCorrelations::~AliAnalysisTaskDStarCorrelations() { Info("AliAnalysisTaskDStarCorrelations","Calling Destructor"); if(fhandler) {delete fhandler; fhandler = 0;} - if(fPoolMgr) {delete fPoolMgr; fPoolMgr = 0;} + //if(fPoolMgr) {delete fPoolMgr; fPoolMgr = 0;} if(fmcArray) {delete fmcArray; fmcArray = 0;} if(fCounter) {delete fCounter; fCounter = 0;} + if(fCorrelator) {delete fCorrelator; fCorrelator = 0;} if(fOutput) {delete fOutput; fOutput = 0;} if(fCuts) {delete fCuts; fCuts = 0;} if(fCuts2) {delete fCuts2; fCuts2=0;} @@ -132,6 +138,10 @@ void AliAnalysisTaskDStarCorrelations::Init(){ if(fDebug > 1) printf("AliAnalysisTaskDStarCorrelations::Init() \n"); AliRDHFCutsDStartoKpipi* copyfCuts=new AliRDHFCutsDStartoKpipi(*fCuts); + + + + // Post the D* cuts PostData(2,copyfCuts); @@ -153,39 +163,43 @@ void AliAnalysisTaskDStarCorrelations::UserCreateOutputObjects(){ // define histograms DefineHistoForAnalysis(); - fCounter = new AliNormalizationCounter(Form("%s",GetOutputSlot(3)->GetContainer()->GetName())); fCounter->Init(); - - // definition of the Pool Manager for Event Mixing + Double_t Pi = TMath::Pi(); + fCorrelator = new AliHFCorrelator("Correlator",fCuts2,fSystem); // fCuts2 is the hadron cut object, fSystem to switch between pp or PbPb + fCorrelator->SetDeltaPhiInterval((-0.5-1./32)*Pi,(1.5-1./32)*Pi); // set correct phi interval + //fCorrelator->SetDeltaPhiInterval(-1.57,4.71); + fCorrelator->SetEventMixing(fmixing); //set kFALSE/kTRUE for mixing Off/On + fCorrelator->SetAssociatedParticleType(fselect); // set 1/2/3 for hadron/kaons/kzeros + fCorrelator->SetApplyDisplacementCut(fDisplacement); //set kFALSE/kTRUE for using the displacement cut + fCorrelator->SetUseMC(fmontecarlo); + Bool_t pooldef = fCorrelator->DefineEventPool(); - Int_t MaxNofEvents = 200; - Int_t MinNofTracks = 1000; - - Int_t NofMultiplicityBins = 5; - Double_t MBins[]={0,20,40,60,80,500}; - Double_t * MultiplicityBins = MBins; - - Int_t NofZVrtxBins = 5; - Double_t ZBins[]={-10,-5,-2.5,2.5,5,10}; - Double_t *ZVrtxBins = ZBins; - - - fPoolMgr = new AliEventPoolManager(MaxNofEvents, MinNofTracks, NofMultiplicityBins, MultiplicityBins, NofZVrtxBins, ZVrtxBins); + if(!pooldef) AliInfo("Warning:: Event pool not defined properly"); + } //_________________________________________________ void AliAnalysisTaskDStarCorrelations::UserExec(Option_t *){ - cout << " " << endl; - cout << "=================================================================================" << endl; + cout << "USER EXEC CHECKPOINT 1" << endl; - if(fselect==1) cout << "TASK::Correlation with hadrons "<< endl; - if(fselect==2) cout << "TASK::Correlation with kaons "<< endl; - if(fselect==3) cout << "TASK::Correlation with kzeros "<< endl; + //Double_t pi = TMath::Pi(); + cout << " " << endl; + cout << "=================================================================================" << endl; + if(!fmixing){ + if(fselect==1) cout << "TASK::Correlation with hadrons on SE "<< endl; + if(fselect==2) cout << "TASK::Correlation with kaons on SE "<< endl; + if(fselect==3) cout << "TASK::Correlation with kzeros on SE "<< endl; + } + if(fmixing){ + if(fselect==1) cout << "TASK::Correlation with hadrons on ME "<< endl; + if(fselect==2) cout << "TASK::Correlation with kaons on ME "<< endl; + if(fselect==3) cout << "TASK::Correlation with kzeros on ME "<< endl; + } if (!fInputEvent) { Error("UserExec","NO EVENT FOUND!"); return; @@ -196,35 +210,30 @@ void AliAnalysisTaskDStarCorrelations::UserExec(Option_t *){ AliError("AOD event not found!"); return; } - Double_t pi = TMath::Pi(); + + fCorrelator->SetAODEvent(aodEvent); // set the event to be processed fEvents++; // event counter ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(0); fmcArray = dynamic_cast(aodEvent->FindListObject(AliAODMCParticle::StdBranchName())); + if(fmontecarlo && !fmcArray){ AliError("Array of MC particles not found"); return; } - - // initialize the pool for event mixing - Int_t multiplicity = aodEvent->GetNTracks(); - AliAODVertex *vtx = aodEvent->GetPrimaryVertex(); - Double_t zvertex = vtx->GetZ(); - Double_t multip = multiplicity; - - if(TMath::Abs(zvertex)>=10 || multip>500 || multip == 0) { - AliInfo(Form("Event with Zvertex = %.2f cm and multiplicity = %.0f out of pool bounds, SKIPPING",zvertex,multip)); + Bool_t isEvSel=fCuts->IsEventSelected(aodEvent); + if(!isEvSel) return; + ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(1); + // + Bool_t correlatorON = fCorrelator->Initialize(); //define the pool for mixing + if(!correlatorON) { + AliInfo("AliHFCorrelator didn't initialize the pool correctly or processed a bad event"); return; } - - AliEventPool* pool = fPoolMgr->GetEventPool(multip, zvertex); - if (!pool) AliFatal(Form("No pool found for multiplicity = %f, zVtx = %f cm", multip, zvertex)); - - - + if(fmontecarlo) fCorrelator->SetMCArray(fmcArray); // D* reconstruction - + cout << "USER EXEC CHECKPOINT 2" << endl; TClonesArray *arrayDStartoD0pi=0; @@ -257,20 +266,21 @@ void AliAnalysisTaskDStarCorrelations::UserExec(Option_t *){ Double_t deltainvMDStar; - Double_t mPDGD0=TDatabasePDG::Instance()->GetParticle(421)->Mass(); - Double_t mPDGDstar=TDatabasePDG::Instance()->GetParticle(413)->Mass(); + Double_t mPDGD0=1.8648;//TDatabasePDG::Instance()->GetParticle(421)->Mass(); + Double_t mPDGDstar=2.01022;//TDatabasePDG::Instance()->GetParticle(413)->Mass(); + cout << "USER EXEC CHECKPOINT 3" << endl; - - if(fselect ==3){// check the K0 invariant mass - TObjArray * fillkzerospectra = AcceptAndReduceKZero(aodEvent, 0,0); - delete fillkzerospectra; - } + //if(fselect ==3){// check the K0 invariant mass + //TObjArray * fillkzerospectra = AcceptAndReduceKZero(aodEvent, 0,0); + //delete fillkzerospectra; + //} //MC tagging for DStar //D* and D0 prongs needed to MatchToMC method Int_t pdgDgDStartoD0pi[2]={421,211}; Int_t pdgDgD0toKpi[2]={321,211}; + //loop on D* candidates for (Int_t iDStartoD0pi = 0; iDStartoD0piGetEntriesFast(); iDStartoD0pi++) { @@ -283,26 +293,24 @@ void AliAnalysisTaskDStarCorrelations::UserExec(Option_t *){ invMassDZero = - 999; deltainvMDStar = -998; - - AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)arrayDStartoD0pi->At(iDStartoD0pi); if(!dstarD0pi->GetSecondaryVtx()) continue; AliAODRecoDecayHF2Prong* theD0particle = (AliAODRecoDecayHF2Prong*)dstarD0pi->Get2Prong(); if (!theD0particle) continue; + // track quality cuts Int_t isTkSelected = fCuts->IsSelected(dstarD0pi,AliRDHFCuts::kTracks); // quality cuts on tracks // region of interest + topological cuts + PID Int_t isSelected=fCuts->IsSelected(dstarD0pi,AliRDHFCuts::kCandidate); //selected - //apply selections if(!isTkSelected) continue; if(!isSelected) continue; if(!fCuts->IsInFiducialAcceptance(dstarD0pi->Pt(),dstarD0pi->YDstar())) continue; - + Int_t mcLabelDStar = -999; if(fmontecarlo){ // find associated MC particle for D* ->D0toKpi - Int_t mcLabelDStar = dstarD0pi->MatchToMC(413,421,pdgDgDStartoD0pi,pdgDgD0toKpi,fmcArray); + mcLabelDStar = dstarD0pi->MatchToMC(413,421,pdgDgDStartoD0pi,pdgDgD0toKpi,fmcArray,kFALSE); if(mcLabelDStar>=0) isDStarMCtag = kTRUE; } @@ -310,31 +318,45 @@ void AliAnalysisTaskDStarCorrelations::UserExec(Option_t *){ phiDStar = dstarD0pi->Phi(); etaDStar = dstarD0pi->Eta(); + phiDStar = fCorrelator->SetCorrectPhiRange(phiDStar); // set the phi of the D meson in the correct range + Int_t ptbin=fCuts->PtBin(dstarD0pi->Pt()); Double_t dmDStarWindow =0.0019;// 0.0019 = 3 sigma Double_t mD0Window=0.074; - if (ptbin==1) mD0Window = 0.026; //0.5-1 - if (ptbin==2) mD0Window = 0.022; //1-2 - if (ptbin==3) mD0Window = 0.024; //2-3 - if (ptbin==4) mD0Window = 0.032; - if (ptbin==5) mD0Window = 0.032; - if (ptbin==6) mD0Window = 0.036; - if (ptbin==7) mD0Window = 0.036; - if (ptbin==8) mD0Window = 0.036; - if (ptbin==9) mD0Window = 0.058; - if (ptbin==10) mD0Window = 0.058; - if (ptbin>10) mD0Window = 0.074; - + if (!fSystem){ // pp + if (ptbin==1) mD0Window = 0.026; //0.5-1 + if (ptbin==2) mD0Window = 0.022; //1-2 + if (ptbin==3) mD0Window = 0.024; //2-3 + if (ptbin==4) mD0Window = 0.032; + if (ptbin==5) mD0Window = 0.032; + if (ptbin==6) mD0Window = 0.036; + if (ptbin==7) mD0Window = 0.036; + if (ptbin==8) mD0Window = 0.036; + if (ptbin==9) mD0Window = 0.058; + if (ptbin==10) mD0Window = 0.058; + if (ptbin>10) mD0Window = 0.074; + } + if(fSystem){// PbPb + if (ptbin==0) mD0Window = 0.032; //1-1 + if (ptbin==1) mD0Window = 0.032; //2-3 + if (ptbin==2) mD0Window = 0.032; //3-4 + if (ptbin==3) mD0Window = 0.032; //4-5 + if (ptbin==4) mD0Window = 0.036; //5-6 + if (ptbin==5) mD0Window = 0.036; //6-8 + if (ptbin==6) mD0Window = 0.055; //8-12 + if (ptbin==7) mD0Window = 0.074; //12-16 + if (ptbin==8) mD0Window = 0.074; //16-24 + if (ptbin==9) mD0Window = 0.074; //24-35 + } invMassDZero = dstarD0pi->InvMassD0(); ((TH2F*)fOutput->FindObject("D0InvMass"))->Fill(ptDStar,invMassDZero); deltainvMDStar = dstarD0pi->DeltaInvMass(); - - + //good candidates if (TMath::Abs(invMassDZero-mPDGD0)FindObject("RecoPtDStar"))->Fill(ptDStar); isInPeak = kTRUE; + ((TH2F*)fOutput->FindObject("PhiEtaTrigger"))->Fill(phiDStar,etaDStar); } }// end if good candidates @@ -354,155 +377,113 @@ void AliAnalysisTaskDStarCorrelations::UserExec(Option_t *){ if(TMath::Abs(deltainvMDStar-(mPDGDstar-mPDGD0))FindObject("RecoPtBkg"))->Fill(ptDStar); isInSideBand = kTRUE; + ((TH2F*)fOutput->FindObject("PhiEtaSideBand"))->Fill(phiDStar,etaDStar); } }//end if sidebands - + // getting the number of triggers in the MCtag D* case + + if(fmontecarlo && isDStarMCtag) ((TH1F*)fOutput->FindObject("MCtagPtDStar"))->Fill(ptDStar); + if(!isInPeak && !isInSideBand) continue; // skip if it is not side band or peak event - SAVE CPU TIME - + + fCorrelator->SetTriggerParticleProperties(ptDStar,phiDStar,etaDStar); // pass to the object the necessary trigger part parameters + + + Int_t trackiddaugh0 = ((AliAODTrack*)theD0particle->GetDaughter(0))->GetID(); Int_t trackiddaugh1 = ((AliAODTrack*)theD0particle->GetDaughter(1))->GetID(); Int_t trackidsoftPi = ((AliAODTrack*)dstarD0pi->GetBachelor())->GetID(); - ptDStar = dstarD0pi->Pt(); - phiDStar = dstarD0pi->Phi(); - etaDStar = dstarD0pi->Eta(); + Bool_t execPool = fCorrelator->ProcessEventPool(); + if(fmixing && !execPool) { + AliInfo("Mixed event analysis: pool is not ready"); + continue; + } + + Int_t NofEventsinPool = 1; + if(fmixing) NofEventsinPool = fCorrelator->GetNofEventsInPool(); - if(!fmixing){ // analysis on Single Event - - TObjArray* selectedtracks = new TObjArray(); - if(fselect==1 || fselect ==2) selectedtracks = AcceptAndReduceTracks(aodEvent); - if(fselect==3) {selectedtracks = AcceptAndReduceKZero(aodEvent,iDStartoD0pi,1);} - Int_t noftracks = selectedtracks->GetEntriesFast(); - Int_t counterPeak =0; - Int_t counterSB = 0; - - for(Int_t i =0; iAt(i); - Double_t phiHad=redpart->Phi(); - - Double_t ptHad=redpart->Pt(); - Double_t etaHad=redpart->Eta(); - - Int_t label = redpart->GetLabel(); - - Int_t trackid = redpart->GetID(); - - - // check that the track is not a D* daughter - if(trackid == trackiddaugh0) continue; - if(trackid == trackiddaugh1) continue; - if(trackid == trackidsoftPi) continue; + for (Int_t jMix =0; jMix < NofEventsinPool; jMix++){// loop on events in the pool; if it is SE analysis, stops at one + + Bool_t analyzetracks = fCorrelator->ProcessAssociatedTracks(jMix); + + if(!analyzetracks) { + AliInfo("AliHFCorrelator::Cannot process the track array"); + continue; + } + + + Int_t NofTracks = fCorrelator->GetNofTracks(); + + for(Int_t iTrack = 0; iTrackCorrelate(iTrack); + if(!runcorrelation) continue; + + Double_t DeltaPhi = fCorrelator->GetDeltaPhi(); + Double_t DeltaEta = fCorrelator->GetDeltaEta(); + + AliReducedParticle * hadron = fCorrelator->GetAssociatedParticle(); + + Double_t ptHad = hadron->Pt(); + Double_t phiHad = hadron->Phi(); + Double_t etaHad = hadron->Eta(); + Double_t label = hadron->GetLabel(); + Double_t trackid = hadron->GetID(); + + phiHad = fCorrelator->SetCorrectPhiRange(phiHad); + + if(!fmixing){ // skip D* Daughetrs + if(trackid == trackiddaugh0) continue; + if(trackid == trackiddaugh1) continue; + if(trackid == trackidsoftPi) continue; + } + + // from here on it is up to the user to decide what object to fill if(fmontecarlo && isDStarMCtag){ // check correlations of MC tagged DStars in MonteCarlo - + Int_t PartSource = fCuts2->IsMCpartFromHF(label,fmcArray); // check source of associated particle (hadron/kaon/K0) - - FillMCTagCorrelations(ptDStar,phiDStar,etaDStar,ptHad,phiHad,etaHad,PartSource); + + FillMCTagCorrelations(ptDStar,DeltaPhi,DeltaEta,ptHad,PartSource); + } + if(isInPeak) { - FillCorrelations(ptDStar,phiDStar,etaDStar,phiHad,etaHad);// function for correlations - counterPeak++; - if (phiDStar > 1.5*pi) phiDStar = phiDStar - 2*pi; - if (phiDStar < -0.5*pi) phiDStar = phiDStar + 2*pi; - - ((TH1F*)fOutput->FindObject("PhiTrigger"))->Fill(phiDStar); - - - - if (phiHad > 1.5*pi) phiHad = phiHad - 2*pi; - if (phiHad < -0.5*pi) phiHad = phiHad + 2*pi; - ((TH1F*)fOutput->FindObject("PhiPart"))->Fill(phiHad); - - } - if(isInSideBand) { - - FillSideBandCorrelations(ptDStar,phiDStar,etaDStar,phiHad,etaHad); // function for sidebands - if (phiDStar > 1.5*pi) phiDStar = phiDStar - 2*pi; - if (phiDStar < -0.5*pi) phiDStar = phiDStar + 2*pi; - ((TH1F*)fOutput->FindObject("PhiSideBand"))->Fill(phiDStar); - - counterSB++; - } + if(fselect==1) ((TH3D*)fOutput->FindObject("DPhiDStarHadron"))->Fill(DeltaPhi,ptDStar,DeltaEta); + if(fselect==2) ((TH3D*)fOutput->FindObject("DPhiDStarKaon"))->Fill(DeltaPhi,ptDStar,DeltaEta); + if(fselect==3) ((TH3D*)fOutput->FindObject("DPhiDStarKZero"))->Fill(DeltaPhi,ptDStar,DeltaEta); + ((TH2F*)fOutput->FindObject("PhiEtaPart"))->Fill(phiHad,etaHad); + //counterPeak++; // count tracks per peak per event - } // end loop on tracks - - if(counterPeak) ((TH1D*)fOutput->FindObject("NofTracksInPeak"))->Fill(counterPeak); - if(counterSB) ((TH1D*)fOutput->FindObject("NofTracksInSB"))->Fill(counterSB); - - - + } - } // end if SE Analysis - - if(fmixing) { // analysis on Mixed Events - if (pool->IsReady()|| pool->GetCurrentNEvents()>=5){ // check if the pool is ready - - pool->PrintInfo(); - - Int_t multbinflag = pool->MultBinIndex(); - Int_t zvtxflag = pool->ZvtxBinIndex(); - - if(isInPeak) ((TH2I*)fOutput->FindObject("EventMixingCheck"))->Fill(multbinflag,zvtxflag); + if(isInSideBand) { + + if(fselect==1) ((TH3D*)fOutput->FindObject("bkgDPhiDStarHadron"))->Fill(DeltaPhi,ptDStar,DeltaEta); + if(fselect==2) ((TH3D*)fOutput->FindObject("bkgDPhiDStarKaon"))->Fill(DeltaPhi,ptDStar,DeltaEta); + if(fselect==3) ((TH3D*)fOutput->FindObject("bkgDPhiDStarKZero"))->Fill(DeltaPhi,ptDStar,DeltaEta); + - - TObjArray* mixedtracks = 0x0; - for (Int_t jMix=0; jMixGetCurrentNEvents(); jMix++) {//loop over the events in the pool - mixedtracks = pool->GetEvent(jMix); - if(!mixedtracks) cout << "No Mixed tracks " << endl; - Int_t jMax = mixedtracks->GetEntriesFast(); - - for(Int_t iMix =0; iMixAt(iMix); - Double_t phiHad=redpart->Phi(); - Double_t etaHad=redpart->Eta(); - Double_t ptHad=redpart->Pt(); - Int_t label = redpart->GetLabel(); - - - if(fmontecarlo && isDStarMCtag){ // check correlations of MC tagged DStars in MonteCarlo - - Int_t PartSource = fCuts2->IsMCpartFromHF(label,fmcArray); // check source of associated particle (hadron/kaon/K0) - - FillMCTagCorrelations(ptDStar,phiDStar,etaDStar,ptHad,phiHad,etaHad,PartSource); - } - - if(isInPeak) { - FillCorrelations(ptDStar,phiDStar,etaDStar,phiHad,etaHad);// function for correlations - - if (phiDStar > 1.5*pi) phiDStar = phiDStar - 2*pi; - if (phiDStar < -0.5*pi) phiDStar = phiDStar + 2*pi; - - ((TH1F*)fOutput->FindObject("PhiTrigger"))->Fill(phiDStar); - if (phiHad > 1.5*pi) phiHad = phiHad - 2*pi; - if (phiHad < -0.5*pi) phiHad = phiHad + 2*pi; - ((TH1F*)fOutput->FindObject("PhiPart"))->Fill(phiHad); - - - - } - - if(isInSideBand) FillSideBandCorrelations(ptDStar,phiDStar,etaDStar,phiHad,etaHad); // function for sidebands - - } // end loop on tracks - }// end loop on events in pool - } // end if pool is ready + //counterSB++; + } - } // end ME analysis - + + } // end loop on track candidates + } // end loop on events in the pool + }// end loop on D* candidates + cout << "USER EXEC CHECKPOINT 4" << endl; + Bool_t updated = fCorrelator->PoolUpdate(); + if(!updated) AliInfo("Pool was not updated"); - if(fmixing) { // update the pool for Event Mixing - if(fselect==1 || fselect==2)pool->UpdatePool(AcceptAndReduceTracks(aodEvent)); // updating the pool for hadrons - if(fselect==3) pool->UpdatePool(AcceptAndReduceKZero(aodEvent,0,0)); // updating the pool for K0s - } - //cout << ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> END OF THE EVENT <<<<<<<<<<<<<<<<<<<<<<<<<<<<<" << endl; + //cout << ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> END OF THE EVENT <<<<<<<<<<<<<<<<<<<<<<<<<<<<<" << endl; PostData(1,fOutput); // set the outputs PostData(3,fCounter); // set the outputs @@ -528,93 +509,6 @@ void AliAnalysisTaskDStarCorrelations::Terminate(Option_t*) } -//_____________________________________________________ -TObjArray* AliAnalysisTaskDStarCorrelations::AcceptAndReduceTracks(AliAODEvent* inputEvent){ - - Int_t nTracks = inputEvent->GetNTracks(); - AliAODVertex * vtx = inputEvent->GetPrimaryVertex(); - Double_t Bz = inputEvent->GetMagneticField(); - - - TObjArray* tracksClone = new TObjArray; - tracksClone->SetOwner(kTRUE); - for (Int_t iTrack=0; iTrackGetTrack(iTrack); - if (! track) continue; - - if(!fCuts2->IsHadronSelected(track,vtx,Bz)) continue; // apply selection of hadrons - - - if(fselect ==2){ - if(!fCuts2->CheckKaonCompatibility(track,fmontecarlo,fmcArray)) continue; // check if it is a Kaon - data and MC - } - - AliVParticle * part = (AliVParticle*)track; - tracksClone->Add(new AliReducedParticle(part->Eta(), part->Phi(), part->Pt(),track->GetLabel(),track->GetID())); - - } - return tracksClone; -} - -//_____________________________________________________ -TObjArray* AliAnalysisTaskDStarCorrelations::AcceptAndReduceKZero(AliAODEvent* inputEvent, Int_t loopindex, Int_t plotassociation){ - - Int_t nOfVZeros = inputEvent->GetNumberOfV0s(); - TObjArray* KZeroClone = new TObjArray; - AliAODVertex *vertex1 = (AliAODVertex*)inputEvent->GetPrimaryVertex(); - Int_t v0label = -1; - Int_t pdgDgK0toPipi[2] = {211,211}; - Double_t mPDGK0=TDatabasePDG::Instance()->GetParticle(310)->Mass(); - const Int_t size = inputEvent->GetNumberOfV0s(); - Int_t prevnegID[size]; - Int_t prevposID[size]; - for(Int_t iv0 =0; iv0< nOfVZeros; iv0++){// loop on all v0 candidates - AliAODv0 *v0 = (static_cast (inputEvent))->GetV0(iv0); - if(!v0) {cout << "is not a v0 at step " << iv0 << endl; continue;} - - if(!fCuts2->IsKZeroSelected(v0,vertex1)) continue; // check if it is a k0 - - // checks to avoid double counting - Int_t negID = -999; - Int_t posID = -998; - //Int_t a = 0; - prevnegID[iv0] = -997; - prevposID[iv0]= -996; - negID = v0->GetNegID(); - posID = v0->GetPosID(); - Bool_t DoubleCounting = kFALSE; - - for(Int_t k=0; kMatchToMC(310,fmcArray, 0, pdgDgK0toPipi); //match a K0 short - Double_t k0pt = v0->Pt(); - Double_t k0eta = v0->Eta(); - Double_t k0Phi = v0->Phi(); - Double_t k0InvMass = v0->MassK0Short(); - - if (loopindex == 0) { - if(!plotassociation) ((TH2F*)fOutput->FindObject("KZeroSpectra"))->Fill(k0InvMass,k0pt); // spectra for all k0 - if(plotassociation) ((TH2F*)fOutput->FindObject("KZeroSpectraifHF"))->Fill(k0InvMass,k0pt); // spectra for k0 in association with a D* - } - // if there are more D* candidates per event, loopindex == 0 makes sure you fill the mass spectra only once! - - if(TMath::Abs(k0InvMass-mPDGK0)>3*0.004) continue; // select candidates within 3 sigma - KZeroClone->Add(new AliReducedParticle(k0eta , k0Phi, k0pt,v0label,0)); - - } - - return KZeroClone; -} - //_____________________________________________________ void AliAnalysisTaskDStarCorrelations::DefineHistoForAnalysis(){ @@ -647,6 +541,9 @@ void AliAnalysisTaskDStarCorrelations::DefineHistoForAnalysis(){ TH1F *RecoPtBkg = new TH1F("RecoPtBkg","RECO pt distribution side bands",30,0,30); fOutput->Add(RecoPtBkg); + TH1F *MCtagPtDStar = new TH1F("MCtagPtDStar","RECO pt of MCtagged DStars side bands",30,0,30); + fOutput->Add(MCtagPtDStar); + TH2F *KZeroSpectra = new TH2F("KZeroSpectra","Spectra of K0s",500,0.3,0.8,250,0,25); if(fselect==3) fOutput->Add(KZeroSpectra); @@ -672,14 +569,14 @@ void AliAnalysisTaskDStarCorrelations::DefineHistoForAnalysis(){ MCSources->GetXaxis()->SetBinLabel(5," from b->B"); if(fmontecarlo) fOutput->Add(MCSources); - TH1F * PhiTrigger = new TH1F("PhiTrigger","#phi distribution of the trigger particle",36,-0.5*Pi,1.5*Pi); - fOutput->Add(PhiTrigger); + TH2F * PhiEtaTrigger = new TH2F("PhiEtaTrigger","#phi distribution of the trigger particle",36,-0.5*Pi,1.5*Pi,18,-0.9,0.9); + fOutput->Add(PhiEtaTrigger); - TH1F * PhiSideBand = new TH1F("PhiSideBand","#phi distribution of the sideband particle",36,-0.5*Pi,1.5*Pi); - fOutput->Add(PhiSideBand); + TH2F * PhiEtaSideBand = new TH2F("PhiEtaSideBand","#phi distribution of the sideband particle",36,-0.5*Pi,1.5*Pi,18,-0.9,0.9); + fOutput->Add(PhiEtaSideBand); - TH1F * PhiPart = new TH1F("PhiPart","#phi distribution of the associated particle",36,-0.5*Pi,1.5*Pi); - fOutput->Add(PhiPart); + TH2F * PhiEtaPart = new TH2F("PhiEtaPart","#phi distribution of the associated particle",36,-0.5*Pi,1.5*Pi,18,-0.9,0.9); + fOutput->Add(PhiEtaPart); //correlations histograms @@ -745,68 +642,23 @@ void AliAnalysisTaskDStarCorrelations::DefineHistoForAnalysis(){ } -//____________________________ Function for correlations ___________________________________________________ -void AliAnalysisTaskDStarCorrelations::FillCorrelations(Double_t ptTrig, Double_t phiTrig, Double_t etaTrig, Double_t phiTrack, Double_t etaTrack){ - Double_t pi = TMath::Pi(); - Double_t deltaPhi, deltaEta; - deltaPhi = phiTrig - phiTrack; - deltaEta = etaTrig - etaTrack; - // set correct Delta Phi range - if (deltaPhi > 1.5*pi -pi/32) deltaPhi = deltaPhi - 2*pi; - if (deltaPhi < -0.5*pi -pi/32) deltaPhi = deltaPhi + 2*pi; - - if(fselect==1) ((TH3D*)fOutput->FindObject("DPhiDStarHadron"))->Fill(deltaPhi,ptTrig,deltaEta); - if(fselect==2) ((TH3D*)fOutput->FindObject("DPhiDStarKaon"))->Fill(deltaPhi,ptTrig,deltaEta); - if(fselect==3) ((TH3D*)fOutput->FindObject("DPhiDStarKZero"))->Fill(deltaPhi,ptTrig,deltaEta); - - - return; -} -//____________________________ Function for sidebands ___________________________________________________ -void AliAnalysisTaskDStarCorrelations::FillSideBandCorrelations(Double_t ptTrig, Double_t phiTrig, Double_t etaTrig, Double_t phiTrack, Double_t etaTrack){ - - Double_t pi = TMath::Pi(); - Double_t deltaPhi, deltaEta; - deltaPhi = phiTrig - phiTrack; - deltaEta = etaTrig - etaTrack; - // set correct Delta Phi range - if (deltaPhi > 1.5*pi -pi/32) deltaPhi = deltaPhi - 2*pi; - if (deltaPhi < -0.5*pi -pi/32) deltaPhi = deltaPhi + 2*pi; - - if(fselect==1) ((TH3D*)fOutput->FindObject("bkgDPhiDStarHadron"))->Fill(deltaPhi,ptTrig,deltaEta); - if(fselect==2) ((TH3D*)fOutput->FindObject("bkgDPhiDStarKaon"))->Fill(deltaPhi,ptTrig,deltaEta); - if(fselect==3) ((TH3D*)fOutput->FindObject("bkgDPhiDStarKZero"))->Fill(deltaPhi,ptTrig,deltaEta); - - - return; - - -} +//____________________________ Function for MC correlations ___________________________________________________ +void AliAnalysisTaskDStarCorrelations::FillMCTagCorrelations(Double_t ptTrig, Double_t DelPhi, Double_t DelEta, Double_t ptTrack, Int_t mcSource){ -//____________________________ Function for sidebands ___________________________________________________ -void AliAnalysisTaskDStarCorrelations::FillMCTagCorrelations(Double_t ptTrig, Double_t phiTrig, Double_t etaTrig, Double_t ptTrack, Double_t phiTrack, Double_t etaTrack, Int_t mcSource){ -Double_t deltaPhi = phiTrig - phiTrack; -Double_t deltaEta = etaTrig - etaTrack; -// set correct Delta Phi range - - Double_t pi = TMath::Pi(); - -if (deltaPhi > 1.5*pi -pi/32) deltaPhi = deltaPhi - 2*pi; -if (deltaPhi < -0.5*pi -pi/32) deltaPhi = deltaPhi + 2*pi; - if(fselect==1) ((TH3D*)fOutput->FindObject("MCTagDPhiDStarHadron"))->Fill(deltaPhi,ptTrig,deltaEta); - if(fselect==2 && ptTrack <1.5) ((TH3D*)fOutput->FindObject("MCTagDPhiDStarKaon"))->Fill(deltaPhi,ptTrig,deltaEta); - if(fselect==3) ((TH3D*)fOutput->FindObject("MCTagDPhiDStarKZero"))->Fill(deltaPhi,ptTrig,deltaEta); + if(fselect==1) ((TH3D*)fOutput->FindObject("MCTagDPhiDStarHadron"))->Fill(DelPhi,ptTrig,DelEta); + if(fselect==2 && ptTrack <1.5) ((TH3D*)fOutput->FindObject("MCTagDPhiDStarKaon"))->Fill(DelPhi,ptTrig,DelEta); + if(fselect==3) ((TH3D*)fOutput->FindObject("MCTagDPhiDStarKZero"))->Fill(DelPhi,ptTrig,DelEta); ((TH1F*)fOutput->FindObject("MCSources"))->Fill(0); if (mcSource==44){ // is from charm ->D - if(fselect==1) ((TH3D*)fOutput->FindObject("CharmDOriginDPhiDStarHadronMC"))->Fill(deltaPhi,ptTrig,deltaEta); - if(fselect==2 && ptTrack <1.5) ((TH3D*)fOutput->FindObject("CharmDOriginDPhiDStarKaonMC"))->Fill(deltaPhi,ptTrig,deltaEta); - if(fselect==3) ((TH3D*)fOutput->FindObject("CharmDOriginDPhiDStarKZeroMC"))->Fill(deltaPhi,ptTrig,deltaEta); + if(fselect==1) ((TH3D*)fOutput->FindObject("CharmDOriginDPhiDStarHadronMC"))->Fill(DelPhi,ptTrig,DelEta); + if(fselect==2 && ptTrack <1.5) ((TH3D*)fOutput->FindObject("CharmDOriginDPhiDStarKaonMC"))->Fill(DelPhi,ptTrig,DelEta); + if(fselect==3) ((TH3D*)fOutput->FindObject("CharmDOriginDPhiDStarKZeroMC"))->Fill(DelPhi,ptTrig,DelEta); ((TH1F*)fOutput->FindObject("MCSources"))->Fill(1); @@ -814,17 +666,17 @@ if (mcSource==44){ // is from charm ->D } if (mcSource==54){ // is from beauty -> D - if(fselect==1) ((TH3D*)fOutput->FindObject("BeautyDOriginDPhiDStarHadronMC"))->Fill(deltaPhi,ptTrig,deltaEta); - if(fselect==2 && ptTrack <1.5) ((TH3D*)fOutput->FindObject("BeautyDOriginDPhiDStarKaonMC"))->Fill(deltaPhi,ptTrig,deltaEta); - if(fselect==3) ((TH3D*)fOutput->FindObject("BeautyDOriginDPhiDStarKZeroMC"))->Fill(deltaPhi,ptTrig,deltaEta); + if(fselect==1) ((TH3D*)fOutput->FindObject("BeautyDOriginDPhiDStarHadronMC"))->Fill(DelPhi,ptTrig,DelEta); + if(fselect==2 && ptTrack <1.5) ((TH3D*)fOutput->FindObject("BeautyDOriginDPhiDStarKaonMC"))->Fill(DelPhi,ptTrig,DelEta); + if(fselect==3) ((TH3D*)fOutput->FindObject("BeautyDOriginDPhiDStarKZeroMC"))->Fill(DelPhi,ptTrig,DelEta); if(fselect==3) ((TH1F*)fOutput->FindObject("MCSources"))->Fill(1); if(fselect==3) ((TH1F*)fOutput->FindObject("MCSources"))->Fill(3); } if (mcSource==55){ // is from beauty ->B - if(fselect==1) ((TH3D*)fOutput->FindObject("BeautyBOriginDPhiDStarHadronMC"))->Fill(deltaPhi,ptTrig,deltaEta); - if(fselect==2 && ptTrack <1.5) ((TH3D*)fOutput->FindObject("BeautyBOriginDPhiDStarKaonMC"))->Fill(deltaPhi,ptTrig,deltaEta); - if(fselect==3) ((TH3D*)fOutput->FindObject("BeautyOriginBDPhiDStarKZeroMC"))->Fill(deltaPhi,ptTrig,deltaEta); + if(fselect==1) ((TH3D*)fOutput->FindObject("BeautyBOriginDPhiDStarHadronMC"))->Fill(DelPhi,ptTrig,DelEta); + if(fselect==2 && ptTrack <1.5) ((TH3D*)fOutput->FindObject("BeautyBOriginDPhiDStarKaonMC"))->Fill(DelPhi,ptTrig,DelEta); + if(fselect==3) ((TH3D*)fOutput->FindObject("BeautyOriginBDPhiDStarKZeroMC"))->Fill(DelPhi,ptTrig,DelEta); if(fselect==3) ((TH1F*)fOutput->FindObject("MCSources"))->Fill(1); if(fselect==3) ((TH1F*)fOutput->FindObject("MCSources"))->Fill(4); } diff --git a/PWGHF/correlationHF/AliAnalysisTaskDStarCorrelations.h b/PWGHF/correlationHF/AliAnalysisTaskDStarCorrelations.h index 3930e95a437..9a4ec13e401 100644 --- a/PWGHF/correlationHF/AliAnalysisTaskDStarCorrelations.h +++ b/PWGHF/correlationHF/AliAnalysisTaskDStarCorrelations.h @@ -37,6 +37,7 @@ #include "AliAODRecoCascadeHF.h" #include "AliHFAssociatedTrackCuts.h" #include "AliNormalizationCounter.h" +#include "AliHFCorrelator.h" class TParticle ; class TClonesArray ; @@ -66,22 +67,17 @@ class AliAnalysisTaskDStarCorrelations : public AliAnalysisTaskSE virtual void UserExec(Option_t *option); virtual void Terminate(Option_t *option); - - // methods to get the tracks to correlate - TObjArray* AcceptAndReduceTracks(AliAODEvent* inputEvent); - TObjArray* AcceptAndReduceKZero(AliAODEvent* inputEvent, Int_t loopindex, Int_t plotassociation); void DefineHistoForAnalysis(); // correlators - void FillCorrelations(Double_t ptTrig, Double_t phiTrig, Double_t etaTrig, Double_t phiTrack, Double_t etaTrack); - void FillSideBandCorrelations(Double_t ptTrig, Double_t phiTrig, Double_t etaTrig, Double_t phiTrack, Double_t etaTrack); - void FillMCTagCorrelations(Double_t ptTrig, Double_t phiTrig, Double_t etaTrig, Double_t ptTrack, Double_t phiTrack, Double_t etaTrack, Int_t mcSource); - + void FillMCTagCorrelations(Double_t ptTrig, Double_t DelPhi, Double_t DelEta, Double_t ptTrack, Int_t mcSource); // setters void SetMonteCarlo(Bool_t k) {fmontecarlo = k;} void SetUseMixing (Bool_t j) {fmixing = j;} - void SetCorrelator(Int_t l) {fselect = l;} + void SetCorrelator(Int_t l) {fselect = l;} // select 1 for hadrons, 2 for Kaons, 3 for Kzeros + void SetUseDisplacement(Int_t m) {fDisplacement=m;} // select 0 for no displ, 1 for abs displ, 2 for d0/sigma_d0 + void SetRunPbPb(Bool_t system){fSystem=system;} // select between pp (kFALSE) or PbPb (kTRUE) @@ -91,23 +87,27 @@ class AliAnalysisTaskDStarCorrelations : public AliAnalysisTaskSE AliAnalysisTaskDStarCorrelations& operator=(const AliAnalysisTaskDStarCorrelations& source); TObject* fhandler; //! Analysis Handler - AliEventPoolManager* fPoolMgr; //! event pool manager + // AliEventPoolManager* fPoolMgr; //! event pool manager TClonesArray* fmcArray; //mcarray AliNormalizationCounter *fCounter; // counter + AliHFCorrelator * fCorrelator; // object for correlations + Int_t fselect; // select what to correlate with a D* 1-chargedtracks,2-chargedkaons,3-k0s Bool_t fmontecarlo;//switch for MC Bool_t fmixing;// switch for event mixing + Bool_t fSystem; // pp or PbPb Int_t fEvents; //! number of event Int_t fDebug; //! debug level + Int_t fDisplacement; // set 0 for no displacement cut, 1 for absolute d0, 2 for d0/sigma_d0 TList *fOutput; //! user output data AliRDHFCutsDStartoKpipi *fCuts; // Cuts D* AliHFAssociatedTrackCuts *fCuts2; // cuts for associated - ClassDef(AliAnalysisTaskDStarCorrelations,1); // class for D meson correlations + ClassDef(AliAnalysisTaskDStarCorrelations,2); // class for D meson correlations }; #endif diff --git a/PWGHF/correlationHF/AliAnalysisTaskSED0Correlations.cxx b/PWGHF/correlationHF/AliAnalysisTaskSED0Correlations.cxx index bd33b81e593..3087627bf46 100644 --- a/PWGHF/correlationHF/AliAnalysisTaskSED0Correlations.cxx +++ b/PWGHF/correlationHF/AliAnalysisTaskSED0Correlations.cxx @@ -74,7 +74,11 @@ AliAnalysisTaskSE(), fNentries(0), fCutsD0(0), fCutsTracks(0), + fCorrelatorTr(0), + fCorrelatorKc(0), + fCorrelatorK0(0), fReadMC(0), + fMixing(kFALSE), fCounter(0), fNPtBins(1), fFillOnlyD0D0bar(0), @@ -102,7 +106,11 @@ AliAnalysisTaskSED0Correlations::AliAnalysisTaskSED0Correlations(const char *nam fNentries(0), fCutsD0(0), fCutsTracks(cutsTrk), + fCorrelatorTr(0), + fCorrelatorKc(0), + fCorrelatorK0(0), fReadMC(0), + fMixing(kFALSE), fCounter(0), fNPtBins(1), fFillOnlyD0D0bar(0), @@ -148,7 +156,11 @@ AliAnalysisTaskSED0Correlations::AliAnalysisTaskSED0Correlations(const AliAnalys fNentries(source.fNentries), fCutsD0(source.fCutsD0), fCutsTracks(source.fCutsTracks), + fCorrelatorTr(source.fCorrelatorTr), + fCorrelatorKc(source.fCorrelatorKc), + fCorrelatorK0(source.fCorrelatorK0), fReadMC(source.fReadMC), + fMixing(source.fMixing), fCounter(source.fCounter), fNPtBins(source.fNPtBins), fFillOnlyD0D0bar(source.fFillOnlyD0D0bar), @@ -183,7 +195,19 @@ AliAnalysisTaskSED0Correlations::~AliAnalysisTaskSED0Correlations() delete fNentries; fNentries = 0; } - if(fCounter){ + if (fCorrelatorTr) { + delete fCorrelatorTr; + fCorrelatorTr = 0; + } + if (fCorrelatorKc) { + delete fCorrelatorKc; + fCorrelatorKc = 0; + } + if (fCorrelatorK0) { + delete fCorrelatorK0; + fCorrelatorK0 = 0; + } + if (fCounter){ delete fCounter; fCounter=0; } @@ -208,7 +232,11 @@ AliAnalysisTaskSED0Correlations& AliAnalysisTaskSED0Correlations::operator=(cons fNentries = orig.fNentries; fCutsD0 = orig.fCutsD0; fCutsTracks = orig.fCutsTracks; + fCorrelatorTr = orig.fCorrelatorTr; + fCorrelatorKc = orig.fCorrelatorKc; + fCorrelatorK0 = orig.fCorrelatorK0; fReadMC = orig.fReadMC; + fMixing = orig.fMixing; fCounter = orig.fCounter; fNPtBins = orig.fNPtBins; fFillOnlyD0D0bar = orig.fFillOnlyD0D0bar; @@ -232,6 +260,7 @@ void AliAnalysisTaskSED0Correlations::Init() const char* nameoutput=GetOutputSlot(3)->GetContainer()->GetName(); copyfCutsD0->SetName(nameoutput); + //needer to clear completely the objects inside with Clear() method // Post the data PostData(3,copyfCutsD0); PostData(7,fCutsTracks); @@ -247,6 +276,33 @@ void AliAnalysisTaskSED0Correlations::UserCreateOutputObjects() // if(fDebug > 1) printf("AnalysisTaskSED0Correlations::UserCreateOutputObjects() \n"); + //HFCorrelator creation and definition + fCorrelatorTr = new AliHFCorrelator("CorrelatorTr",fCutsTracks,fSys); + fCorrelatorKc = new AliHFCorrelator("CorrelatorKc",fCutsTracks,fSys); + fCorrelatorK0 = new AliHFCorrelator("CorrelatorK0",fCutsTracks,fSys); + fCorrelatorTr->SetDeltaPhiInterval(-1.57,4.71);// set the Delta Phi Interval you want (in this case -0.5Pi to 1.5 Pi) + fCorrelatorKc->SetDeltaPhiInterval(-1.57,4.71); + fCorrelatorK0->SetDeltaPhiInterval(-1.57,4.71); + fCorrelatorTr->SetEventMixing(fMixing);// sets the analysis on a single event (kFALSE) or mixed events (kTRUE) + fCorrelatorKc->SetEventMixing(fMixing); + fCorrelatorK0->SetEventMixing(fMixing); + fCorrelatorTr->SetAssociatedParticleType(1);// set 1 for correlations with hadrons, 2 with kaons, 3 with KZeros + fCorrelatorKc->SetAssociatedParticleType(2);// set 1 for correlations with hadrons, 2 with kaons, 3 with KZeros + fCorrelatorK0->SetAssociatedParticleType(3);// set 1 for correlations with hadrons, 2 with kaons, 3 with KZeros + fCorrelatorTr->SetApplyDisplacementCut(2); //0: don't calculate d0; 1: return d0; 2: return d0/d0err + fCorrelatorKc->SetApplyDisplacementCut(2); + fCorrelatorK0->SetApplyDisplacementCut(0); + fCorrelatorTr->SetUseMC(fReadMC);// sets Montecarlo flag + fCorrelatorKc->SetUseMC(fReadMC); + fCorrelatorK0->SetUseMC(fReadMC); + fCorrelatorKc->SetPIDmode(2); //switch for K+/- PID option + Bool_t pooldefTr = fCorrelatorTr->DefineEventPool();// method that defines the properties ot the event mixing (zVtx and Multipl. bins) + Bool_t pooldefKc = fCorrelatorKc->DefineEventPool();// method that defines the properties ot the event mixing (zVtx and Multipl. bins) + Bool_t pooldefK0 = fCorrelatorK0->DefineEventPool();// method that defines the properties ot the event mixing (zVtx and Multipl. bins) + if(!pooldefTr) AliInfo("Warning:: Event pool not defined properly"); + if(!pooldefKc) AliInfo("Warning:: Event pool not defined properly"); + if(!pooldefK0) AliInfo("Warning:: Event pool not defined properly"); + // Several histograms are more conveniently managed in a TList fOutputMass = new TList(); fOutputMass->SetOwner(); @@ -405,7 +461,7 @@ void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/) return; } } - + //histogram filled with 1 for every AOD fNentries->Fill(0); fCounter->StoreEvent(aod,fCutsD0,fReadMC); @@ -439,9 +495,25 @@ void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/) if (skipEvent) return; } + //HFCorrelators initialization (for this event) + fCorrelatorTr->SetAODEvent(aod); // set the AOD event from which you are processing + fCorrelatorKc->SetAODEvent(aod); + fCorrelatorK0->SetAODEvent(aod); + Bool_t correlatorONTr = fCorrelatorTr->Initialize(); // initialize the pool for event mixing + Bool_t correlatorONKc = fCorrelatorKc->Initialize(); + Bool_t correlatorONK0 = fCorrelatorK0->Initialize(); + if(!correlatorONTr) {AliInfo("AliHFCorrelator (tracks) didn't initialize the pool correctly or processed a bad event"); return;} + if(!correlatorONKc) {AliInfo("AliHFCorrelator (charged K) didn't initialize the pool correctly or processed a bad event"); return;} + if(!correlatorONK0) {AliInfo("AliHFCorrelator (K0) didn't initialize the pool correctly or processed a bad event"); return;} + + if(fReadMC) { + fCorrelatorTr->SetMCArray(mcArray); // set the TClonesArray *fmcArray for analysis on monte carlo + fCorrelatorKc->SetMCArray(mcArray); + fCorrelatorK0->SetMCArray(mcArray); + } + // AOD primary vertex AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex(); - Bool_t isGoodVtx=kFALSE; //vtx1->Print(); @@ -454,7 +526,7 @@ void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/) //Reset flag for tracks distributions fill fAlreadyFilled=kFALSE; - // loop over candidates + //***** Loop over D0 candidates ***** Int_t nInD0toKpi = inputArray->GetEntriesFast(); if(fDebug>2) printf("Number of D0->Kpi: %d\n",nInD0toKpi); @@ -467,7 +539,7 @@ void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/) for(int j=0; j<2; j++) {idArrayV0[iV0][j]=-2;} AliAODv0 *v0 = (AliAODv0*)v0array->UncheckedAt(iV0); if(SelectV0(v0,vtx1,2,idArrayV0)) { //option 2 = for mass inv plots only - if(fReadMC && (v0->MatchToMC(311,mcArray,2,pdgCodes)<0)) continue; + if(fReadMC && (v0->MatchToMC(310,mcArray,2,pdgCodes)<0)) continue; //310 = K0s, 311 = K0 generico!! ((TH2F*)fOutputStudy->FindObject("hK0MassInv"))->Fill(v0->MassK0Short(),v0->Pt()); //invariant mass plot ((TH1F*)fOutputStudy->FindObject("hist_Pt_K0_AllEv"))->Fill(v0->Pt()); //pT distribution (in all events), K0 case } @@ -478,7 +550,7 @@ void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/) //rejection of tracks if(track->GetID() < 0) continue; //discard negative ID tracks if(track->Pt() < fPtThreshLow.at(0) || track->Pt() > fPtThreshUp.at(0)) continue; //discard tracks outside pt range for hadrons/K - if(!fCutsTracks->IsHadronSelected(track,vtx1,aod->GetMagneticField())) continue; + if(!fCutsTracks->IsHadronSelected(track) || !fCutsTracks->CheckHadronKinematic(track->Pt(),0.1)) continue; //0.1 = dummy (d0 value, no cut on it for me) //pT distribution (in all events), charg and hadr cases ((TH1F*)fOutputStudy->FindObject("hist_Pt_Charg_AllEv"))->Fill(track->Pt()); if(fCutsTracks->CheckKaonCompatibility(track,kFALSE,0,2)) ((TH1F*)fOutputStudy->FindObject("hist_Pt_Kcharg_AllEv"))->Fill(track->Pt()); @@ -504,20 +576,42 @@ void AliAnalysisTaskSED0Correlations::UserExec(Option_t */*option*/) Int_t ptbin=fCutsD0->PtBin(d->Pt()); if(ptbin==-1) {fNentries->Fill(4); continue;} //out of bounds - fIsSelectedCandidate=fCutsD0->IsSelected(d,AliRDHFCuts::kAll,aod); //selected + fIsSelectedCandidate=fCutsD0->IsSelected(d,AliRDHFCuts::kAll,aod); //D0 selected if(!fIsSelectedCandidate) continue; - if(!fReadMC) CalculateCorrelations(aod,d); //correlations on real data - else { //correlations on MC -> association of selected D0 to MCinfo with MCtruth + //D0 infos + Double_t phiD0 = fCorrelatorTr->SetCorrectPhiRange(d->Phi()); + phiD0 = fCorrelatorKc->SetCorrectPhiRange(d->Phi()); //bad usage, but returns a Double_t... + phiD0 = fCorrelatorK0->SetCorrectPhiRange(d->Phi()); + fCorrelatorTr->SetTriggerParticleProperties(d->Pt(),phiD0,d->Eta()); // sets the parameters of the trigger particles that are needed + fCorrelatorKc->SetTriggerParticleProperties(d->Pt(),phiD0,d->Eta()); + fCorrelatorK0->SetTriggerParticleProperties(d->Pt(),phiD0,d->Eta()); + fCorrelatorTr->SetD0Properties(d,fIsSelectedCandidate); //sets special properties for D0 + fCorrelatorKc->SetD0Properties(d,fIsSelectedCandidate); + fCorrelatorK0->SetD0Properties(d,fIsSelectedCandidate); + + if(!fReadMC) { + CalculateCorrelations(d); //correlations on real data + } else { //correlations on MC -> association of selected D0 to MCinfo with MCtruth Int_t pdgDgD0toKpi[2]={321,211}; Int_t labD0 = d->MatchToMC(421,mcArray,2,pdgDgD0toKpi); //return MC particle label if the array corresponds to a D0, -1 if not - if (labD0>-1) CalculateCorrelations(aod,d,labD0,mcArray); + if (labD0>-1) { + CalculateCorrelations(d,labD0,mcArray); + } } FillMassHists(d,mcArray,fCutsD0,fOutputMass); } } //end for prongs + + if(fMixing /* && fAlreadyFilled*/) { // update the pool for Event Mixing, if: enabled, event is ok, at least a SelD0 found! (fAlreadyFilled's role!) + Bool_t updatedTr = fCorrelatorTr->PoolUpdate(); + Bool_t updatedKc = fCorrelatorKc->PoolUpdate(); + Bool_t updatedK0 = fCorrelatorK0->PoolUpdate(); + if(!updatedTr || !updatedKc || !updatedK0) AliInfo("Pool was not updated"); + } + fCounter->StoreCandidates(aod,nSelectedloose,kTRUE); fCounter->StoreCandidates(aod,nSelectedtight,kFALSE); @@ -666,7 +760,7 @@ Int_t AliAnalysisTaskSED0Correlations::CheckD0Origin(TClonesArray* arrayMC, AliA // // checking whether the mother of the particles come from a charm or a bottom quark // - printf(" AliAnalysisTaskSED0Correlations::CheckD0Origin() \n"); + printf("AliAnalysisTaskSED0Correlations::CheckD0Origin() \n"); Int_t pdgGranma = 0; Int_t mother = 0; @@ -705,276 +799,201 @@ void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() { TString namePlot = ""; //These for limits in THnSparse (one per bin, same limits). - //Vars: DeltaPhi, InvMass, PtTrack + //Vars: DeltaPhi, InvMass, PtTrack, Displacement, DeltaEta Int_t nBinsPhi[5] = {32,150,30,3,16}; Double_t binMinPhi[5] = {-1.6,1.6,0.,0.,-1.6}; //is the minimum for all the bins Double_t binMaxPhi[5] = {4.8,2.2,3.0,3.,1.6}; //is the maximum for all the bins - for(Int_t i=0;iSumw2(); - fOutputCorr->Add(hPhiK); - - namePlot="hPhi_Kcharg_Bin"; - namePlot+=i; - - THnSparseI *hPhiH = new THnSparseI(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiH->Sumw2(); - fOutputCorr->Add(hPhiH); - - namePlot="hPhi_Charg_Bin"; - namePlot+=i; - - THnSparseI *hPhiC = new THnSparseI(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiC->Sumw2(); - fOutputCorr->Add(hPhiC); - - //histos for c/b origin for D0 (MC only) - if (fReadMC) { + for(Int_t i=0;iSumw2(); - fOutputCorr->Add(hPhiK_c); + THnSparseI *hPhiK = new THnSparseI(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiK->Sumw2(); + fOutputCorr->Add(hPhiK); - namePlot="hPhi_Kcharg_From_c_Bin"; + namePlot="hPhi_Kcharg_Bin"; namePlot+=i; - THnSparseI *hPhiH_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiH_c->Sumw2(); - fOutputCorr->Add(hPhiH_c); + THnSparseI *hPhiH = new THnSparseI(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiH->Sumw2(); + fOutputCorr->Add(hPhiH); - namePlot="hPhi_Charg_From_c_Bin"; + namePlot="hPhi_Charg_Bin"; namePlot+=i; - THnSparseI *hPhiC_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiC_c->Sumw2(); - fOutputCorr->Add(hPhiC_c); + THnSparseI *hPhiC = new THnSparseI(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiC->Sumw2(); + fOutputCorr->Add(hPhiC); - namePlot="hPhi_K0_From_b_Bin"; - namePlot+=i; - - THnSparseI *hPhiK_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiK_b->Sumw2(); - fOutputCorr->Add(hPhiK_b); - - namePlot="hPhi_Kcharg_From_b_Bin"; - namePlot+=i; - - THnSparseI *hPhiH_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiH_b->Sumw2(); - fOutputCorr->Add(hPhiH_b); - - namePlot="hPhi_Charg_From_b_Bin"; - namePlot+=i; - - THnSparseI *hPhiC_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiC_b->Sumw2(); - fOutputCorr->Add(hPhiC_b); - - //HF-only tracks (c for c->D0, b for b->D0) - namePlot="hPhi_K0_HF_From_c_Bin"; - namePlot+=i; - - THnSparseI *hPhiK_HF_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiK_HF_c->Sumw2(); - fOutputCorr->Add(hPhiK_HF_c); - - namePlot="hPhi_Kcharg_HF_From_c_Bin"; - namePlot+=i; - - THnSparseI *hPhiH_HF_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiH_HF_c->Sumw2(); - fOutputCorr->Add(hPhiH_HF_c); - - namePlot="hPhi_Charg_HF_From_c_Bin"; - namePlot+=i; - THnSparseI *hPhiC_HF_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiC_HF_c->Sumw2(); - fOutputCorr->Add(hPhiC_HF_c); - - namePlot="hPhi_K0_HF_From_b_Bin"; - namePlot+=i; - - THnSparseI *hPhiK_HF_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiK_HF_b->Sumw2(); - fOutputCorr->Add(hPhiK_HF_b); - - namePlot="hPhi_Kcharg_HF_From_b_Bin"; - namePlot+=i; - - THnSparseI *hPhiH_HF_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiH_HF_b->Sumw2(); - fOutputCorr->Add(hPhiH_HF_b); - - namePlot="hPhi_Charg_HF_From_b_Bin"; - namePlot+=i; - - THnSparseI *hPhiC_HF_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); - hPhiC_HF_b->Sumw2(); - fOutputCorr->Add(hPhiC_HF_b); - } - - //leading hadron correlations - namePlot="hPhi_Lead_Bin"; - namePlot+=i; - - TH2F *hCorrLead = new TH2F(namePlot.Data(), "Leading particle correlation; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - hCorrLead->Sumw2(); - fOutputCorr->Add(hCorrLead); - - if (fReadMC) { - namePlot="hPhi_Lead_From_c_Bin"; - namePlot+=i; - - TH2F *hCorrLead_c = new TH2F(namePlot.Data(), "Leading particle correlation - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - hCorrLead_c->Sumw2(); - fOutputCorr->Add(hCorrLead_c); - - namePlot="hPhi_Lead_From_b_Bin"; - namePlot+=i; - - TH2F *hCorrLead_b = new TH2F(namePlot.Data(), "Leading particle correlation - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - hCorrLead_b->Sumw2(); - fOutputCorr->Add(hCorrLead_b); - - namePlot="hPhi_Lead_HF_From_c_Bin"; - namePlot+=i; - - TH2F *hCorrLead_HF_c = new TH2F(namePlot.Data(), "Leading particle correlation HF - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - hCorrLead_HF_c->Sumw2(); - fOutputCorr->Add(hCorrLead_HF_c); - - namePlot="hPhi_Lead_HF_From_b_Bin"; - namePlot+=i; - - TH2F *hCorrLead_HF_b = new TH2F(namePlot.Data(), "Leading particle correlation HF - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - hCorrLead_HF_b->Sumw2(); - fOutputCorr->Add(hCorrLead_HF_b); - } - - //pT weighted correlations - namePlot="hPhi_Weig_Bin"; - namePlot+=i; + //histos for c/b origin for D0 (MC only) + if (fReadMC) { - TH2F *hCorrWeig = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted); #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeig); + //generic origin for tracks + namePlot="hPhi_K0_From_c_Bin"; + namePlot+=i; - if (fReadMC) { - namePlot="hPhi_Weig_From_c_Bin"; - namePlot+=i; + THnSparseI *hPhiK_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiK_c->Sumw2(); + fOutputCorr->Add(hPhiK_c); - TH2F *hCorrWeig_c = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted) - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeig_c); + namePlot="hPhi_Kcharg_From_c_Bin"; + namePlot+=i; - namePlot="hPhi_Weig_From_b_Bin"; - namePlot+=i; + THnSparseI *hPhiH_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiH_c->Sumw2(); + fOutputCorr->Add(hPhiH_c); - TH2F *hCorrWeig_b = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted) - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeig_b); + namePlot="hPhi_Charg_From_c_Bin"; + namePlot+=i; - namePlot="hPhi_Weig_HF_From_c_Bin"; - namePlot+=i; + THnSparseI *hPhiC_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiC_c->Sumw2(); + fOutputCorr->Add(hPhiC_c); + + namePlot="hPhi_K0_From_b_Bin"; + namePlot+=i; - TH2F *hCorrWeig_HF_c = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted) HF - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeig_HF_c); + THnSparseI *hPhiK_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiK_b->Sumw2(); + fOutputCorr->Add(hPhiK_b); - namePlot="hPhi_Weig_HF_From_b_Bin"; - namePlot+=i; + namePlot="hPhi_Kcharg_From_b_Bin"; + namePlot+=i; - TH2F *hCorrWeig_HF_b = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted) HF - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeig_HF_b); - } + THnSparseI *hPhiH_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiH_b->Sumw2(); + fOutputCorr->Add(hPhiH_b); - //pT weighted correlations (squared weights) - namePlot="hPhi_WeigSqr_Bin"; - namePlot+=i; + namePlot="hPhi_Charg_From_b_Bin"; + namePlot+=i; - TH2F *hCorrWeigSqr = new TH2F(namePlot.Data(), "Charged particle correlation (pT Weighted - squared); #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeigSqr); + THnSparseI *hPhiC_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiC_b->Sumw2(); + fOutputCorr->Add(hPhiC_b); - if (fReadMC) { - namePlot="hPhi_WeigSqr_From_c_Bin"; - namePlot+=i; + //HF-only tracks (c for c->D0, b for b->D0) + namePlot="hPhi_K0_HF_From_c_Bin"; + namePlot+=i; - TH2F *hCorrWeigSqr_c = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted - squared) - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeigSqr_c); + THnSparseI *hPhiK_HF_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiK_HF_c->Sumw2(); + fOutputCorr->Add(hPhiK_HF_c); - namePlot="hPhi_WeigSqr_From_b_Bin"; - namePlot+=i; + namePlot="hPhi_Kcharg_HF_From_c_Bin"; + namePlot+=i; - TH2F *hCorrWeigSqr_b = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted - squared) - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeigSqr_b); + THnSparseI *hPhiH_HF_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiH_HF_c->Sumw2(); + fOutputCorr->Add(hPhiH_HF_c); - namePlot="hPhi_WeigSqr_HF_From_c_Bin"; - namePlot+=i; + namePlot="hPhi_Charg_HF_From_c_Bin"; + namePlot+=i; - TH2F *hCorrWeigSqr_HF_c = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted - squared) HF - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeigSqr_HF_c); + THnSparseI *hPhiC_HF_c = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiC_HF_c->Sumw2(); + fOutputCorr->Add(hPhiC_HF_c); - namePlot="hPhi_WeigSqr_HF_From_b_Bin"; - namePlot+=i; + namePlot="hPhi_K0_HF_From_b_Bin"; + namePlot+=i; - TH2F *hCorrWeigSqr_HF_b = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted - squared) HF - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeigSqr_HF_b); - } + THnSparseI *hPhiK_HF_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiK_HF_b->Sumw2(); + fOutputCorr->Add(hPhiK_HF_b); - //pT weighted correlations (inverse of pT distribution weights) - namePlot="hPhi_WeigDist_Bin"; - namePlot+=i; + namePlot="hPhi_Kcharg_HF_From_b_Bin"; + namePlot+=i; - TH2F *hCorrWeigDist = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted - inv.dis.); #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeigDist); + THnSparseI *hPhiH_HF_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiH_HF_b->Sumw2(); + fOutputCorr->Add(hPhiH_HF_b); - if (fReadMC) { - namePlot="hPhi_WeigDist_From_c_Bin"; - namePlot+=i; + namePlot="hPhi_Charg_HF_From_b_Bin"; + namePlot+=i; - TH2F *hCorrWeigDist_c = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted - inv.dis.) - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeigDist_c); + THnSparseI *hPhiC_HF_b = new THnSparseI(namePlot.Data(), "Azimuthal correlation HF - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",5,nBinsPhi,binMinPhi,binMaxPhi); + hPhiC_HF_b->Sumw2(); + fOutputCorr->Add(hPhiC_HF_b); + } - namePlot="hPhi_WeigDist_From_b_Bin"; + //leading hadron correlations + namePlot="hPhi_Lead_Bin"; namePlot+=i; - TH2F *hCorrWeigDist_b = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted - inv.dis.) - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeigDist_b); + TH2F *hCorrLead = new TH2F(namePlot.Data(), "Leading particle correlation; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); + hCorrLead->Sumw2(); + fOutputCorr->Add(hCorrLead); - namePlot="hPhi_WeigDist_HF_From_c_Bin"; - namePlot+=i; - - TH2F *hCorrWeigDist_HF_c = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted - inv.dis.) HF - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeigDist_HF_c); + if (fReadMC) { + namePlot="hPhi_Lead_From_c_Bin"; + namePlot+=i; - namePlot="hPhi_WeigDist_HF_From_b_Bin"; + TH2F *hCorrLead_c = new TH2F(namePlot.Data(), "Leading particle correlation - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); + hCorrLead_c->Sumw2(); + fOutputCorr->Add(hCorrLead_c); + + namePlot="hPhi_Lead_From_b_Bin"; + namePlot+=i; + + TH2F *hCorrLead_b = new TH2F(namePlot.Data(), "Leading particle correlation - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); + hCorrLead_b->Sumw2(); + fOutputCorr->Add(hCorrLead_b); + + namePlot="hPhi_Lead_HF_From_c_Bin"; + namePlot+=i; + + TH2F *hCorrLead_HF_c = new TH2F(namePlot.Data(), "Leading particle correlation HF - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); + hCorrLead_HF_c->Sumw2(); + fOutputCorr->Add(hCorrLead_HF_c); + + namePlot="hPhi_Lead_HF_From_b_Bin"; + namePlot+=i; + + TH2F *hCorrLead_HF_b = new TH2F(namePlot.Data(), "Leading particle correlation HF - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); + hCorrLead_HF_b->Sumw2(); + fOutputCorr->Add(hCorrLead_HF_b); + } + + //pT weighted correlations + namePlot="hPhi_Weig_Bin"; namePlot+=i; - - TH2F *hCorrWeigDist_HF_b = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted - inv.dis.) HF - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); - fOutputCorr->Add(hCorrWeigDist_HF_b); - } - - //counters - namePlot = "hist_Count_Charg_Bin"; namePlot+=i; - TH1F *hCountC = new TH1F(namePlot.Data(), "Charged track counter; # Tracks",100,0.,100.); - hCountC->SetMinimum(0); - fOutputStudy->Add(hCountC); - - namePlot = "hist_Count_Kcharg_Bin"; namePlot+=i; - TH1F *hCountH = new TH1F(namePlot.Data(), "Hadrons counter; # Tracks",100,0.,100.); - hCountH->SetMinimum(0); - fOutputStudy->Add(hCountH); - - namePlot = "hist_Count_K0_Bin"; namePlot+=i; - TH1F *hCountK = new TH1F(namePlot.Data(), "Kaons counter; # Tracks",10,0.,10.); - hCountK->SetMinimum(0); - fOutputStudy->Add(hCountK); + + TH2F *hCorrWeig = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted); #Delta#phi",32,-1.6,4.8,300,1.6,2.2); + fOutputCorr->Add(hCorrWeig); + + if (fReadMC) { + namePlot="hPhi_Weig_From_c_Bin"; + namePlot+=i; + + TH2F *hCorrWeig_c = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted) - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); + fOutputCorr->Add(hCorrWeig_c); + + namePlot="hPhi_Weig_From_b_Bin"; + namePlot+=i; + + TH2F *hCorrWeig_b = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted) - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); + fOutputCorr->Add(hCorrWeig_b); + + namePlot="hPhi_Weig_HF_From_c_Bin"; + namePlot+=i; + + TH2F *hCorrWeig_HF_c = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted) HF - c origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); + fOutputCorr->Add(hCorrWeig_HF_c); + + namePlot="hPhi_Weig_HF_From_b_Bin"; + namePlot+=i; + + TH2F *hCorrWeig_HF_b = new TH2F(namePlot.Data(), "Charged particle correlation (pT weighted) HF - b origin; #Delta#phi",32,-1.6,4.8,300,1.6,2.2); + fOutputCorr->Add(hCorrWeig_HF_b); + } //pT distribution histos namePlot = "hist_Pt_Charg_Bin"; namePlot+=i; @@ -998,7 +1017,33 @@ void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() { hDstarPions->GetXaxis()->SetBinLabel(1,"Not rejected"); hDstarPions->GetXaxis()->SetBinLabel(2,"Rejected"); hDstarPions->SetMinimum(0); - fOutputStudy->Add(hDstarPions); + fOutputStudy->Add(hDstarPions); + + } + + if(fMixing) { + //THnSparse plots for event mixing! + namePlot="hPhi_K0_Bin"; + namePlot+=i;namePlot+="_EvMix"; + + THnSparseI *hPhiK_EvMix = new THnSparseI(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + hPhiK_EvMix->Sumw2(); + fOutputCorr->Add(hPhiK_EvMix); + + namePlot="hPhi_Kcharg_Bin"; + namePlot+=i;namePlot+="_EvMix"; + + THnSparseI *hPhiH_EvMix = new THnSparseI(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + hPhiH_EvMix->Sumw2(); + fOutputCorr->Add(hPhiH_EvMix); + + namePlot="hPhi_Charg_Bin"; + namePlot+=i;namePlot+="_EvMix"; + + THnSparseI *hPhiC_EvMix = new THnSparseI(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix); + hPhiC_EvMix->Sumw2(); + fOutputCorr->Add(hPhiC_EvMix); + } } //out of bin loop @@ -1008,6 +1053,18 @@ void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() { hRejTracks->GetXaxis()->SetBinLabel(2,"Rejected"); fOutputStudy->Add(hRejTracks); + TH1F *hCountC = new TH1F("hist_Count_Charg", "Charged track counter; # Tracks",100,0.,100.); + hCountC->SetMinimum(0); + fOutputStudy->Add(hCountC); + + TH1F *hCountH = new TH1F("hist_Count_Kcharg", "Hadrons counter; # Tracks",20,0.,20.); + hCountH->SetMinimum(0); + fOutputStudy->Add(hCountH); + + TH1F *hCountK = new TH1F("hist_Count_K0", "Kaons counter; # Tracks",20,0.,20.); + hCountK->SetMinimum(0); + fOutputStudy->Add(hCountK); + if (fFillGlobal) { //all-events plots //pt distributions TH1F *hPtCAll = new TH1F("hist_Pt_Charg_AllEv", "Charged track pT (All); p_{T} (GeV/c)",240,0.,12.); @@ -1018,10 +1075,27 @@ void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() { hPtHAll->SetMinimum(0); fOutputStudy->Add(hPtHAll); - TH1F *hPtKAll = new TH1F("hist_Pt_K0_AllEv", "Kaons pT (All); p_{T} (GeV/c)",240,0.,12.); + TH1F *hPtKAll = new TH1F("hist_Pt_K0_AllEv", "Kaons pT (All); p_{T} (GeV/c)",240,0.,12.); hPtKAll->SetMinimum(0); fOutputStudy->Add(hPtKAll); + //phi distributions + TH1F *hPhiDistCAll = new TH1F("hist_PhiDistr_Charg", "Charged track phi distr. (All); p_{T} (GeV/c)",32,-1.6,4.8); + hPhiDistCAll->SetMinimum(0); + fOutputStudy->Add(hPhiDistCAll); + + TH1F *hPhiDistHAll = new TH1F("hist_PhiDistr_Kcharg", "Hadrons phi distr. (All); p_{T} (GeV/c)",32,-1.6,4.8); + hPhiDistHAll->SetMinimum(0); + fOutputStudy->Add(hPhiDistHAll); + + TH1F *hPhiDistKAll = new TH1F("hist_PhiDistr_K0", "Kaons phi distr. (All); p_{T} (GeV/c)",32,-1.6,4.8); + hPhiDistKAll->SetMinimum(0); + fOutputStudy->Add(hPhiDistKAll); + + TH1F *hPhiDistDAll = new TH1F("hist_PhiDistr_D0", "D^{0} phi distr. (All); p_{T} (GeV/c)",32,-1.6,4.8); + hPhiDistDAll->SetMinimum(0); + fOutputStudy->Add(hPhiDistDAll); + //K0 Invariant Mass plots TH2F *hK0MassInv = new TH2F("hK0MassInv", "K0 invariant mass; Invariant mass (MeV/c^{2}); pT (GeV/c)",200,0.4,0.6,100,0.,10.); hK0MassInv->SetMinimum(0); @@ -1178,345 +1252,317 @@ void AliAnalysisTaskSED0Correlations::CreateCorrelationsObjs() { } } + if(fMixing) { //check for event mixing (events in bins) + TH2I* EventMixingCheckTr = new TH2I("EventMixingCheckTr","EventMixingCheckTr",5,-0.5,4.5,7,-0.5,6.5); + fOutputStudy->Add(EventMixingCheckTr); + + TH2I* EventMixingCheckK0 = new TH2I("EventMixingCheckK0","EventMixingCheckK0",5,-0.5,4.5,7,-0.5,6.5); + fOutputStudy->Add(EventMixingCheckK0); + } + } //________________________________________________________________________ -void AliAnalysisTaskSED0Correlations::CalculateCorrelations(AliAODEvent* aod, AliAODRecoDecayHF2Prong* d, Int_t labD0, TClonesArray* mcArray) { +void AliAnalysisTaskSED0Correlations::CalculateCorrelations(AliAODRecoDecayHF2Prong* d, Int_t labD0, TClonesArray* mcArray) { // // Method for correlations D0-hadrons study // - - Double_t mD0, mD0bar, deltaphi, d0, d0err, deltaeta; + Int_t N_Charg = 0, N_KCharg = 0, N_Kaons = 0; + Double_t mD0, mD0bar; + Int_t origD0 = 0, PDGD0 = 0; d->InvMassD0(mD0,mD0bar); Int_t ptbin = PtBinCorr(d->Pt()); if(ptbin < 0) return; - AliAODVertex *vtx = (AliAODVertex*)aod->GetPrimaryVertex(); - Int_t N_Charg = 0, N_Kcharg = 0, N_Kaons = 0; + //Origin of D0 + TString orig=""; + if(fReadMC) { + origD0=CheckD0Origin(mcArray,(AliAODMCParticle*)mcArray->At(labD0)); + PDGD0 = ((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode(); + switch (CheckD0Origin(mcArray,(AliAODMCParticle*)mcArray->At(labD0))) { + case 4: + orig = "_From_c"; + ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(0.); + break; + case 5: + orig = "_From_b"; + ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(1.); + break; + default: + return; + } + } + + Double_t highPt = 0; Double_t lead[3] = {0,0,0}; //infos for leading particle (pt,deltaphi) - if(fReadMC == 0) { - Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped - Double_t highPt = 0; Double_t lead[2] = {0,0}; //infos for leading particle (pt,deltaphi) + //loop over the tracks in the pool + Bool_t execPool = fCorrelatorTr->ProcessEventPool(); //pool is ready? (only in ME, in SE returns kFALSE) + + if(fMixing && !execPool) { + AliInfo("Mixed event analysis: pool is not ready"); + return; + } + + Int_t NofEventsinPool = 1; + if(fMixing) NofEventsinPool = fCorrelatorTr->GetNofEventsInPool(); + + for (Int_t jMix =0; jMix < NofEventsinPool; jMix++) {// loop on events in the pool; if it is SE analysis, stops at one (index not needed there) + Bool_t analyzetracksTr = fCorrelatorTr->ProcessAssociatedTracks(jMix);// process all the tracks in the aodEvent, by applying the selection cuts + Bool_t analyzetracksKc = fCorrelatorKc->ProcessAssociatedTracks(jMix); + Bool_t analyzetracksK0 = fCorrelatorK0->ProcessAssociatedTracks(jMix); + if(!analyzetracksTr || !analyzetracksKc || !analyzetracksK0) { + AliInfo("AliHFCorrelator::Cannot process the track array"); + continue; + } - for(Int_t itrack=0; itrackGetNTracks(); itrack++) { // loop on tracks - AliAODTrack * track = aod->GetTrack(itrack); - - if(!TrackSelectedInLoop(track,d,aod,ptbin,idDaughs)) continue; //rejection of track (daughter of D0/quality cuts not passed/soft pion/negative ID - - GetImpParameter(track,aod,d0,d0err); //evaluates d0 and sigma_d0 - - ((TH1F*)fOutputStudy->FindObject("hRejTracks"))->Fill(0); //track accepted by quality cuts - deltaphi = d->Phi()-track->Phi(); - deltaeta = d->Eta()-track->Eta(); - if (deltaphi < -1.571) deltaphi+=6.283; - if (deltaphi > 4.712) deltaphi-=6.283; - Double_t pttrack = track->Pt(); - - if (pttrack > highPt) {highPt = pttrack; lead[0] = pttrack; lead[1] = deltaphi;} //for leading particle - - //Lines needed to include overflow into THnSparse projections! - Double_t ptLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(2)->GetXmax(); //all plots have same axes... - Double_t displLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(3)->GetXmax(); - Double_t EtaLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(4)->GetXmax(); - if(pttrack > ptLim_Sparse) pttrack = ptLim_Sparse-0.01; - if(d0/d0err > displLim_Sparse) d0 = (displLim_Sparse-0.001)*d0err; - if(deltaeta > EtaLim_Sparse) deltaeta = EtaLim_Sparse-0.01; - if(deltaeta < -EtaLim_Sparse) deltaeta = -EtaLim_Sparse+0.01; - - //variables for filling histos - Double_t fillSpPhiD0[5] = {deltaphi,mD0,pttrack,d0/d0err,deltaeta}; - Double_t fillSpPhiD0bar[5] = {deltaphi,mD0bar,pttrack,d0/d0err,deltaeta}; - - //generic charged tracks (NO PID selection) - if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) { //D0 - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->Fill(fillSpPhiD0); + //Charged tracks + for(Int_t iTrack = 0; iTrackGetNofTracks(); iTrack++){ // looping on track candidates + + Bool_t runcorrelation = fCorrelatorTr->Correlate(iTrack); + if(!runcorrelation) continue; + + AliReducedParticle* track = fCorrelatorTr->GetAssociatedParticle(); + + if(!fMixing) { + if(!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; + } else { //not a soft pion + if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0); + if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0bar); + } } - if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) { //D0bar - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->Fill(fillSpPhiD0bar); + if(track->Pt() < fPtThreshLow.at(ptbin) || track->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K + + Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped + if(track->GetID() == idDaughs[0] || track->GetID() == idDaughs[1]) continue; //discards daughters of candidate + + FillSparsePlots(mcArray,d,origD0,PDGD0,track,kTrack); //fills for charged tracks + + if(!fMixing) N_Charg++; + + //retrieving leading info... + if(track->Pt() > highPt) { + lead[0] = fCorrelatorTr->GetDeltaPhi(); + lead[1] = fCorrelatorTr->GetDeltaEta(); + lead[2] = fReadMC ? CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(track->GetLabel())) : 0; + highPt = track->Pt(); } - if(!fAlreadyFilled) ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_Charg_Bin%d",ptbin)))->Fill(track->Pt()); - N_Charg++; - //pT_weighted track correlations fill - if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) { //D0 - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0,pttrack); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigSqr_Bin%d",ptbin)))->Fill(deltaphi,mD0,pow(pttrack,2.)); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigDist_Bin%d",ptbin)))->Fill(deltaphi,mD0,PtWeig(ptbin,pttrack)); - } - if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) { //D0bar - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0bar,pttrack); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigSqr_Bin%d",ptbin)))->Fill(deltaphi,mD0bar,pow(pttrack,2.)); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigDist_Bin%d",ptbin)))->Fill(deltaphi,mD0bar,PtWeig(ptbin,pttrack)); - } - - //Charged Kaon identification - if(fCutsTracks->CheckKaonCompatibility(track,kFALSE,0,2)) { - if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) { //D0 - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Kcharg_Bin%d",ptbin)))->Fill(fillSpPhiD0); - } - if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) { //D0bar - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Kcharg_Bin%d",ptbin)))->Fill(fillSpPhiD0bar); + + } // end of tracks loop + + //Charged Kaons loop + for(Int_t iTrack = 0; iTrackGetNofTracks(); iTrack++){ // looping on charged kaons candidates + + Bool_t runcorrelation = fCorrelatorKc->Correlate(iTrack); + if(!runcorrelation) continue; + + AliReducedParticle* kCharg = fCorrelatorKc->GetAssociatedParticle(); + + if(!fMixing) { + if(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; + } else { + if (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0); + if (fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0bar); } - if(!fAlreadyFilled) ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_Kcharg_Bin%d",ptbin)))->Fill(track->Pt()); - N_Kcharg++; - } //end hadron case + } + if(kCharg->Pt() < fPtThreshLow.at(ptbin) || kCharg->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K - } //end tracks loop + Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped + if(kCharg->GetID() == idDaughs[0] || kCharg->GetID() == idDaughs[1]) continue; //discards daughters of candidate + + FillSparsePlots(mcArray,d,origD0,PDGD0,kCharg,kKCharg); //fills for charged tracks - //kaon identification - TClonesArray *v0array = (TClonesArray*)aod->GetList()->FindObject("v0s"); - Int_t idArrayV0[v0array->GetEntriesFast()][2]; //array for skipping K0 double-counting - for(int iV0=0; iV0GetEntriesFast(); iV0++) { - for(int j=0; j<2; j++) {idArrayV0[iV0][j]=-2;} + if(!fMixing) N_KCharg++; - AliAODv0 *v0 = (AliAODv0*)v0array->UncheckedAt(iV0); - if(v0->Pt() < fPtThreshLow.at(ptbin) || v0->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K - if(!SelectV0(v0,vtx,1,idArrayV0)) continue; //option 1 = for correlations - - Double_t ptV0=v0->Pt(), deltaphiV0, deltaetaV0; - Double_t ptLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(2)->GetXmax(); //all plots have same axes... - Double_t EtaLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(4)->GetXmax(); - if (ptV0 > ptLim_Sparse) ptV0 = ptLim_Sparse-0.01; - deltaphiV0 = d->Phi()-v0->Phi(); - deltaetaV0 = d->Eta()-v0->Eta(); - if (deltaphiV0 < -1.571) deltaphiV0+=6.283; - if (deltaphiV0 > 4.712) deltaphiV0-=6.283; - if(deltaetaV0 > EtaLim_Sparse) deltaetaV0 = EtaLim_Sparse-0.01; - if(deltaetaV0 < -EtaLim_Sparse) deltaetaV0 = -EtaLim_Sparse+0.01; - - Double_t fillSpPhiD0K0[5] = {deltaphiV0,mD0,ptV0,0.,deltaetaV0}; - Double_t fillSpPhiD0barK0[5] = {deltaphiV0,mD0bar,ptV0,0.,deltaetaV0}; + } // end of charged kaons loop - if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) { //D0 - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_K0_Bin%d",ptbin)))->Fill(fillSpPhiD0K0); - } - if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) { //D0bar - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_K0_Bin%d",ptbin)))->Fill(fillSpPhiD0barK0); - } - if(!fAlreadyFilled) ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_K0_Bin%d",ptbin)))->Fill(v0->Pt()); - N_Kaons++; - } // end kaon case + //K0 loop + for(Int_t iTrack = 0; iTrackGetNofTracks(); iTrack++){ // looping on k0 candidates + + Bool_t runcorrelation = fCorrelatorK0->Correlate(iTrack); + if(!runcorrelation) continue; + + AliReducedParticle* k0 = fCorrelatorK0->GetAssociatedParticle(); - //Leading track correlations fill - if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) { //D0 - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[1],mD0); + if(k0->Pt() < fPtThreshLow.at(ptbin) || k0->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K + + Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped + if(k0->GetID() == idDaughs[0] || k0->GetID() == idDaughs[1]) continue; //discards daughters of candidate + + FillSparsePlots(mcArray,d,origD0,PDGD0,k0,kK0); //fills for charged tracks + + if(!fMixing) N_Kaons++; + + } // end of charged kaons loop + + } //end of events loop + + //leading track correlations fill + if(fReadMC) { + if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==421) { //D0 + ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[1],mD0); //c and b D0 + ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0); //c or b D0 + if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0); + if(origD0==5&&(int)lead[2]>=3&&(int)lead[2]<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0); } - if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) { //D0bar - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[1],mD0bar); + if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==-421) { //D0bar + ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[1],mD0bar); + ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0bar); //c or b D0 + if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0bar); + if(origD0==5&&(int)lead[2]>=3&&(int)lead[2]<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[0],mD0bar); } + } else { + if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[0],mD0); //c and b D0 + if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[0],mD0bar); + } - } //end fReadMC = 0 + //Fill of count histograms + if (!fAlreadyFilled) { + ((TH1F*)fOutputStudy->FindObject("hist_Count_Charg"))->Fill(N_Charg); + ((TH1F*)fOutputStudy->FindObject("hist_Count_Kcharg"))->Fill(N_KCharg); + ((TH1F*)fOutputStudy->FindObject("hist_Count_K0"))->Fill(N_Kaons); + } - if(fReadMC == 1) { - Int_t idDaughs[2] = {((AliVTrack*)d->GetDaughter(0))->GetID(),((AliVTrack*)d->GetDaughter(1))->GetID()}; //IDs of daughters to be skipped - Double_t highPt = 0; Double_t lead[3] = {0,0,0}; //infos for leading particle (pt,deltaphi,origtrack) + fAlreadyFilled=kTRUE; //at least a D0 analyzed in the event; distribution plots already filled - //Origin of D0 - TString orig=""; Int_t origD0=CheckD0Origin(mcArray,(AliAODMCParticle*)mcArray->At(labD0)); - switch (CheckD0Origin(mcArray,(AliAODMCParticle*)mcArray->At(labD0))) - { - case 4: - orig = "_From_c"; - ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(0.); - break; - case 5: - orig = "_From_b"; - ((TH1F*)fOutputStudy->FindObject(Form("histOrig_D0_Bin%d",ptbin)))->Fill(1.); - break; - default: - return; +} + +//________________________________________________________________________ +void AliAnalysisTaskSED0Correlations::FillSparsePlots(TClonesArray* mcArray, AliAODRecoDecayHF2Prong *d, Int_t origD0, Int_t PdgD0, AliReducedParticle* track, Int_t type) { + // + //fills the THnSparse for correlations, calculating the variables + // + + //Initialization of variables + Int_t ptbin = PtBinCorr(d->Pt()); + if(ptbin < 0) return; + Double_t mD0, mD0bar, deltaphi = 0., deltaeta = 0.; + d->InvMassD0(mD0,mD0bar); + + Double_t ptTrack = track->Pt(); + Double_t d0Track = type!=kK0 ? track->GetImpPar() : 0.; + Double_t phiTr = track->Phi(); + Double_t origTr = fReadMC ? CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(track->GetLabel())) : 0; + + TString part = "", orig = ""; + + switch (type) { + case(kTrack): { + part = "Charg"; + deltaphi = fCorrelatorTr->GetDeltaPhi(); + deltaeta = fCorrelatorTr->GetDeltaEta(); + break; + } + case(kKCharg): { + part = "Kcharg"; + deltaphi = fCorrelatorKc->GetDeltaPhi(); + deltaeta = fCorrelatorKc->GetDeltaEta(); + break; + } + case(kK0): { + part = "K0"; + deltaphi = fCorrelatorK0->GetDeltaPhi(); + deltaeta = fCorrelatorK0->GetDeltaEta(); + break; } + } - for(Int_t itrack=0; itrackGetNTracks(); itrack++) { // loop on tracks - AliAODTrack * track = aod->GetTrack(itrack); + //Fixes limits; needed to include overflow into THnSparse projections! + Double_t pTorig = track->Pt(); + Double_t d0orig = track->GetImpPar(); + Double_t ptLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(2)->GetXmax(); //all plots have same axes... + Double_t displLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(3)->GetXmax(); + Double_t EtaLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(4)->GetXmax(); + if(ptTrack > ptLim_Sparse) ptTrack = ptLim_Sparse-0.01; + if(d0Track > displLim_Sparse) d0Track = (displLim_Sparse-0.001); + if(deltaeta > EtaLim_Sparse) deltaeta = EtaLim_Sparse-0.01; + if(deltaeta < -EtaLim_Sparse) deltaeta = -EtaLim_Sparse+0.01; + + if(fMixing == kSE) { + + //variables for filling histos + Double_t fillSpPhiD0[5] = {deltaphi,mD0,ptTrack,d0Track,deltaeta}; + Double_t fillSpPhiD0bar[5] = {deltaphi,mD0bar,ptTrack,d0Track,deltaeta}; - if(!TrackSelectedInLoop(track,d,aod,ptbin,idDaughs)) continue; //rejection of track (daughter of D0/quality cuts not passed/soft pion/negative ID - - Int_t label = track->GetLabel(); - if (label<0) continue; //discard negative label tracks - Int_t PDGtrack = ((AliAODMCParticle*)mcArray->At(label))->GetPdgCode(); - - GetImpParameter(track,aod,d0,d0err); //evaluates d0 and sigma_d0 - - //Infos on track (origin, phi, eta) - Int_t origTr = CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(label)); - deltaphi = d->Phi()-track->Phi(); - deltaeta = d->Eta()-track->Eta(); - if (deltaphi < -1.571) deltaphi+=6.283; - if (deltaphi > 4.712) deltaphi-=6.283; - Double_t pttrack = track->Pt(); - - if (pttrack > highPt) {highPt = pttrack; lead[0] = pttrack; lead[1] = deltaphi; lead[2] = origTr;} //for leading particle - - //Lines needed to include overflow into THnSparse projections! - Double_t d0orig = d0; - Double_t ptLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(2)->GetXmax(); //all plots have same axes... - Double_t displLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(3)->GetXmax(); - Double_t EtaLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(4)->GetXmax(); - if(pttrack > ptLim_Sparse) pttrack = ptLim_Sparse-0.01; - if(d0/d0err > displLim_Sparse) d0 = (displLim_Sparse-0.001)*d0err; - if(deltaeta > EtaLim_Sparse) deltaeta = EtaLim_Sparse-0.01; - if(deltaeta < -EtaLim_Sparse) deltaeta = -EtaLim_Sparse+0.01; - - //variables for filling histos - Double_t fillSpPhiD0[5] = {deltaphi,mD0,pttrack,d0/d0err,deltaeta}; - Double_t fillSpPhiD0bar[5] = {deltaphi,mD0bar,pttrack,d0/d0err,deltaeta}; - - //generic charged tracks case - if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==421) { //D0 (from MCTruth) - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->Fill(fillSpPhiD0); - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0); - if(origD0==4&&origTr>=1&&origTr<=2) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0);//c tr - if(origD0==5&&origTr>=3&&origTr<=6) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0);//b tr + if(fReadMC == 0) { + //sparse fill for data (tracks, K+-, K0) + weighted + if(fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3) { //D0 + ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0); + ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0,pTorig); } - if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==-421) { //D0bar - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->Fill(fillSpPhiD0bar); - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0bar); - if(origD0==4&&origTr>=1&&origTr<=2) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0bar); - if(origD0==5&&origTr>=3&&origTr<=6) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0bar); + if(fIsSelectedCandidate == 2 || fIsSelectedCandidate == 3) { //D0bar + ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar); + ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0bar,pTorig); } if(!fAlreadyFilled) { - ((TH1F*)fOutputStudy->FindObject(Form("histDispl_Charg_Bin%d",ptbin)))->Fill(d0orig); //Fills displacement histos - if (origTr>=1&&origTr<=6) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_Charg_HF_Bin%d",ptbin)))->Fill(d0orig); - if (origTr>=1&&origTr<=6) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_Charg_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(d0orig); - ((TH1F*)fOutputStudy->FindObject(Form("histDispl_Charg%s_Bin%d",orig.Data(),ptbin)))->Fill(d0orig); //Fills displacement histos - ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_Charg_Bin%d",ptbin)))->Fill(track->Pt()); - ((TH1F*)fOutputStudy->FindObject(Form("histOrig_Charg_Bin%d",ptbin)))->Fill(origTr); + ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_%s_Bin%d",part.Data(),ptbin)))->Fill(pTorig); + ((TH1F*)fOutputStudy->FindObject(Form("hist_PhiDistr_%s",part.Data())))->Fill(phiTr); } - N_Charg++; - //pT_weighted track correlations fill - if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==421) { //D0 - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0,pttrack); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pttrack); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigSqr_Bin%d",ptbin)))->Fill(deltaphi,mD0,pow(pttrack,2.)); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigSqr%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pow(pttrack,2.)); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigDist_Bin%d",ptbin)))->Fill(deltaphi,mD0,PtWeig(ptbin,pttrack)); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigDist%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,PtWeig(ptbin,pttrack)); - if(origD0==4&&origTr>=1&&origTr<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pttrack);//c tr - if(origD0==5&&origTr>=3&&origTr<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pttrack);//b tr - if(origD0==4&&origTr>=1&&origTr<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigSqr_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pow(pttrack,2.));//c tr - if(origD0==5&&origTr>=3&&origTr<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigSqr_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pow(pttrack,2.));//b tr - if(origD0==4&&origTr>=1&&origTr<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigDist_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,PtWeig(ptbin,pttrack));//c tr - if(origD0==5&&origTr>=3&&origTr<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigDist_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,PtWeig(ptbin,pttrack));//b tr - } - if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==-421) { //D0bar - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0bar,pttrack); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pttrack); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigSqr_Bin%d",ptbin)))->Fill(deltaphi,mD0bar,pow(pttrack,2.)); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigSqr%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pow(pttrack,2.)); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigDist_Bin%d",ptbin)))->Fill(deltaphi,mD0bar,PtWeig(ptbin,pttrack)); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigDist%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,PtWeig(ptbin,pttrack)); - if(origD0==4&&origTr>=1&&origTr<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pttrack);//c tr - if(origD0==5&&origTr>=3&&origTr<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pttrack);//b tr - if(origD0==4&&origTr>=1&&origTr<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigSqr_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pow(pttrack,2.));//c tr - if(origD0==5&&origTr>=3&&origTr<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigSqr_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pow(pttrack,2.));//b tr - if(origD0==4&&origTr>=1&&origTr<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigDist_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,PtWeig(ptbin,pttrack));//c tr - if(origD0==5&&origTr>=3&&origTr<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_WeigDist_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,PtWeig(ptbin,pttrack));//b tr - } - - //Charged Kaon identification (K from MCTruth) - if(TMath::Abs(PDGtrack) == 321) { - if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==421) { //D0 (from MCTruth) - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Kcharg_Bin%d",ptbin)))->Fill(fillSpPhiD0); - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Kcharg%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0); - if(origD0==4&&origTr>=1&&origTr<=2) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Kcharg_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0); - if(origD0==5&&origTr>=3&&origTr<=6) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Kcharg_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0); - } - if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==-421) { //D0bar - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Kcharg_Bin%d",ptbin)))->Fill(fillSpPhiD0bar); - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Kcharg%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0bar); - if(origD0==4&&origTr>=1&&origTr<=2) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Kcharg_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0bar); - if(origD0==5&&origTr>=3&&origTr<=6) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Kcharg_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0bar); - } - if(!fAlreadyFilled) { - ((TH1F*)fOutputStudy->FindObject(Form("histDispl_Kcharg_Bin%d",ptbin)))->Fill(d0orig); //Fills displacement histos - if (origTr>=1&&origTr<=6) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_Kcharg_HF_Bin%d",ptbin)))->Fill(d0orig); - if (origTr>=1&&origTr<=6) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_Kcharg_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(d0orig); - ((TH1F*)fOutputStudy->FindObject(Form("histDispl_Kcharg%s_Bin%d",orig.Data(),ptbin)))->Fill(d0orig); //Fills displacement histos - ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_Kcharg_Bin%d",ptbin)))->Fill(track->Pt()); - ((TH1F*)fOutputStudy->FindObject(Form("histOrig_Kcharg_Bin%d",ptbin)))->Fill(origTr); - } - N_Kcharg++; - - } //end hadron case + } - } //end tracks loop + if(fReadMC) { - //Kaon identification - TClonesArray *v0array = (TClonesArray*)aod->GetList()->FindObject("v0s"); - Int_t idArrayV0[v0array->GetEntriesFast()][2]; //array for skipping K0 double-counting - for(int iV0=0; iV0GetEntriesFast();iV0++) { - for(int j=0; j<2; j++) {idArrayV0[iV0][j]=-2;} - AliAODv0 *v0 = (AliAODv0*)v0array->UncheckedAt(iV0); + if(origD0==4) {orig = "_From_c";} else {orig = "_From_b";} - if(v0->Pt() < fPtThreshLow.at(ptbin) || v0->Pt() > fPtThreshUp.at(ptbin)) continue; //discard tracks outside pt range for hadrons/K - if(!SelectV0(v0,vtx,1,idArrayV0)) continue; //option 1 = for correlations - Int_t pdgCodes[2] = {211,211}; - Int_t labV0 = v0->MatchToMC(311,mcArray,2,pdgCodes); - if(labV0<=0) continue; - - Double_t ptV0=v0->Pt(), deltaphiV0, deltaetaV0; - Double_t ptLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(2)->GetXmax(); //all plots have same axes... - Double_t EtaLim_Sparse = ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_Charg_Bin%d",ptbin)))->GetAxis(4)->GetXmax(); - deltaphiV0 = d->Phi()-v0->Phi(); - deltaetaV0 = d->Eta()-v0->Eta(); - if (deltaphiV0 < -1.571) deltaphiV0+=6.283; - if (deltaphiV0 > 4.712) deltaphiV0-=6.283; - if (ptV0 > ptLim_Sparse) ptV0 = ptLim_Sparse-0.01; - if(deltaetaV0 > EtaLim_Sparse) deltaetaV0 = EtaLim_Sparse-0.01; - if(deltaetaV0 < -EtaLim_Sparse) deltaetaV0 = -EtaLim_Sparse+0.01; - - Double_t fillSpPhiD0K0[5] = {deltaphiV0,mD0,ptV0,0.,deltaetaV0}; - Double_t fillSpPhiD0barK0[5] = {deltaphiV0,mD0bar,ptV0,0.,deltaetaV0}; - - Int_t origV0 = CheckTrackOrigin(mcArray,(AliAODMCParticle*)mcArray->At(labV0)); - - if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==421) { //D0 (from MCTruth) - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_K0_Bin%d",ptbin)))->Fill(fillSpPhiD0K0); - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_K0%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0K0); - if(origD0==4&&origV0>=1&&origV0<=2) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_K0_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0K0); - if(origD0==5&&origV0>=3&&origV0<=6) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_K0_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0K0); - } - if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==-421) { //D0bar - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_K0_Bin%d",ptbin)))->Fill(fillSpPhiD0barK0); - ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_K0%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0barK0); - if(origD0==4&&origV0>=1&&origV0<=2) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_K0_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0barK0); - if(origD0==5&&origV0>=3&&origV0<=6) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_K0_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(fillSpPhiD0barK0); + //sparse fill for data (tracks, K+-, K0) + weighted + if(PdgD0==421) { //D0 (from MCTruth) + ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0); + ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0); + if(origD0==4&&origTr>=1&&origTr<=2) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0); + if(origD0==5&&origTr>=3&&origTr<=6) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0); + ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0,pTorig); + ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pTorig); + if(origD0==4&&origTr>=1&&origTr<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pTorig); + if(origD0==5&&origTr>=3&&origTr<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0,pTorig); } + if(PdgD0==-421) { //D0bar + ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d",part.Data(),ptbin)))->Fill(fillSpPhiD0bar); + ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar); + if(origD0==4&&origTr>=1&&origTr<=2) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar); + if(origD0==5&&origTr>=3&&origTr<=6) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(fillSpPhiD0bar); + ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_Bin%d",ptbin)))->Fill(deltaphi,mD0bar,pTorig); + ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pTorig); + if(origD0==4&&origTr>=1&&origTr<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pTorig); + if(origD0==5&&origTr>=3&&origTr<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Weig_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(deltaphi,mD0bar,pTorig); + } if(!fAlreadyFilled) { - ((TH1F*)fOutputStudy->FindObject(Form("histDispl_K0_Bin%d",ptbin)))->Fill(0.); //Fills displacement histos - if (origV0>=1&&origV0<=6) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_K0_HF_Bin%d",ptbin)))->Fill(0.); - if (origV0>=1&&origV0<=6) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_K0_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(0.); - ((TH1F*)fOutputStudy->FindObject(Form("histDispl_K0%s_Bin%d",orig.Data(),ptbin)))->Fill(0.); //Fills displacement histos - ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_K0_Bin%d",ptbin)))->Fill(v0->Pt()); - ((TH1F*)fOutputStudy->FindObject(Form("histOrig_K0_Bin%d",ptbin)))->Fill(origV0); + ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_Bin%d",part.Data(),ptbin)))->Fill(d0orig); //Fills displacement histos + if (origTr>=1&&origTr<=6) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_HF_Bin%d",part.Data(),ptbin)))->Fill(d0orig); + if (origTr>=1&&origTr<=6) ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s_HF%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(d0orig); + ((TH1F*)fOutputStudy->FindObject(Form("histDispl_%s%s_Bin%d",part.Data(),orig.Data(),ptbin)))->Fill(d0orig); //Fills displacement histos + ((TH1F*)fOutputStudy->FindObject(Form("hist_Pt_%s_Bin%d",part.Data(),ptbin)))->Fill(pTorig); + ((TH1F*)fOutputStudy->FindObject(Form("histOrig_%s_Bin%d",part.Data(),ptbin)))->Fill(origTr); + ((TH1F*)fOutputStudy->FindObject(Form("hist_PhiDistr_%s",part.Data())))->Fill(phiTr); } - N_Kaons++; - } // end kaon case + }//end MC case - //leading track correlations fill - if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==421) { //D0 - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[1],mD0); //c and b D0 - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[1],mD0); //c or b D0 - if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[1],mD0); - if(origD0==5&&(int)lead[2]>=3&&(int)lead[2]<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[1],mD0); - } - if(((AliAODMCParticle*)mcArray->At(labD0))->GetPdgCode()==-421) { //D0bar - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_Bin%d",ptbin)))->Fill(lead[1],mD0bar); - ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[1],mD0bar); //c or b D0 - if(origD0==4&&(int)lead[2]>=1&&(int)lead[2]<=2) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[1],mD0bar); - if(origD0==5&&(int)lead[2]>=3&&(int)lead[2]<=6) ((TH2F*)fOutputCorr->FindObject(Form("hPhi_Lead_HF%s_Bin%d",orig.Data(),ptbin)))->Fill(lead[1],mD0bar); - } + } //end of SE fill - } //end fReadMC = 1 + if(fMixing == kME) { - if (!fAlreadyFilled) { - ((TH1F*)fOutputStudy->FindObject(Form("hist_Count_Charg_Bin%d",ptbin)))->Fill(N_Charg); - ((TH1F*)fOutputStudy->FindObject(Form("hist_Count_Kcharg_Bin%d",ptbin)))->Fill(N_Kcharg); - ((TH1F*)fOutputStudy->FindObject(Form("hist_Count_K0_Bin%d",ptbin)))->Fill(N_Kaons); - } + //variables for filling histos + Double_t fillSpPhiD0[3] = {deltaphi,mD0,deltaeta}; + Double_t fillSpPhiD0bar[3] = {deltaphi,mD0bar,deltaeta}; - fAlreadyFilled=kTRUE; //distribution plots for tracks filled + if(fReadMC == 0) { + //sparse fill for data (tracks, K+-, K0) + if(fIsSelectedCandidate == 1||fIsSelectedCandidate == 3) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0); + if(fIsSelectedCandidate == 2||fIsSelectedCandidate == 3) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0bar); + } + if(fReadMC == 1) { + //sparse fill for data (tracks, K+-, K0) + if(PdgD0==421) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0); + if(PdgD0==-421) ((THnSparseI*)fOutputCorr->FindObject(Form("hPhi_%s_Bin%d_EvMix",part.Data(),ptbin)))->Fill(fillSpPhiD0bar); + }//end MC case + } //end of ME fill + + return; } //_________________________________________________________________________________________________ @@ -1533,7 +1579,7 @@ Int_t AliAnalysisTaskSED0Correlations::CheckTrackOrigin(TClonesArray* arrayMC, A // 7) c hadronization // 8) b hadronization // - printf(" AliAnalysisTaskSED0Correlations::CheckD0Origin() \n"); + if(fDebug>2) printf("AliAnalysisTaskSED0Correlations::CheckTrkOrigin() \n"); Int_t pdgGranma = 0; Int_t mother = 0; @@ -1591,7 +1637,6 @@ Int_t AliAnalysisTaskSED0Correlations::CheckTrackOrigin(TClonesArray* arrayMC, A else return 0; //no HF quark } - //________________________________________________________________________ Int_t AliAnalysisTaskSED0Correlations::PtBinCorr(Double_t pt) const { // @@ -1606,62 +1651,6 @@ Int_t AliAnalysisTaskSED0Correlations::PtBinCorr(Double_t pt) const { return ptbin; } -//________________________________________________________________________ -Double_t AliAnalysisTaskSED0Correlations::PtWeig(Int_t ptbin, Double_t x) const { - // - //gives the weight for Weighted Corrs (inverse of pT distribution, from 3 D0 pt bins) - // - if(x>11.5) x=11.5; if(x<0.3) x=0.3; //bounds for fit of distributions! - if(ptbin >= 3 && ptbin <= 5) return 1/(0.22958*(1.507e+03-5.035e+02*x+5.681e+01*x*x-2.186e+00*x*x*x+exp(1.336e+01-2.146e+00*x+1.206e-01*x*x))); - if(ptbin >= 6 && ptbin <= 8) return 1/(1.71828*(1.984e+02-6.113e+01*x+6.444e+00*x*x-2.342e-01*x*x*x+exp(1.023e+01-2.059e+00*x+1.194e-01*x*x))); - if(ptbin >= 9 && ptbin <= 10) return 1/(1.25905*(1.990e+02-6.306e+01*x+6.995e+00*x*x-2.681e-01*x*x*x+exp(9.125e+00-2.053e+00*x+1.276e-01*x*x))); - - return 0; //for other bins plot is disabled! -} - -//________________________________________________________________________ -void AliAnalysisTaskSED0Correlations::GetImpParameter(AliAODTrack *track, AliAODEvent* aod, Double_t &d0, Double_t &d0err) const { - // - //calculates d0 and error on d0 for the track - // - Double_t d0z0[2],covd0z0[3]; - AliESDtrack esdTrack(track); // AliESDTrack conversion of AOD track - esdTrack.PropagateToDCA(aod->GetPrimaryVertex(),aod->GetMagneticField(),10000.,d0z0,covd0z0); - d0 = TMath::Abs(d0z0[0]); // impact parameter - d0err = TMath::Sqrt(covd0z0[0]); // resolution of impact parameter -} - -//________________________________________________________________________ -Bool_t AliAnalysisTaskSED0Correlations::TrackSelectedInLoop(AliAODTrack* track, AliAODRecoDecayHF2Prong *d, AliAODEvent *aod, Int_t ptbin, Int_t idDaughs[]) const { - // - //rejection of tracks in loop for correlations - // - Bool_t output = kTRUE; - AliAODVertex *vtx = (AliAODVertex*)aod->GetPrimaryVertex(); - Double_t bz = aod->GetMagneticField(); - Double_t mD0, mD0bar; - d->InvMassD0(mD0,mD0bar); - - if(track->GetID() == idDaughs[0] || track->GetID() == idDaughs[1] || track->GetID() < 0) output = kFALSE; //discards daughters of candidate - if(track->Pt() < fPtThreshLow.at(ptbin) || track->Pt() > fPtThreshUp.at(ptbin)) output = kFALSE; //discard tracks outside pt range for hadrons/K - if(!fCutsTracks->IsHadronSelected(track,vtx,bz)) { //track discarded by quality cuts - ((TH1F*)fOutputStudy->FindObject("hRejTracks"))->Fill(1); - output = kFALSE; - } else ((TH1F*)fOutputStudy->FindObject("hRejTracks"))->Fill(0); //track accepted by quality cuts - if(!fCutsTracks->InvMassDstarRejection(d,track,fIsSelectedCandidate)) { - if (fReadMC == 0 && (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3)) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0); - if (fReadMC == 0 && fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,mD0bar); - if (fReadMC == 1) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,1.864); - output = kFALSE; - } else { - if (fReadMC == 0 && (fIsSelectedCandidate == 1 || fIsSelectedCandidate == 3)) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0); - if (fReadMC == 0 && fIsSelectedCandidate >= 2) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(0.,mD0bar); - if (fReadMC == 1) ((TH2F*)fOutputStudy->FindObject(Form("hDstarPions_Bin%d",ptbin)))->Fill(1.,1.864); - } - - return output; -} - //--------------------------------------------------------------------------- Bool_t AliAnalysisTaskSED0Correlations::SelectV0(AliAODv0* v0, AliAODVertex *vtx, Int_t opt, Int_t idArrayV0[][2]) const { @@ -1676,9 +1665,7 @@ Bool_t AliAnalysisTaskSED0Correlations::SelectV0(AliAODv0* v0, AliAODVertex *vtx AliAODTrack *v0Daug2 = (AliAODTrack*)v0->GetDaughter(1); if(opt==1) { //additional cuts for correlations (V0 has to be closer than 3 sigma from K0 mass) - if(v0->Pt() < 3. && TMath::Abs(v0->MassK0Short()-0.4976) > 3*0.0040) return kFALSE; - if(v0->Pt() > 3. && v0->Pt() < 6. && TMath::Abs(v0->MassK0Short()-0.4976) > 3*0.0052) return kFALSE; - if(v0->Pt() > 6. && TMath::Abs(v0->MassK0Short()-0.4976) > 3*0.0075) return kFALSE; + if(TMath::Abs(v0->MassK0Short()-0.4976) > 3*0.004) return kFALSE; } //This part removes double counting for swapped tracks! @@ -1715,5 +1702,7 @@ void AliAnalysisTaskSED0Correlations::PrintBinsAndLimits() { cout << "\n--------------------------\n"; cout << "MC Truth = "< #include #include #include @@ -25,10 +24,9 @@ #include "AliAnalysisTaskSE.h" #include "AliRDHFCutsD0toKpi.h" #include "AliHFAssociatedTrackCuts.h" +#include "AliHFCorrelator.h" #include "AliNormalizationCounter.h" -using std::vector; - class AliAODEvent; class AliAnalysisTaskSED0Correlations : public AliAnalysisTaskSE @@ -66,22 +64,23 @@ class AliAnalysisTaskSED0Correlations : public AliAnalysisTaskSE void SetPtTreshUp(Double_t* pttreshup) {for(int i=0;i fBinLimsCorr; // limits of pt bins per correlations @@ -96,7 +95,11 @@ class AliAnalysisTaskSED0Correlations : public AliAnalysisTaskSE TH1F *fNentries; //!histogram with number of events on output slot 2 AliRDHFCutsD0toKpi *fCutsD0; // Cuts for D0, output 3 AliHFAssociatedTrackCuts *fCutsTracks;// Cuts for tracks and K0, output 7 + AliHFCorrelator* fCorrelatorTr; // Correlator for tracks + AliHFCorrelator* fCorrelatorKc; // Correlator for charged K + AliHFCorrelator* fCorrelatorK0; // Correlator for K0 Bool_t fReadMC; // flag for MC array: kTRUE = read it, kFALSE = do not read it + Bool_t fMixing; // flag to enable also event mixing AliNormalizationCounter *fCounter; //!AliNormalizationCounter on output slot 4 Int_t fNPtBins; // Number of pt bins Int_t fFillOnlyD0D0bar; // flag to fill mass histogram with D0/D0bar only (0 = fill with both, 1 = fill with D0 only, 2 = fill with D0bar only) @@ -105,7 +108,7 @@ class AliAnalysisTaskSED0Correlations : public AliAnalysisTaskSE 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) - ClassDef(AliAnalysisTaskSED0Correlations,1); // AliAnalysisTaskSE for D0->Kpi + ClassDef(AliAnalysisTaskSED0Correlations,2); // AliAnalysisTaskSE for D0->Kpi }; #endif diff --git a/PWGHF/correlationHF/AliHFAssociatedTrackCuts.cxx b/PWGHF/correlationHF/AliHFAssociatedTrackCuts.cxx index accb795a23b..aef609c8902 100644 --- a/PWGHF/correlationHF/AliHFAssociatedTrackCuts.cxx +++ b/PWGHF/correlationHF/AliHFAssociatedTrackCuts.cxx @@ -41,6 +41,17 @@ AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(): AliAnalysisCuts(), fESDTrackCuts(0), fPidObj(0), + +fPoolMaxNEvents(0), +fPoolMinNTracks(0), +fMinEventsToMix(0), +fNzVtxBins(0), +fNzVtxBinsDim(0), +fZvtxBins(0), +fNCentBins(0), +fNCentBinsDim(0), +fCentBins(0), + fNTrackCuts(0), fAODTrackCuts(0), fTrackCutsNames(0), @@ -63,6 +74,17 @@ AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(const char* name, const char* AliAnalysisCuts(name,title), fESDTrackCuts(0), fPidObj(0), + +fPoolMaxNEvents(0), +fPoolMinNTracks(0), +fMinEventsToMix(0), +fNzVtxBins(0), +fNzVtxBinsDim(0), +fZvtxBins(0), +fNCentBins(0), +fNCentBinsDim(0), +fCentBins(0), + fNTrackCuts(0), fAODTrackCuts(0), fTrackCutsNames(0), @@ -81,6 +103,17 @@ AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(const AliHFAssociatedTrackCut AliAnalysisCuts(source), fESDTrackCuts(source.fESDTrackCuts), fPidObj(source.fPidObj), + +fPoolMaxNEvents(source.fPoolMaxNEvents), +fPoolMinNTracks(source.fPoolMinNTracks), +fMinEventsToMix(source.fMinEventsToMix), +fNzVtxBins(source.fNzVtxBins), +fNzVtxBinsDim(source.fNzVtxBinsDim), +fZvtxBins(source.fZvtxBins), +fNCentBins(source.fNCentBins), +fNCentBinsDim(source.fNCentBinsDim), +fCentBins(source.fCentBins), + fNTrackCuts(source.fNTrackCuts), fAODTrackCuts(source.fAODTrackCuts), fTrackCutsNames(source.fTrackCutsNames), @@ -127,6 +160,8 @@ AliHFAssociatedTrackCuts::~AliHFAssociatedTrackCuts() { if(fESDTrackCuts) {delete fESDTrackCuts; fESDTrackCuts = 0;} if(fPidObj) {delete fPidObj; fPidObj = 0;} + if(fZvtxBins) {delete[] fZvtxBins; fZvtxBins=0;} + if(fCentBins) {delete[] fCentBins; fCentBins=0;} if(fAODTrackCuts) {delete[] fAODTrackCuts; fAODTrackCuts=0;} if(fTrackCutsNames) {delete[] fTrackCutsNames; fTrackCutsNames=0;} if(fAODvZeroCuts){delete[] fAODvZeroCuts; fAODvZeroCuts=0;} @@ -140,29 +175,24 @@ Bool_t AliHFAssociatedTrackCuts::IsInAcceptance() return kFALSE; } //-------------------------------------------------------------------------- -Bool_t AliHFAssociatedTrackCuts::IsHadronSelected(AliAODTrack * track, AliAODVertex *vtx1, Double_t bz) +Bool_t AliHFAssociatedTrackCuts::IsHadronSelected(AliAODTrack * track) { - - - - Double_t d0 = -999.; - Double_t d0err = -998.; - Double_t VeryBig = 100; - AliESDtrack esdtrack(track); - Double_t d0z0[2],covd0z0[3]; + if(!fESDTrackCuts->IsSelected(&esdtrack)) return kFALSE; + return kTRUE; - esdtrack.PropagateToDCA(vtx1,bz,VeryBig,d0z0,covd0z0); +} + +//-------------------------------------------------------------------------- +Bool_t AliHFAssociatedTrackCuts::CheckHadronKinematic(Double_t pt, Double_t d0) +{ - d0 = TMath::Abs(d0z0[0]); - d0err = TMath::Sqrt(covd0z0[0]); - if(!fESDTrackCuts->IsSelected(&esdtrack)) return kFALSE; - if(track->Pt() < fAODTrackCuts[0]) return kFALSE; - if(track->Pt() > fAODTrackCuts[1]) return kFALSE; - if(d0 < fAODTrackCuts[2]) return kFALSE; - if(d0 > fAODTrackCuts[3]) return kFALSE; + if(pt < fAODTrackCuts[0]) return kFALSE; + if(pt > fAODTrackCuts[1]) return kFALSE; + if(d0 < fAODTrackCuts[2]) return kFALSE; + if(d0 > fAODTrackCuts[3]) return kFALSE; return kTRUE; @@ -273,17 +303,57 @@ Int_t AliHFAssociatedTrackCuts::IsMCpartFromHF(Int_t label, TClonesArray*mcArray if (isBeauty && isB) return 55; return 0; } + +//-------------------------------------------------------------------------- +Bool_t AliHFAssociatedTrackCuts::InvMassDstarRejection(AliAODRecoDecayHF2Prong* d, AliAODTrack *track, Int_t hypD0) const { + // + // Calculates invmass of track+D0 and rejects if compatible with D* + // (to remove pions from D*) + // + Double_t nsigma = 3.; + + Double_t mD0, mD0bar; + d->InvMassD0(mD0,mD0bar); + + Double_t invmassDstar1 = 0, invmassDstar2 = 0; + Double_t e1Pi = d->EProng(0,211), e2K = d->EProng(1,321); //hyp 1 (pi,K) - D0 + Double_t e1K = d->EProng(0,321), e2Pi = d->EProng(1,211); //hyp 2 (K,pi) - D0bar + Double_t psum2 = (d->Px()+track->Px())*(d->Px()+track->Px()) + +(d->Py()+track->Py())*(d->Py()+track->Py()) + +(d->Pz()+track->Pz())*(d->Pz()+track->Pz()); + + switch(hypD0) { + case 1: + invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2); + if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE; + break; + case 2: + invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2); + if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE; + break; + case 3: + invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2); + invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2); + if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE; + if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE; + break; + } + + return kTRUE; +} //________________________________________________________ void AliHFAssociatedTrackCuts::SetAODTrackCuts(Float_t *cutsarray) { - + cout << "Check 1" << endl; if(!fAODTrackCuts) fAODTrackCuts = new Float_t[fNTrackCuts]; - + cout << "Check 2" << endl; for(Int_t i =0; iInvMassD0(mD0,mD0bar); - - Double_t invmassDstar1 = 0, invmassDstar2 = 0; - Double_t e1Pi = d->EProng(0,211), e2K = d->EProng(1,321); //hyp 1 (pi,K) - D0 - Double_t e1K = d->EProng(0,321), e2Pi = d->EProng(1,211); //hyp 2 (K,pi) - D0bar - Double_t psum2 = (d->Px()+track->Px())*(d->Px()+track->Px()) - +(d->Py()+track->Py())*(d->Py()+track->Py()) - +(d->Pz()+track->Pz())*(d->Pz()+track->Pz()); - - switch(hypD0) { - case 1: - invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2); - if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE; - break; - case 2: - invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2); - if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE; - break; - case 3: - invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2); - invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2); - if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE; - if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE; - break; - } - - return kTRUE; +void AliHFAssociatedTrackCuts::PrintPoolParameters() +{ + printf("\nEvent Pool settings: \n \n"); + + printf("Number of zVtx Bins: %d\n", fNzVtxBins); + printf("\nzVtx Bins:\n"); + //Double_t zVtxbinLims[fNzVtxBins+1] = fNzVtxBins; + for(Int_t k=0; k +#include +#include +#include "TROOT.h" +#include "AliHFCorrelator.h" +#include "AliRDHFCutsDStartoKpipi.h" +#include "AliHFAssociatedTrackCuts.h" +#include "AliEventPoolManager.h" +#include "AliReducedParticle.h" +#include "AliCentrality.h" + +//_____________________________________________________ +AliHFCorrelator::AliHFCorrelator() : +// +// default constructor +// +TNamed(), +fPoolMgr(0x0), +fPool(0x0), +fhadcuts(0x0), +fAODEvent(0x0), +fAssociatedTracks(0x0), +fmcArray(0x0), +fReducedPart(0x0), +fD0cand(0x0), +fhypD0(0), + +fmixing(kFALSE), +fmontecarlo(kFALSE), +fsystem(kFALSE), +fselect(0), +fUseImpactParameter(0), +fPIDmode(0), + +fNofTracks(0), +fPoolContent(0), + +fPhiMin(0), +fPhiMax(0), + +fPtTrigger(0), +fPhiTrigger(0), +fEtaTrigger(0), + + +fDeltaPhi(0), +fDeltaEta(0), +fk0InvMass(0) + +{ + // default constructor +} + + + +//_____________________________________________________ +AliHFCorrelator::AliHFCorrelator(const Char_t* name, AliHFAssociatedTrackCuts *cuts, Bool_t ppOrPbPb) : +TNamed(name,"title"), +fPoolMgr(0x0), +fPool(0x0), +fhadcuts(0x0), +fAODEvent(0x0), +fAssociatedTracks(0x0), +fmcArray(0x0), +fReducedPart(0x0), +fD0cand(0x0), +fhypD0(0), + +fmixing(kFALSE), +fmontecarlo(kFALSE), +fsystem(ppOrPbPb), +fselect(0), +fUseImpactParameter(0), +fPIDmode(0), + +fNofTracks(0), +fPoolContent(0), + +fPhiMin(0), +fPhiMax(0), + +fPtTrigger(0), +fPhiTrigger(0), +fEtaTrigger(0), + + +fDeltaPhi(0), +fDeltaEta(0), +fk0InvMass(0) +{ + fhadcuts = cuts; +} + +//_____________________________________________________ +AliHFCorrelator::~AliHFCorrelator() +{ +// +// destructor +// + + if(fPoolMgr) {delete fPoolMgr; fPoolMgr=0;} + if(fPool) {delete fPool; fPool=0;} + if(fhadcuts) {delete fhadcuts; fhadcuts=0;} + if(fAODEvent) {delete fAODEvent; fAODEvent=0;} + if(fAssociatedTracks) {delete fAssociatedTracks; fAssociatedTracks=0;} + if(fmcArray) {delete fmcArray; fmcArray=0;} + if(fReducedPart) {delete fReducedPart; fReducedPart=0;} + if(fD0cand) {delete fD0cand; fD0cand=0;} + + + if(fNofTracks) fNofTracks = 0; + + if(fPhiMin) fPhiMin = 0; + if(fPhiMax) fPhiMax = 0; + + if(fPtTrigger) fPtTrigger=0; + if(fPhiTrigger) fPhiTrigger=0; + if(fEtaTrigger) fEtaTrigger=0; + + if(fDeltaPhi) fDeltaPhi=0; + if(fDeltaEta) fDeltaEta=0; + + if(fk0InvMass) fk0InvMass=0; + +} +//_____________________________________________________ +Bool_t AliHFCorrelator::DefineEventPool(){ + // definition of the Pool Manager for Event Mixing + + + Int_t MaxNofEvents = fhadcuts->GetMaxNEventsInPool(); + Int_t MinNofTracks = fhadcuts->GetMinNTracksInPool(); + Int_t NofCentBins = fhadcuts->GetNCentPoolBins(); + Double_t * CentBins = fhadcuts->GetCentPoolBins(); + Int_t NofZVrtxBins = fhadcuts->GetNZvtxPoolBins(); + Double_t *ZVrtxBins = fhadcuts->GetZvtxPoolBins(); + + + fPoolMgr = new AliEventPoolManager(MaxNofEvents, MinNofTracks, NofCentBins, CentBins, NofZVrtxBins, ZVrtxBins); + if(!fPoolMgr) return kFALSE; + return kTRUE; +} +//_____________________________________________________ +Bool_t AliHFCorrelator::Initialize(){ + + cout << "AliHFCorrelator::Initialize"<< endl; + if(!fAODEvent) cout << "No AOD event" << endl; + + AliCentrality *centralityObj = 0; + Int_t multiplicity = -1; + Double_t MultipOrCent = -1; + + // initialize the pool for event mixing + if(!fsystem){ // pp + multiplicity = fAODEvent->GetNTracks(); + MultipOrCent = multiplicity; // convert from Int_t to Double_t + } + if(fsystem){ // PbPb + + centralityObj = fAODEvent->GetHeader()->GetCentralityP(); + MultipOrCent = centralityObj->GetCentralityPercentileUnchecked("V0M"); + AliInfo(Form("Centrality is %f", MultipOrCent)); + } + + AliAODVertex *vtx = fAODEvent->GetPrimaryVertex(); + Double_t zvertex = vtx->GetZ(); // zvertex + Double_t * CentBins = fhadcuts->GetCentPoolBins(); + Double_t poolmin=CentBins[0]; + Double_t poolmax=CentBins[fhadcuts->GetNCentPoolBins()]; + cout << "pool limits in centr/multipl " << poolmin << " - " << poolmax << endl; +/* for (Int_t iM=0; iM evp; + for (Int_t iZ=0; iZ=10 || multip>500 || multip == 0) { +// AliInfo(Form("Event with Zvertex = %.2f cm and multiplicity = %.0f out of pool bounds, SKIPPING",zvertex,multip)); +// return kFALSE; +// } + + if(TMath::Abs(zvertex)>=10 || MultipOrCent>poolmax || MultipOrCent < poolmin) { + if(!fsystem)AliInfo(Form("pp Event with Zvertex = %.2f cm and multiplicity = %.0f out of pool bounds, SKIPPING",zvertex,MultipOrCent)); + else AliInfo(Form("PbPb Event with Zvertex = %.2f cm and centrality = %.1f out of pool bounds, SKIPPING",zvertex,MultipOrCent)); + + return kFALSE; + } + + fPool = fPoolMgr->GetEventPool(MultipOrCent, zvertex); + + if (!fPool){ + AliInfo(Form("No pool found for multiplicity = %f, zVtx = %f cm", MultipOrCent, zvertex)); + return kFALSE; + } + fPool->PrintInfo(); + return kTRUE; +} + +//_____________________________________________________ +Bool_t AliHFCorrelator::ProcessEventPool(){ + // analysis on Mixed Events + cout << "AliHFCorrelator::ProcessEventPool"<< endl; + if(!fmixing) return kFALSE; + if(!fPool->IsReady()) return kFALSE; + if(fPool->GetCurrentNEvents()GetMinEventsToMix()) return kFALSE; + fPool->PrintInfo(); + fPoolContent = fPool->GetCurrentNEvents(); + + return kTRUE; + +} + +//_____________________________________________________ +Bool_t AliHFCorrelator::ProcessAssociatedTracks(Int_t EventLoopIndex){ + + fAssociatedTracks = new TObjArray(); + + if(!fmixing){ // analysis on Single Event + + + if(fselect==1 || fselect ==2) fAssociatedTracks = AcceptAndReduceTracks(fAODEvent); + if(fselect==3) {fAssociatedTracks = AcceptAndReduceKZero(fAODEvent);} + + + } + + if(fmixing) { // analysis on Mixed Events + + + fAssociatedTracks = fPool->GetEvent(EventLoopIndex); + + + + + } // end if mixing + + if(!fAssociatedTracks) return kFALSE; + + fNofTracks = fAssociatedTracks->GetEntriesFast(); + + return kTRUE; + +} +//_____________________________________________________ +Bool_t AliHFCorrelator::Correlate(Int_t loopindex){ + + if(loopindex >= fNofTracks) return kFALSE; + if(!fAssociatedTracks) return kFALSE; + + fReducedPart = (AliReducedParticle*)fAssociatedTracks->At(loopindex); + + + fDeltaPhi = SetCorrectPhiRange(fPhiTrigger - fReducedPart->Phi()); + + fDeltaEta = fEtaTrigger - fReducedPart->Eta(); + + return kTRUE; + +} + +//_____________________________________________________ +Bool_t AliHFCorrelator::PoolUpdate(){ + + if(!fmixing) return kFALSE; + if(!fPool) return kFALSE; + if(fmixing) { // update the pool for Event Mixing + if(fselect==1 || fselect==2)fPool->UpdatePool(AcceptAndReduceTracks(fAODEvent)); // updating the pool for hadrons + if(fselect==3) fPool->UpdatePool(AcceptAndReduceKZero(fAODEvent)); // updating the pool for K0s + } + + return kTRUE; +} + +//_____________________________________________________ +Double_t AliHFCorrelator::SetCorrectPhiRange(Double_t phi){ + Double_t pi = TMath::Pi(); + + if(phifPhiMax) phi = phi - 2*pi; + + return phi; +} + +//_____________________________________________________ +TObjArray* AliHFCorrelator::AcceptAndReduceTracks(AliAODEvent* inputEvent){ + + Int_t nTracks = inputEvent->GetNTracks(); + AliAODVertex * vtx = inputEvent->GetPrimaryVertex(); + Double_t Bz = inputEvent->GetMagneticField(); + + + TObjArray* tracksClone = new TObjArray; + tracksClone->SetOwner(kTRUE); + for (Int_t iTrack=0; iTrackGetTrack(iTrack); + if (!track) continue; + if(!fhadcuts->IsHadronSelected(track)) continue; // apply ESD level selections + + Double_t pT = track->Pt(); + + //compute impact parameter + Double_t d0z0[2],covd0z0[3]; + Double_t d0=-999999.; + if(fUseImpactParameter) track->PropagateToDCA(vtx,Bz,100,d0z0,covd0z0); + else d0z0[0] = 1. ; // random number - be careful with the cuts you applied + + if(fUseImpactParameter==1) d0 = TMath::Abs(d0z0[0]); // use impact parameter + if(fUseImpactParameter==2) { // use impact parameter over resolution + if(TMath::Abs(covd0z0[0])>0.00000001) d0 = TMath::Abs(d0z0[0])/TMath::Abs(covd0z0[0]); + else d0 = -1.; // if the resoultion is Zero, rejects the track - to be on the safe side + + } + + + + + if(!fhadcuts->CheckHadronKinematic(pT,d0))continue; // apply kinematic cuts + Bool_t rejectsoftpi = kFALSE; + if(fD0cand) rejectsoftpi = fhadcuts->InvMassDstarRejection(fD0cand,track,fhypD0); + + + if(fselect ==2){ + if(!fhadcuts->CheckKaonCompatibility(track,fmontecarlo,fmcArray,fPIDmode)) continue; // check if it is a Kaon - data and MC + } + + //AliVParticle * part = (AliVParticle*)track; + tracksClone->Add(new AliReducedParticle(track->Eta(), track->Phi(), pT,track->GetLabel(),track->GetID(),d0,rejectsoftpi)); + } + + + return tracksClone; +} + +//_____________________________________________________ +TObjArray* AliHFCorrelator::AcceptAndReduceKZero(AliAODEvent* inputEvent){ + + Int_t nOfVZeros = inputEvent->GetNumberOfV0s(); + TObjArray* KZeroClone = new TObjArray; + AliAODVertex *vertex1 = (AliAODVertex*)inputEvent->GetPrimaryVertex(); + Int_t v0label = -1; + Int_t pdgDgK0toPipi[2] = {211,211}; + Double_t mPDGK0=0.497614;//TDatabasePDG::Instance()->GetParticle(310)->Mass(); + const Int_t size = inputEvent->GetNumberOfV0s(); + Int_t prevnegID[size]; + Int_t prevposID[size]; + for(Int_t iv0 =0; iv0< nOfVZeros; iv0++){// loop on all v0 candidates + AliAODv0 *v0 = (static_cast (inputEvent))->GetV0(iv0); + if(!v0) {cout << "is not a v0 at step " << iv0 << endl; continue;} + + if(!fhadcuts->IsKZeroSelected(v0,vertex1)) continue; // check if it is a k0 + + // checks to avoid double counting + Int_t negID = -999; + Int_t posID = -998; + //Int_t a = 0; + prevnegID[iv0] = -997; + prevposID[iv0]= -996; + negID = v0->GetNegID(); + posID = v0->GetPosID(); + Bool_t DoubleCounting = kFALSE; + + for(Int_t k=0; kMatchToMC(310,fmcArray, 0, pdgDgK0toPipi); //match a K0 short + Double_t k0pt = v0->Pt(); + Double_t k0eta = v0->Eta(); + Double_t k0Phi = v0->Phi(); + fk0InvMass = v0->MassK0Short(); + + //if (loopindex == 0) { + // if(!plotassociation) ((TH2F*)fOutput->FindObject("KZeroSpectra"))->Fill(k0InvMass,k0pt); // spectra for all k0 + // if(plotassociation) ((TH2F*)fOutput->FindObject("KZeroSpectraifHF"))->Fill(k0InvMass,k0pt); // spectra for k0 in association with a D* + //} + // if there are more D* candidates per event, loopindex == 0 makes sure you fill the mass spectra only once! + + if(TMath::Abs(fk0InvMass-mPDGK0)>3*0.004) continue; // select candidates within 3 sigma + KZeroClone->Add(new AliReducedParticle(k0eta,k0Phi,k0pt,v0label)); + + } + + return KZeroClone; +} diff --git a/PWGHF/correlationHF/AliHFCorrelator.h b/PWGHF/correlationHF/AliHFCorrelator.h new file mode 100644 index 00000000000..0744edbfce9 --- /dev/null +++ b/PWGHF/correlationHF/AliHFCorrelator.h @@ -0,0 +1,160 @@ +#ifndef AliHFCorrelator_H +#define AliHFCorrelator_H + +/************************************************************************** + * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. * + * * + * Author: The ALICE Off-line Project. * + * Contributors are mentioned in the code where appropriate. * + * * + * Permission to use, copy, modify and distribute this software and its * + * documentation strictly for non-commercial purposes is hereby granted * + * without fee, provided that the above copyright notice appears in all * + * copies and that both the copyright notice and this permission notice * + * appear in the supporting documentation. The authors make no claims * + * about the suitability of this software for any purpose. It is * + * provided "as is" without express or implied warranty. * + **************************************************************************/ + +// +// Base class for Heavy Flavour Correlations Analysis +// Single Event and Mixed Event Analysis are implemented +//----------------------------------------------------------------------- +// +// +// Author S.Bjelogrlic +// Utrecht University +// sandro.bjelogrlic@cern.ch +// +//----------------------------------------------------------------------- + +/* $Id$ */ + +#include "AliHFAssociatedTrackCuts.h" +#include "AliEventPoolManager.h" +#include "AliVParticle.h" +#include "AliReducedParticle.h" + + +class AliHFCorrelator : public TNamed +{ + + public: + + AliHFCorrelator(); + AliHFCorrelator(const Char_t* name, AliHFAssociatedTrackCuts *cuts, Bool_t ppOrPbPb); + virtual ~AliHFCorrelator(); + + //setters + void SetDeltaPhiInterval (Double_t min, Double_t max){ + fPhiMin = min; fPhiMax = max; + if(TMath::Abs(fPhiMin-fPhiMax) != 2*TMath::Pi()) AliInfo("AliHFCorrelator::Warning: the delta phi interval is not set to 2 Pi"); + } + void SetEventMixing(Bool_t mixON){fmixing=mixON;} + void SetTriggerParticleProperties(Double_t ptTrig, Double_t phiTrig, Double_t etaTrig) + {fPtTrigger = ptTrig; fPhiTrigger = phiTrig; fEtaTrigger = etaTrig;} + + + void SetAssociatedParticleType(Int_t type){fselect = type;} + void SetAODEvent(AliAODEvent* inputevent){fAODEvent = inputevent;} + void SetMCArray(TClonesArray* mcArray){fmcArray = mcArray;} + void SetUseMC(Bool_t useMC){fmontecarlo = useMC;} + void SetApplyDisplacementCut(Int_t applycut){fUseImpactParameter = applycut;} + void SetPIDmode(Int_t mode){fPIDmode = mode;} + + void SetD0Properties(AliAODRecoDecayHF2Prong* d, Int_t D0hyp) + {fD0cand = d; fhypD0 = D0hyp;} + + + + Bool_t DefineEventPool(); // Definition of the Event pool parameters + Bool_t Initialize(); // function that initlize everything for the analysis + Bool_t ProcessEventPool(); // processes the event pool + Bool_t ProcessAssociatedTracks(Int_t EventLoopIndex); // + Bool_t Correlate(Int_t loopindex); // function that computes the correlations between the trigger particle and the track n. loopindex + Bool_t PoolUpdate();// updates the event pool + Double_t SetCorrectPhiRange(Double_t phi); // sets all the angles in the correct range + + //getters + AliEventPool* GetPool() {return fPool;} + TObjArray * GetTrackArray(){return fAssociatedTracks;} + AliHFAssociatedTrackCuts* GetSelectionCuts() {return fhadcuts;} + AliReducedParticle* GetAssociatedParticle() {return fReducedPart;} + + Int_t GetNofTracks(){return fNofTracks;} + Int_t GetNofEventsInPool(){return fPoolContent;} + + /*Double_t GetAssociatedTrackPt() {return fPtAssoc;} // pt of the associated track + Double_t GetAssociatedTrackPhi() {return fPhiAssoc;} // phi of the associated track + Double_t GetAssociatedTrackEta() {return fEtaAssoc;} // Eta of the associated track + Int_t GetAssociatedTrackLabel() {return fTrackLabel;} // tracklabel associated track + Int_t GetAssociatedTrackID() {return fTrackID;} // trackID associated track + Double_t GetAssociatedTrackImpPar() {return fTrackImpPar;} // impact parameter of the track + Bool_t GetAssociatedTrackSoftPiCompatibility() {return fIsSoftPiCand;} + */ + Double_t GetDeltaPhi(){return fDeltaPhi;} // Delta Phi, needs to be called after the method correlate + Double_t GetDeltaEta(){return fDeltaEta;} // Delta Eta + + Double_t GetAssociatedKZeroInvariantmass(){return fk0InvMass;} + + + + // methods to reduce the tracks to correlate with track selection cuts applied here + TObjArray* AcceptAndReduceTracks(AliAODEvent* inputEvent); // selecting hadrons and kaons + TObjArray* AcceptAndReduceKZero(AliAODEvent* inputEvent); // selecting kzeros + + + private: + + AliHFCorrelator(const AliHFCorrelator& vtxr); + AliHFCorrelator& operator=(const AliHFCorrelator& vtxr ); + + AliEventPoolManager* fPoolMgr; //! event pool manager + AliEventPool * fPool; //! Pool for event mixing + AliHFAssociatedTrackCuts* fhadcuts;//! hadron cuts + AliAODEvent * fAODEvent;//! AOD Event + TObjArray* fAssociatedTracks; // Array of associated tracks + TClonesArray* fmcArray; //mcarray + AliReducedParticle * fReducedPart; // reduced AOD particle; + AliAODRecoDecayHF2Prong* fD0cand; //D0 candidate + Int_t fhypD0; //hypothesis necessary for + + Bool_t fmixing;// switch for event mixing + Bool_t fmontecarlo; // switch for MonteCarlo + Bool_t fsystem; // select pp (kFALSE) or PbPb (kTRUE) + Int_t fselect; // 1 for hadrons, 2 for kaons, 3 for KZeros + Int_t fUseImpactParameter; // switch to use the impact parameter cut + Int_t fPIDmode; // set the PID mode for Kaon identification + + Int_t fNofTracks; // number of tracks in track array + Int_t fPoolContent; // n of events in pool + + Double_t fPhiMin; // min for phi + Double_t fPhiMax; // max for phi + + Double_t fPtTrigger; // pt of the trigger D meson + Double_t fPhiTrigger; // phi of the trigger D meson + Double_t fEtaTrigger; // Eta of the trigger D meson + +// Double_t fPtAssoc; // pt of the trigger D meson +// Double_t fPhiAssoc; // phi of the trigger D meson +// Double_t fEtaAssoc; // Eta of the trigger D meson +// Int_t fTrackLabel; // tracklabel +// Int_t fTrackID; // trackID +// Double_t fTrackImpPar; // Impact parameter of the track in respect to the vertex +// Bool_t fIsSoftPiCand; // is kTRUE if the track is compatible with a soft pion hypothesis (necessary for D0 analysis) + + Double_t fDeltaPhi; // delta phi between D meson and associated track + Double_t fDeltaEta; // delta eta between D meson and associated track + + Double_t fk0InvMass; // KZero invariant mass + + + ClassDef(AliHFCorrelator,1); // class for HF correlations +}; + + + + + +#endif diff --git a/PWGHF/correlationHF/AliReducedParticle.cxx b/PWGHF/correlationHF/AliReducedParticle.cxx index f420e4b53f4..3df0170d4df 100644 --- a/PWGHF/correlationHF/AliReducedParticle.cxx +++ b/PWGHF/correlationHF/AliReducedParticle.cxx @@ -30,12 +30,28 @@ #include "AliVParticle.h" #include "AliReducedParticle.h" -AliReducedParticle::AliReducedParticle(Double_t eta, Double_t phi, Double_t pt, Int_t McLabel, Int_t trackid) : +AliReducedParticle::AliReducedParticle(Double_t eta, Double_t phi, Double_t pt, Int_t mcLabel, Int_t trackid, Double_t impPar, Bool_t checkSoftPi) : fEta(eta), fPhi(phi), fpT(pt), -fMcLabel(McLabel), -fid(trackid) +fMcLabel(mcLabel), +fid(trackid), +fImpPar(impPar), +fCheckSoftPi(checkSoftPi) +{ + // + // default constructor + // +} + +AliReducedParticle::AliReducedParticle(Double_t eta, Double_t phi, Double_t pt, Int_t McLabel) : +fEta(eta), +fPhi(phi), +fpT(pt), +fMcLabel(McLabel), +fid(0), +fImpPar(0.), +fCheckSoftPi(kFALSE) { // // default constructor @@ -48,5 +64,5 @@ AliReducedParticle::~AliReducedParticle() { // destructor // - Info("AliReducedParticle","Calling Destructor"); + //Info("AliReducedParticle","Calling Destructor"); } diff --git a/PWGHF/correlationHF/AliReducedParticle.h b/PWGHF/correlationHF/AliReducedParticle.h index 9c525dbd247..7edf23a6c29 100644 --- a/PWGHF/correlationHF/AliReducedParticle.h +++ b/PWGHF/correlationHF/AliReducedParticle.h @@ -38,20 +38,22 @@ // class to get the reduced hadron candidate class AliReducedParticle : public AliVParticle -//class AliReducedParticle : public AliAODTrack { public: - AliReducedParticle(Double_t eta, Double_t phi, Double_t pt, Int_t McLabel, Int_t trackid); + AliReducedParticle(Double_t eta, Double_t phi, Double_t pt, Int_t mcLabel, Int_t trackid, Double_t impPar, Bool_t checkSoftPi); + AliReducedParticle(Double_t eta, Double_t phi, Double_t pt, Int_t mcLabel); ~AliReducedParticle(); // kinematics virtual Double_t Pt() const { return fpT; } - virtual Double_t Phi() const { return fPhi; } - virtual Double_t Eta() const { return fEta; } - virtual Int_t GetLabel() const { return fMcLabel; } - virtual Int_t GetID() const{return fid;} + virtual Double_t Phi() const { return fPhi; } + virtual Double_t Eta() const { return fEta; } + virtual Int_t GetLabel() const { return fMcLabel; } + virtual Int_t GetID() const{return fid;} + virtual Double_t GetImpPar() const{return fImpPar;} + virtual Bool_t CheckSoftPi() const{return fCheckSoftPi;} // kinematics virtual Double_t Px() const { AliFatal("Not implemented"); return 0; } @@ -89,10 +91,11 @@ private: Double_t fEta; // eta Double_t fPhi; // phi Double_t fpT; // pT - Int_t fMcLabel; //mclabel - Int_t fid; - + Int_t fMcLabel; //mclabel + Int_t fid; // track ID + Double_t fImpPar; // impact parameter + Bool_t fCheckSoftPi; // check if the track is compatible with a softpion from D* - ClassDef(AliReducedParticle, 1); // class which contains only quantities requires for this analysis to reduce memory consumption for event mixing + ClassDef(AliReducedParticle, 2); // class which contains only quantities requires for this analysis to reduce memory consumption for event mixing }; #endif diff --git a/PWGHF/correlationHF/macros/AddTaskD0Correlations.C b/PWGHF/correlationHF/macros/AddTaskD0Correlations.C index 39e68c91228..cf55c31ec3d 100644 --- a/PWGHF/correlationHF/macros/AddTaskD0Correlations.C +++ b/PWGHF/correlationHF/macros/AddTaskD0Correlations.C @@ -1,4 +1,4 @@ -AliAnalysisTaskSED0Correlations *AddTaskD0Correlations(Bool_t readMC=kFALSE, Int_t system=0/*0=pp,1=PbPb*/, Int_t flagD0D0bar=0, +AliAnalysisTaskSED0Correlations *AddTaskD0Correlations(Bool_t readMC=kFALSE, Bool_t mixing=kFALSE, 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.root", TString cutsD0name="D0toKpiCuts", TString cutsTrkname="AssociatedTrkCuts", Bool_t flagAOD049=kFALSE, Int_t standardbins=1, Bool_t stdcuts=kFALSE) { @@ -93,8 +93,8 @@ AliAnalysisTaskSED0Correlations *AddTaskD0Correlations(Bool_t readMC=kFALSE, Int //Cuts for D0 AliRDHFCutsD0toKpi* RDHFD0Corrs=new AliRDHFCutsD0toKpi(); //RDHFD0Corrs->SetUsePhysicsSelection(kFALSE); + RDHFD0Corrs->SetLowPt(kFALSE); //low-pt special PID disabled if(stdcuts) { - RDHFD0Corrs->SetLowPt(kFALSE); //low-pt special PID disabled if(system==0) RDHFD0Corrs->SetStandardCutsPP2010(); else { RDHFD0Corrs->SetStandardCutsPbPb2010(); @@ -145,6 +145,7 @@ AliAnalysisTaskSED0Correlations *AddTaskD0Correlations(Bool_t readMC=kFALSE, Int AliAnalysisTaskSED0Correlations *massD0Task = new AliAnalysisTaskSED0Correlations(taskname.Data(),RDHFD0Corrs,corrCuts); massD0Task->SetDebugLevel(2); massD0Task->SetReadMC(readMC); + massD0Task->SetEvMixing(mixing); massD0Task->SetFillOnlyD0D0bar(flagD0D0bar); massD0Task->SetSystem(system); //0=pp, 1=PbPb diff --git a/PWGHF/correlationHF/macros/AddTaskDStarCorrelations.C b/PWGHF/correlationHF/macros/AddTaskDStarCorrelations.C index b657814e98f..f9603b10900 100644 --- a/PWGHF/correlationHF/macros/AddTaskDStarCorrelations.C +++ b/PWGHF/correlationHF/macros/AddTaskDStarCorrelations.C @@ -1,119 +1,121 @@ - //DEFINITION OF A FEW CONSTANTS -//---------------------------------------------------- - -/* $Id$ */ - -AliAnalysisTaskDStarCorrelations *AddTaskDStarCorrelations(Int_t theMCon =5, Int_t mixing=4, Int_t trackselect =1) -{ - - AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); - if (!mgr) { - ::Error("AliAnalysisTaskDStarCorrelations", "No analysis manager to connect to."); - return NULL; - } - - TFile* filecuts=new TFile("DStartoKpipiCuts_corr.root"); - if(!filecuts->IsOpen()){ - cout<<"Input file not found: exit"<IsOpen()){ - cout<<"Input file2 not found: exit"<Get("DStartoKpipiCuts"); - RDHFDStartoKpipi->SetName("DStartoKpipiCuts"); - - - AliHFAssociatedTrackCuts* corrCuts=new AliHFAssociatedTrackCuts(); - corrCuts = (AliHFAssociatedTrackCuts*)filecuts2->Get("AssociatedCuts"); - corrCuts->SetName("AssociatedCuts"); - corrCuts->PrintAll(); - - // mm let's see if everything is ok - if(!RDHFDStartoKpipi){ - cout<<"Specific AliRDHFCuts not found"<SetMonteCarlo(theMCon); - task->SetUseMixing(mixing); - task->SetCorrelator(trackselect) ; - //task->SetDebugLevel(0); - - - if(trackselect == 1) Info("AliAnalysisTaskDStarCorrelations","Correlating D* with charged hadrons \n"); - else if(trackselect == 2) Info("AliAnalysisTaskDStarCorrelations","Correlating D* with charged kaons \n"); - else if(trackselect == 3) Info("AliAnalysisTaskDStarCorrelations","Correlating D* with reconstructed K0s \n"); - else Fatal("AliAnalysisTaskDStarCorrelations","Nothing to correlate with!"); - if(mixing) Info ("AliAnalysisTaskDStarCorrelations","Event Mixing Analysis\n"); - if(!mixing) Info ("AliAnalysisTaskDStarCorrelations","Single Event Analysis \n"); - - // Create and connect containers for input/output - //TString dcavalue = " "; - if(!theMCon) TString contname = "Data"; - if(theMCon) TString contname = "MonteCarlo"; - if(trackselect ==1) TString particle = "Hadron"; - if(trackselect ==2) TString particle = "Kaon"; - if(trackselect ==3) TString particle = "KZero"; - - TString cutname = "cuts" ; - TString cutname2 = "hadroncuts" ; - TString outputfile = AliAnalysisManager::GetCommonFileName(); - TString counter = "NormCounter"; - outputfile += ":PWGHF_D2H_"; - if(!mixing) { - outputfile += "SE"; - contname += "SE"; - cutname += "SE"; - cutname2 += "SE"; - counter+= "SE"; - } - if(mixing){ - outputfile += "ME"; - contname += "ME"; - cutname += "ME"; - cutname2 += "ME"; - counter+= "ME"; - } - outputfile += "Dphi_DStar"; - outputfile += particle; - cutname += particle; - cutname2 += particle; - contname += particle; - counter+= particle; - - - //cout << "Contname = " << contname << endl; - mgr->AddTask(task); - // ------ input data ------ - AliAnalysisDataContainer *cinput0 = mgr->GetCommonInputContainer(); - - // ----- output data ----- - - // output TH1I for event counting - AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(contname, TList::Class(),AliAnalysisManager::kOutputContainer,outputfile.Data()); - AliAnalysisDataContainer *coutput2 = mgr->CreateContainer(cutname,AliRDHFCutsDStartoKpipi::Class(),AliAnalysisManager::kOutputContainer, outputfile.Data()); //cuts - AliAnalysisDataContainer *coutputDstarNorm = mgr->CreateContainer(counter,AliNormalizationCounter::Class(),AliAnalysisManager::kOutputContainer, outputfile.Data()); - AliAnalysisDataContainer *coutput4 = mgr->CreateContainer(cutname2,AliHFAssociatedTrackCuts::Class(),AliAnalysisManager::kOutputContainer, outputfile.Data()); //cuts - - mgr->ConnectInput(task,0,mgr->GetCommonInputContainer()); - mgr->ConnectOutput(task,1,coutput1); - mgr->ConnectOutput(task,2,coutput2); - mgr->ConnectOutput(task,3,coutputDstarNorm); - mgr->ConnectOutput(task,4,coutput4); - - return task ; - -} - + //DEFINITION OF A FEW CONSTANTS +//---------------------------------------------------- + +/* $Id$ */ + +AliAnalysisTaskDStarCorrelations *AddTaskDStarCorrelations(Bool_t runOnPbPb,Bool_t theMCon, Bool_t mixing, Int_t trackselect =1, Int_t usedispl =0, TString DCutObjName = "DStartoKpipiCuts_corr.root", TString TrackCutObjName = "AssocPartCuts.root") +{ + + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + if (!mgr) { + ::Error(" AliAnalysisTaskDStarCorrelations", "No analysis manager to connect to."); + return NULL; + } + + TFile* filecuts=new TFile(DCutObjName.Data()); + if(!filecuts->IsOpen()){ + cout<<"Input file not found: exit"<IsOpen()){ + cout<<"Input file2 not found: exit"<Get("DStartoKpipiCuts"); + RDHFDStartoKpipi->SetName("DStartoKpipiCuts"); + + + AliHFAssociatedTrackCuts* corrCuts=new AliHFAssociatedTrackCuts(); + corrCuts = (AliHFAssociatedTrackCuts*)filecuts2->Get("AssociatedCuts"); + corrCuts->SetName("AssociatedCuts"); + corrCuts->PrintAll(); + + // mm let's see if everything is ok + if(!RDHFDStartoKpipi){ + cout<<"Specific AliRDHFCuts not found"<SetMonteCarlo(theMCon); + task->SetUseMixing(mixing); + task->SetCorrelator(trackselect) ; + task->SetUseDisplacement(usedispl); + task->SetRunPbPb(runOnPbPb); + //task->SetDebugLevel(0); + + + if(trackselect == 1) Info(" AliAnalysisTaskDStarCorrelations","Correlating D* with charged hadrons \n"); + else if(trackselect == 2) Info(" AliAnalysisTaskDStarCorrelations","Correlating D* with charged kaons \n"); + else if(trackselect == 3) Info(" AliAnalysisTaskDStarCorrelations","Correlating D* with reconstructed K0s \n"); + else Fatal(" AliAnalysisTaskDStarCorrelations","Nothing to correlate with!"); + if(mixing) Info (" AliAnalysisTaskDStarCorrelations","Event Mixing Analysis\n"); + if(!mixing) Info (" AliAnalysisTaskDStarCorrelations","Single Event Analysis \n"); + + // Create and connect containers for input/output + //TString dcavalue = " "; + if(!theMCon) TString contname = "Data"; + if(theMCon) TString contname = "MonteCarlo"; + if(trackselect ==1) TString particle = "Hadron"; + if(trackselect ==2) TString particle = "Kaon"; + if(trackselect ==3) TString particle = "KZero"; + + TString cutname = "cuts" ; + TString cutname2 = "hadroncuts" ; + TString outputfile = AliAnalysisManager::GetCommonFileName(); + TString counter = "NormCounter"; + outputfile += ":PWGHF_D2H_"; + if(!mixing) { + outputfile += "SE"; + contname += "SE"; + cutname += "SE"; + cutname2 += "SE"; + counter+= "SE"; + } + if(mixing){ + outputfile += "ME"; + contname += "ME"; + cutname += "ME"; + cutname2 += "ME"; + counter+= "ME"; + } + outputfile += "Dphi_DStar"; + outputfile += particle; + cutname += particle; + cutname2 += particle; + contname += particle; + counter+= particle; + + + cout << "Contname = " << contname << endl; + mgr->AddTask(task); + // ------ input data ------ + AliAnalysisDataContainer *cinput0 = mgr->GetCommonInputContainer(); + + // ----- output data ----- + + // output TH1I for event counting + AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(contname, TList::Class(),AliAnalysisManager::kOutputContainer,outputfile.Data()); + AliAnalysisDataContainer *coutput2 = mgr->CreateContainer(cutname,AliRDHFCutsDStartoKpipi::Class(),AliAnalysisManager::kOutputContainer, outputfile.Data()); //cuts + AliAnalysisDataContainer *coutputDstarNorm = mgr->CreateContainer(counter,AliNormalizationCounter::Class(),AliAnalysisManager::kOutputContainer, outputfile.Data()); + AliAnalysisDataContainer *coutput4 = mgr->CreateContainer(cutname2,AliHFAssociatedTrackCuts::Class(),AliAnalysisManager::kOutputContainer, outputfile.Data()); //cuts + + mgr->ConnectInput(task,0,mgr->GetCommonInputContainer()); + mgr->ConnectOutput(task,1,coutput1); + mgr->ConnectOutput(task,2,coutput2); + mgr->ConnectOutput(task,3,coutputDstarNorm); + mgr->ConnectOutput(task,4,coutput4); + + return task ; + +} + diff --git a/PWGHF/correlationHF/macros/makeTFileAssociatedTrackCuts.C b/PWGHF/correlationHF/macros/makeTFileAssociatedTrackCuts.C index 7e194618de9..f9c946faf4b 100644 --- a/PWGHF/correlationHF/macros/makeTFileAssociatedTrackCuts.C +++ b/PWGHF/correlationHF/macros/makeTFileAssociatedTrackCuts.C @@ -1,91 +1,106 @@ -#include -#include -#include "AliHFAssociatedTrackCuts.h" -#include -#include -//#include "AliAODPidHF.h" - -/* $Id$ */ - -void makeInputHFCorrelation(){ - - - AliHFAssociatedTrackCuts* HFCorrelationCuts=new AliHFAssociatedTrackCuts(); - HFCorrelationCuts->SetName("AssociatedCuts"); - HFCorrelationCuts->SetTitle("Cuts for associated track"); - Float_t eta = 0.9; - - //______________________________ set ESD track cuts - AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts(); - - - esdTrackCuts->SetRequireSigmaToVertex(kFALSE); - //esdTrackCuts->SetDCAToVertex2D(kFALSE); - - esdTrackCuts->SetRequireTPCRefit(kTRUE); - esdTrackCuts->SetRequireITSRefit(kTRUE); - esdTrackCuts->SetMinNClustersTPC(80); - esdTrackCuts->SetMinNClustersITS(2); - //esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny); // set requirement on Pixels - // default is kBoth, otherwise kAny - esdTrackCuts->SetPtRange(0.3,1.e10); - esdTrackCuts->SetEtaRange(-eta,eta); - - - HFCorrelationCuts->AddTrackCuts(esdTrackCuts); - - //______________________________ set kinematics cuts for AOD track - const int nofcuts = 4; - Float_t* trackcutsarray; - trackcutsarray=new Float_t[nofcuts]; - trackcutsarray[0] = 0.3;//track min pt - trackcutsarray[1] = 100.;//track max pt - trackcutsarray[2] = 0.;//track min impact parameter - trackcutsarray[3] = 100.;//track max impact parameter - - HFCorrelationCuts->SetNVarsTrack(nofcuts); - HFCorrelationCuts->SetAODTrackCuts(trackcutsarray); - - - - //______________________________ set kinematics cuts for AOD v0 - - const int nofcuts2 = 7; - - Float_t* vzerocutsarray; - vzerocutsarray=new Float_t[nofcuts2]; - vzerocutsarray[0] = 0.2; // max dca between two daugters (cm) - vzerocutsarray[1] = 2; // max chi square - vzerocutsarray[2] = 2.; // min decay length (cm) - vzerocutsarray[3] = 15; // max decay length (cm) - vzerocutsarray[4] = 0.2; // min opening angle between two daugters - vzerocutsarray[5] = 0; // min pt of k0 (GeV/c) - vzerocutsarray[6] = 0.9; // set eta acceptance - - HFCorrelationCuts->SetNVarsVzero(nofcuts2); - HFCorrelationCuts->SetAODvZeroCuts(vzerocutsarray); - - //______________________________ set PID - - Int_t mode =1; - AliAODPidHF* pidObj=new AliAODPidHF(); - pidObj->SetMatch(mode); - pidObj->SetSigma(0,2); // TPC - pidObj->SetSigma(3,3); // TOF - pidObj->SetTPC(kTRUE); - pidObj->SetTOF(kTRUE); - pidObj->SetCompat(kTRUE); - - HFCorrelationCuts->SetPidHF(pidObj); - - - //______________________________ save to *.root file - HFCorrelationCuts->PrintAll(); - - TFile* fout=new TFile("AssocPartCuts.root","recreate"); //set this!! - fout->cd(); - HFCorrelationCuts->Write(); - fout->Close(); - - -} +#include +#include +#include "AliHFAssociatedTrackCuts.h" +#include +#include +#include "AliAODPidHF.h" + +/* $Id$ */ + +void makeInputHFCorrelation(){ + + AliHFAssociatedTrackCuts* HFCorrelationCuts=new AliHFAssociatedTrackCuts(); + HFCorrelationCuts->SetName("AssociatedCuts"); + HFCorrelationCuts->SetTitle("Cuts for associated track"); + Float_t eta = 0.9; + //______________________________ set ESD track cuts + AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts(); + + esdTrackCuts->SetRequireSigmaToVertex(kFALSE); + esdTrackCuts->SetDCAToVertex2D(kFALSE); + + esdTrackCuts->SetRequireTPCRefit(kTRUE); + esdTrackCuts->SetRequireITSRefit(kTRUE); + esdTrackCuts->SetMinNClustersTPC(80); + esdTrackCuts->SetMinNClustersITS(2); + //esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny); // set requirement on Pixels + // default is kBoth, otherwise kAny + esdTrackCuts->SetPtRange(0.3,1.e10); + esdTrackCuts->SetEtaRange(-eta,eta); + HFCorrelationCuts->AddTrackCuts(esdTrackCuts); + //______________________________ set Pool for event mixing + + + + HFCorrelationCuts->SetMaxNEventsInPool(200); + HFCorrelationCuts->SetMinNTracksInPool(1000); + HFCorrelationCuts->SetMinEventsToMix(8); + + HFCorrelationCuts->SetNofPoolBins(3,5); + + //Double_t MBins[]={0,20,40,60,80,500}; + Double_t MBins[]={0,1,2,3,4,500}; + Double_t * MultiplicityBins = MBins; + + //Double_t ZBins[]={-10,-5,-2.5,2.5,5,10}; + Double_t ZBins[]={-10,-2.5,2.5,10}; + Double_t *ZVrtxBins = ZBins; + + HFCorrelationCuts->SetPoolBins(ZVrtxBins,MultiplicityBins); + cout << "Crash 1 " << endl; + //______________________________ set kinematics cuts for AOD track + const int nofcuts = 4; + Float_t* trackcutsarray; + trackcutsarray=new Float_t[nofcuts]; + cout << "Crash 1.1 " << endl; + trackcutsarray[0] = 0.3;//track min pt + trackcutsarray[1] = 100.;//track max pt + trackcutsarray[2] = 0.;//track min impact parameter + trackcutsarray[3] = 100.;//track max impact parameter + cout << "Crash 1.2 " << endl; + HFCorrelationCuts->SetNVarsTrack(nofcuts); + cout << "Crash 1.3 " << endl; + HFCorrelationCuts->SetAODTrackCuts(trackcutsarray); + + cout << "Crash 2 " << endl; + //______________________________ set kinematics cuts for AOD v0 + + const int nofcuts2 = 7; + + Float_t* vzerocutsarray; + vzerocutsarray=new Float_t[nofcuts2]; + vzerocutsarray[0] = 0.2; // max dca between two daugters (cm) + vzerocutsarray[1] = 2; // max chi square + vzerocutsarray[2] = 2.; // min decay length (cm) + vzerocutsarray[3] = 15; // max decay length (cm) + vzerocutsarray[4] = 0.2; // min opening angle between two daugters + vzerocutsarray[5] = 0; // min pt of k0 (GeV/c) + vzerocutsarray[6] = 0.9; // set eta acceptance + + HFCorrelationCuts->SetNVarsVzero(nofcuts2); + HFCorrelationCuts->SetAODvZeroCuts(vzerocutsarray); + //______________________________ set PID + cout << "Crash 3 " << endl; + + Int_t mode =1; + AliAODPidHF* pidObj=new AliAODPidHF(); + pidObj->SetMatch(mode); + pidObj->SetSigma(0,2); // TPC + pidObj->SetSigma(3,3); // TOF + pidObj->SetTPC(kTRUE); + pidObj->SetTOF(kTRUE); + pidObj->SetCompat(kTRUE); + + HFCorrelationCuts->SetPidHF(pidObj); + cout << "Crash 4 " << endl; + + //______________________________ save to *.root file + HFCorrelationCuts->PrintAll(); + + TFile* fout=new TFile("AssocPartCuts.root","recreate"); //set this!! + fout->cd(); + HFCorrelationCuts->Write(); + fout->Close(); + cout << "Crash 5 " << endl; + +} -- 2.43.0