// new in this version (use with caution): vzero event plane, event mixing
#include "TChain.h"
-#include "TTree.h"
#include "TH1.h"
#include "TH1F.h"
#include "TH2F.h"
-#include "TCanvas.h"
#include "TMath.h"
#include "TObjArray.h"
-#include "AliAnalysisTask.h"
#include "AliAnalysisTaskSE.h"
#include "AliAnalysisManager.h"
-#include "AliESDEvent.h"
#include "AliAODEvent.h"
-#include "AliESDInputHandler.h"
#include "AliAODInputHandler.h"
#include "AliCentrality.h"
-#include "AliVEvent.h"
#include "AliAnalysisTaskPhiFlow.h"
#include "AliFlowBayesianPID.h"
-#include "AliStack.h"
-#include "AliESDtrackCuts.h"
+#include "AliPIDCombined.h"
#include "AliMCEvent.h"
#include "TProfile.h"
#include "AliFlowCandidateTrack.h"
#include "AliFlowCommonConstants.h"
#include "AliFlowEvent.h"
#include "TVector3.h"
-#include "AliESDVZERO.h"
#include "AliAODVZERO.h"
-#include "AliPID.h"
#include "AliPIDResponse.h"
#include "AliAODMCParticle.h"
#include "AliAnalysisTaskVnV0.h"
ClassImp(AliPhiMesonHelperTrack)
AliAnalysisTaskPhiFlow::AliAnalysisTaskPhiFlow() : AliAnalysisTaskSE(),
- fDebug(0), fIsMC(0), fEventMixing(0), fTypeMixing(0), fQA(0), fV0(0), fAODAnalysis(0), fMassBins(1), fMinMass(-1.), fMaxMass(0.), fCutsRP(NULL), fNullCuts(0), fPIDResponse(0), fFlowEvent(0), fBayesianResponse(0), fCandidates(0), fOldTrackParam(0), fRequireTPCStandAlone(0), fCandidateEtaPtCut(0), fCandidateMinEta(0), fCandidateMaxEta(0), fCandidateMinPt(0), fCandidateMaxPt(0), fNPtBins(18), fCentrality(999), fVertex(999), fESD(0), fAOD(0), fPoolManager(0), fOutputList(0), fEventStats(0), fCentralityPass(0), fCentralityNoPass(0), fNOPID(0), fPIDk(0),fNOPIDTOF(0), fPIDTOF(0), fPtP(0), fPtN(0), fPtKP(0), fPtKN(0), fCentralityMin(0), fCentralityMax(100), fkCentralityMethod(0), fPOICuts(0), fVertexRange(0), fPhi(0), fPt(0), fEta(0), fVZEROA(0), fVZEROC(0), fTPCM(0), fDeltaDipAngle(0), fDeltaDipPt(0), fApplyDeltaDipCut(0), fDCAAll(0), fDCAXYQA(0), fDCAZQA(0), fDCAPrim(0), fDCASecondaryWeak(0), fDCAMaterial(0), fSubEventDPhiv2(0)
+ fDebug(0), fIsMC(0), fEventMixing(0), fTypeMixing(0), fQA(0), fV0(0), fMassBins(1), fMinMass(-1.), fMaxMass(0.), fCutsRP(NULL), fNullCuts(0), fPIDResponse(0), fFlowEvent(0), fBayesianResponse(0), fCandidates(0), fCandidateEtaPtCut(0), fCandidateMinEta(0), fCandidateMaxEta(0), fCandidateMinPt(0), fCandidateMaxPt(0), fNPtBins(18), fCentrality(999), fVertex(999), fAOD(0), fPoolManager(0), fOutputList(0), fEventStats(0), fCentralityPass(0), fCentralityNoPass(0), fNOPID(0), fPIDk(0),fNOPIDTOF(0), fPIDTOF(0), fPtP(0), fPtN(0), fPtKP(0), fPtKN(0), fMultCorAfterCuts(0), fMultvsCentr(0), fCentralityMin(0), fCentralityMax(100), fkCentralityMethodA(0), fkCentralityMethodB(0), fCentralityCut2010(0), fCentralityCut2011(0), fPOICuts(0), fVertexRange(0), fPhi(0), fPt(0), fEta(0), fVZEROA(0), fVZEROC(0), fTPCM(0)/*, fDeltaDipAngle(0), fDeltaDipPt(0), fApplyDeltaDipCut(0)*/, fDCAAll(0), fDCAXYQA(0), fDCAZQA(0), fDCAPrim(0), fDCASecondaryWeak(0), fDCAMaterial(0), fSubEventDPhiv2(0), fSkipEventSelection(0), fUsePidResponse(0), fPIDCombined(0)
{
// Default constructor
for(Int_t i(0); i < 7; i++) fPIDConfig[i] = 0.;
}
//_____________________________________________________________________________
AliAnalysisTaskPhiFlow::AliAnalysisTaskPhiFlow(const char *name) : AliAnalysisTaskSE(name),
- fDebug(0), fIsMC(0), fEventMixing(0), fTypeMixing(0), fQA(0), fV0(0), fAODAnalysis(0), fMassBins(1), fMinMass(-1.), fMaxMass(0.), fCutsRP(NULL), fNullCuts(0), fPIDResponse(0), fFlowEvent(0), fBayesianResponse(0), fCandidates(0), fOldTrackParam(0), fRequireTPCStandAlone(0), fCandidateEtaPtCut(0), fCandidateMinEta(0), fCandidateMaxEta(0), fCandidateMinPt(0), fCandidateMaxPt(0), fNPtBins(18), fCentrality(999), fVertex(999), fESD(0), fAOD(0), fPoolManager(0), fOutputList(0), fEventStats(0), fCentralityPass(0), fCentralityNoPass(0), fNOPID(0), fPIDk(0), fNOPIDTOF(0), fPIDTOF(0), fPtP(0), fPtN(0), fPtKP(0), fPtKN(0), fCentralityMin(0), fCentralityMax(100), fkCentralityMethod(0), fPOICuts(0), fVertexRange(0), fPhi(0), fPt(0), fEta(0), fVZEROA(0), fVZEROC(0), fTPCM(0), fDeltaDipAngle(0), fDeltaDipPt(0), fApplyDeltaDipCut(0), fDCAAll(0), fDCAXYQA(0), fDCAZQA(0), fDCAPrim(0), fDCASecondaryWeak(0), fDCAMaterial(0), fSubEventDPhiv2(0)
+ fDebug(0), fIsMC(0), fEventMixing(0), fTypeMixing(0), fQA(0), fV0(0), fMassBins(1), fMinMass(-1.), fMaxMass(0.), fCutsRP(NULL), fNullCuts(0), fPIDResponse(0), fFlowEvent(0), fBayesianResponse(0), fCandidates(0), fCandidateEtaPtCut(0), fCandidateMinEta(0), fCandidateMaxEta(0), fCandidateMinPt(0), fCandidateMaxPt(0), fNPtBins(18), fCentrality(999), fVertex(999), fAOD(0), fPoolManager(0), fOutputList(0), fEventStats(0), fCentralityPass(0), fCentralityNoPass(0), fNOPID(0), fPIDk(0), fNOPIDTOF(0), fPIDTOF(0), fPtP(0), fPtN(0), fPtKP(0), fPtKN(0), fMultCorAfterCuts(0), fMultvsCentr(0), fCentralityMin(0), fCentralityMax(100), fkCentralityMethodA(0), fkCentralityMethodB(0), fCentralityCut2010(0), fCentralityCut2011(0), fPOICuts(0), fVertexRange(0), fPhi(0), fPt(0), fEta(0), fVZEROA(0), fVZEROC(0), fTPCM(0)/*, fDeltaDipAngle(0), fDeltaDipPt(0), fApplyDeltaDipCut(0)*/, fDCAAll(0), fDCAXYQA(0), fDCAZQA(0), fDCAPrim(0), fDCASecondaryWeak(0), fDCAMaterial(0), fSubEventDPhiv2(0), fSkipEventSelection(0), fUsePidResponse(0), fPIDCombined(0)
{
// Constructor
for(Int_t i(0); i < 7; i++) fPIDConfig[i] = 0.;
if (fFlowEvent) delete fFlowEvent;
if (fBayesianResponse) delete fBayesianResponse;
if (fEventMixing) delete fPoolManager;
+ if (fPIDCombined) delete fPIDCombined;
if (fDebug > 0) cout << " === Deleted instance of AliAnalysisTaskPhiFlow === " << endl;
}
//_____________________________________________________________________________
fOutputList->Add(fDCAXYQA);
fDCAZQA = new TH1F("fDCAZQA", "fDCAZQA", 1000, -5, 5);
fOutputList->Add(fDCAZQA);
+ if(fCentralityCut2010 || fCentralityCut2011) {
+ fMultCorAfterCuts = new TH2F("fMultCorAfterCuts", "TPC vs Global multiplicity (After cuts); Global multiplicity; TPC multiplicity", 100, 0, 3000, 100, 0, 3000);
+ fOutputList->Add(fMultCorAfterCuts);
+ fMultvsCentr = new TH2F("fMultvsCentr", "Multiplicity vs centrality; centrality; Multiplicity", 9, -0.5, 100.5, 101, 0, 3000);
+ fOutputList->Add(fMultvsCentr);
+ }
}
if(fIsMC || fQA) {
fDCAAll = new TH2F("fDCAAll", "fDCAAll", 1000, 0, 10, 1000, -5, 5);
if(fDebug > 0) cout << " *** UserCreateOutputObjects() *** " << endl;
fNullCuts = new AliFlowTrackCuts("null_cuts");
fBayesianResponse = new AliFlowBayesianPID();
+ // combined pid
+ fPIDCombined = new AliPIDCombined;
+ fPIDCombined->SetDefaultTPCPriors();
+ fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC|AliPIDResponse::kDetTOF);
+
+ // flag to mc
+ if(fSkipEventSelection || fIsMC) fBayesianResponse->SetMC(kTRUE);
Double_t t(0);
for(Int_t i = 0; i < 7; i++) t+=TMath::Abs(fPIDConfig[i]);
if(t < 0.1) AliFatal("No valid PID procedure recognized -- terminating analysis !!!");
cc->SetNbinsMass(fMassBins); cc->SetNbinsEta(200); (fMassBins == 1) ? cc->SetNbinsPt(15) : cc->SetNbinsPt(100); // high pt
cc->SetMassMin(fMinMass); cc->SetEtaMin(-5.0); cc->SetPtMin(0);
cc->SetMassMax(fMaxMass); cc->SetEtaMax(+5.0); (fMassBins == 1) ? cc->SetPtMax(15) : cc->SetPtMax(10); // high pt
- if (!fOldTrackParam) fBayesianResponse->SetNewTrackParam();
+ fBayesianResponse->SetNewTrackParam();
AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
if (man) {
AliInputEventHandler* inputHandler = (AliInputEventHandler*)(man->GetInputEventHandler());
return TMath::Sqrt((es - (pxs + pys + pzs)));
}
//_____________________________________________________________________________
+/*
template <typename T> Double_t AliAnalysisTaskPhiFlow::DeltaDipAngle(const T* track1, const T* track2) const
{
// Calculate the delta dip angle between two particles
if ((TMath::Abs(DeltaDipAngle(track1, track2)) < fDeltaDipAngle) && (PhiPt(track1, track2) < fDeltaDipPt)) return kFALSE;
return kTRUE;
}
+*/
//_____________________________________________________________________________
template <typename T> Bool_t AliAnalysisTaskPhiFlow::CheckCandidateEtaPtCut(const T* track1, const T* track2) const
{
// Impose event cuts
if(fDebug > 0) cout << " *** EventCut() *** " << endl;
if (!event) return kFALSE;
+ if (fSkipEventSelection) return kTRUE;
if (!CheckVertex(event)) return kFALSE;
if (!CheckCentrality(event)) return kFALSE;
if(fQA) PlotMultiplcities(event);
{
// Check if event is within the set centrality range. Falls back to V0 centrality determination if no method is set
if(fDebug > 0) cout << " *** CheckCentrality() *** " << endl;
- if (!fkCentralityMethod) AliFatal("No centrality method set! FATAL ERROR!");
- fCentrality = event->GetCentrality()->GetCentralityPercentile(fkCentralityMethod);
- if ((fCentrality <= fCentralityMin) || (fCentrality > fCentralityMax)) {
+ if (!fkCentralityMethodA) AliFatal("No centrality method set! FATAL ERROR!");
+ fCentrality = event->GetCentrality()->GetCentralityPercentile(fkCentralityMethodA);
+ Double_t cenB(-999);
+ // if a second centrality estimator is requited, set it
+ (fkCentralityMethodB) ? cenB = event->GetCentrality()->GetCentralityPercentile(fkCentralityMethodB) : cenB = fCentrality;
+ if (TMath::Abs(fCentrality-cenB) > 5 || cenB >= 80 || cenB < 0 || fCentrality <= fCentralityMin || fCentrality > fCentralityMax) {
if(fQA) fCentralityNoPass->Fill(fCentrality) ;
return kFALSE;
}
+ const Int_t nGoodTracks = event->GetNumberOfTracks();
+ if(fCentralityCut2010) { // cut on outliers
+ Float_t multTPC(0.); // tpc mult estimate
+ Float_t multGlob(0.); // global multiplicity
+ for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) { // fill tpc mult
+ AliAODTrack* trackAOD = event->GetTrack(iTracks);
+ if (!trackAOD) continue;
+ if (!(trackAOD->TestFilterBit(1))) continue;
+ if ((trackAOD->Pt() < .2) || (trackAOD->Pt() > 5.0) || (TMath::Abs(trackAOD->Eta()) > .8) || (trackAOD->GetTPCNcls() < 70) || (trackAOD->GetDetPid()->GetTPCsignal() < 10.0) || (trackAOD->Chi2perNDF() < 0.2)) continue;
+ multTPC++;
+ }
+ for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) { // fill global mult
+ AliAODTrack* trackAOD = event->GetTrack(iTracks);
+ if (!trackAOD) continue;
+ if (!(trackAOD->TestFilterBit(16))) continue;
+ if ((trackAOD->Pt() < .2) || (trackAOD->Pt() > 5.0) || (TMath::Abs(trackAOD->Eta()) > .8) || (trackAOD->GetTPCNcls() < 70) || (trackAOD->GetDetPid()->GetTPCsignal() < 10.0) || (trackAOD->Chi2perNDF() < 0.1)) continue;
+ Double_t b[2] = {-99., -99.};
+ Double_t bCov[3] = {-99., -99., -99.};
+ AliAODTrack copy(*trackAOD);
+ if (!(copy.PropagateToDCA(event->GetPrimaryVertex(), event->GetMagneticField(), 100., b, bCov))) continue;
+ if ((TMath::Abs(b[0]) > 0.3) || (TMath::Abs(b[1]) > 0.3)) continue;
+ multGlob++;
+ } //track loop
+ // printf(" mult TPC %.2f, mult Glob %.2f \n", multTPC, multGlob);
+ if(! (multTPC > (-40.3+1.22*multGlob) && multTPC < (32.1+1.59*multGlob))) return kFALSE;
+ if(fQA) {
+ fMultCorAfterCuts->Fill(multGlob, multTPC);
+ fMultvsCentr->Fill(fCentrality, multTPC);
+ }
+ }
+
+ if(fCentralityCut2011) { // cut on outliers
+ Float_t multTPC(0.); // tpc mult estimate
+ Float_t multGlob(0.); // global multiplicity
+ for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) { // fill tpc mult
+ AliAODTrack* trackAOD = event->GetTrack(iTracks);
+ if (!trackAOD) continue;
+ if (!(trackAOD->TestFilterBit(1))) continue;
+ if ((trackAOD->Pt() < .2) || (trackAOD->Pt() > 5.0) || (TMath::Abs(trackAOD->Eta()) > .8) || (trackAOD->GetTPCNcls() < 70) || (trackAOD->GetDetPid()->GetTPCsignal() < 10.0) || (trackAOD->Chi2perNDF() < 0.2)) continue;
+ multTPC++;
+ }
+ for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) { // fill global mult
+ AliAODTrack* trackAOD = event->GetTrack(iTracks);
+ if (!trackAOD) continue;
+ if (!(trackAOD->TestFilterBit(16))) continue;
+ if ((trackAOD->Pt() < .2) || (trackAOD->Pt() > 5.0) || (TMath::Abs(trackAOD->Eta()) > .8) || (trackAOD->GetTPCNcls() < 70) || (trackAOD->GetDetPid()->GetTPCsignal() < 10.0) || (trackAOD->Chi2perNDF() < 0.1)) continue;
+ Double_t b[2] = {-99., -99.};
+ Double_t bCov[3] = {-99., -99., -99.};
+ AliAODTrack copy(*trackAOD);
+ if (!(copy.PropagateToDCA(event->GetPrimaryVertex(), event->GetMagneticField(), 100., b, bCov))) continue;
+ if ((TMath::Abs(b[0]) > 0.3) || (TMath::Abs(b[1]) > 0.3)) continue;
+ multGlob++;
+ } //track loop
+ //printf(" mult TPC %.2f, mult Glob %.2f \n", multTPC, multGlob);
+ if(! (multTPC > (-36.73 + 1.48*multGlob) && multTPC < (62.87 + 1.78*multGlob))) return kFALSE;
+ if(fQA) {
+ fMultCorAfterCuts->Fill(multGlob, multTPC);
+ fMultvsCentr->Fill(fCentrality, multTPC);
+ }
+ }
+
fCentralityPass->Fill(fCentrality);
return kTRUE;
}
//_____________________________________________________________________________
-void AliAnalysisTaskPhiFlow::InitializeBayesianPID(AliESDEvent* event)
-{
- // Initialize the Bayesian PID object for ESD analysis
- if(fDebug > 0) cout << " *** InitializeBayesianPID() *** " << endl;
- fBayesianResponse->SetDetResponse(event, fCentrality, AliESDpid::kTOF_T0, kTRUE);
-}
-//_____________________________________________________________________________
void AliAnalysisTaskPhiFlow::InitializeBayesianPID(AliAODEvent* event)
{
// Initialize the Bayesian PID object for AOD
// fDCAConfig[0] > 1 pt dependent dca cut
if(fDebug > 1) cout << " *** PassesDCACut() *** " << endl;
if(fIsMC) return kTRUE;
+ AliAODTrack copy(*track);
if( (fDCAConfig[0] < 0.1) && (fDCAConfig[0] > -0.1) ) {
if(fQA) {
Double_t b[2] = { -99., -99.};
Double_t bCov[3] = { -99., -99., -99.};
- track->PropagateToDCA(fAOD->GetPrimaryVertex(), fAOD->GetMagneticField(), 100., b, bCov);
- fDCAXYQA->Fill(b[0]);
- fDCAZQA->Fill(b[1]);
- fDCAPrim->Fill(track->Pt(), b[0]);
+ if(copy.PropagateToDCA(fAOD->GetPrimaryVertex(), fAOD->GetMagneticField(), 100., b, bCov)) {
+ fDCAXYQA->Fill(b[0]);
+ fDCAZQA->Fill(b[1]);
+ fDCAPrim->Fill(track->Pt(), b[0]);
+ }
}
return kTRUE;
}
Double_t b[2] = { -99., -99.};
Double_t bCov[3] = { -99., -99., -99.};
- track->PropagateToDCA(fAOD->GetPrimaryVertex(), fAOD->GetMagneticField(), 100., b, bCov);
+ if(!copy.PropagateToDCA(fAOD->GetPrimaryVertex(), fAOD->GetMagneticField(), 100., b, bCov)) return kFALSE;
if((!fIsMC)&&fQA) fDCAMaterial->Fill(track->Pt(), b[0]);
if( (fDCAConfig[0] < -.9) && ( (TMath::Abs(b[0]) > fDCAConfig[1]) || (TMath::Abs(b[1]) > fDCAConfig[2])) ) return kFALSE;
if(fDCAConfig[0] > .9) {
return kTRUE;
}
//_____________________________________________________________________________
-Bool_t AliAnalysisTaskPhiFlow::IsKaon(AliESDtrack* track) const
-{
- // Check if particle is a kaon according to method set in steering macro
- if(fDebug > 1) cout << " *** IsKaon() *** " << endl;
- if (fRequireTPCStandAlone && (track->GetStatus()&AliESDtrack::kTPCin) == 0) return kFALSE;
- if (fPIDConfig[1]!=0 || fPIDConfig[4]!=0) AliFatal("TPC || ITS PID not available in ESD anlaysis -- terminating analysis !!!");
- if (PassesTPCbayesianCut(track))
- {
- if(fQA) fPIDk->Fill(track->P(), track->GetTPCsignal());
- return kTRUE;
- }
- return kFALSE;
-}
-//_____________________________________________________________________________
Bool_t AliAnalysisTaskPhiFlow::IsKaon(AliAODTrack* track) const
{
// Kaon identification routine, based on multiple detectors and approaches
if(fDebug > 1) cout << " *** IsKaon() *** " << endl;
- if(fRequireTPCStandAlone && (!track->TestFilterBit(1))) return kFALSE;
+ if(fUsePidResponse) {
+ Double_t prob[10] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
+ fPIDCombined->ComputeProbabilities(track, fPIDResponse, prob);
+ if(prob[3] > fPIDConfig[6]) return kTRUE;
+ }
if(!PassesDCACut(track)) return kFALSE;
if(fQA) {fNOPID->Fill(track->P(), track->GetTPCsignal());fNOPIDTOF->Fill(track->P(), track->GetTOFsignal());}
if(track->Pt() < fPIDConfig[1]) {
}
fAOD = dynamic_cast<AliAODEvent*>(InputEvent()); // check for aod data type
if (fAOD) {
- fAODAnalysis = kTRUE;
if (!EventCut(fAOD)) return; // check for event cuts
InitializeBayesianPID(fAOD); // init event objects
SetNullCuts(fAOD);
}
for (Int_t pTracks = 0; pTracks < unp ; pTracks++) { // perform analysis
for (Int_t nTracks = 0; nTracks < unn ; nTracks++) {
- if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(up[pTracks], un[nTracks]))) continue;
+// if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(up[pTracks], un[nTracks]))) continue;
if (fCandidateEtaPtCut && (!CheckCandidateEtaPtCut(up[pTracks], un[nTracks]))) continue;
PtSelector(0, up[pTracks], un[nTracks]);
Double_t pt = PhiPt(up[pTracks], un[nTracks]);
if(!fEventMixing) { // combinatorial background
for (Int_t pTracks = 0; pTracks < unp ; pTracks++) {
for (Int_t nTracks = pTracks + 1; nTracks < unp ; nTracks++) {
- if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(up[pTracks], up[nTracks]))) continue;
+// if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(up[pTracks], up[nTracks]))) continue;
if (fCandidateEtaPtCut && (!CheckCandidateEtaPtCut(up[pTracks], up[nTracks]))) continue;
PtSelector(1, up[pTracks], up[nTracks]);
}
}
for (Int_t nTracks = 0; nTracks < unn ; nTracks++) {
for (Int_t pTracks = nTracks + 1; pTracks < unn ; pTracks++) {
- if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(un[nTracks], un[pTracks]))) continue;
+// if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(un[nTracks], un[pTracks]))) continue;
if (fCandidateEtaPtCut && (!CheckCandidateEtaPtCut(un[nTracks], un[pTracks]))) continue;
PtSelector(2, un[nTracks], un[pTracks]);
}
PostData(1, fOutputList);
PostData(2, fFlowEvent);
}
- fESD = dynamic_cast<AliESDEvent*>(InputEvent()); // check for esd data type (obsolete)
- if (fESD) {
- if(fV0 || fEventMixing) {
- if (fDebug > 0 ) cout << " FATAL ERROR: AOD only analysis called on ESD data, skipping input !!! " << endl;
- return;
- }
- fAODAnalysis = kFALSE;
- if (!EventCut(fESD)) return; // check event cuts
- InitializeBayesianPID(fESD);
- SetNullCuts(fESD);
- PrepareFlowEvent(fESD->GetNumberOfTracks());
- if(fQA) fEventStats->Fill(0);
- Int_t unTracks = fESD->GetNumberOfTracks();
- AliESDtrack* un[unTracks];
- AliESDtrack* up[unTracks];
- Int_t unp(0);
- Int_t unn(0);
- for (Int_t iTracks = 0; iTracks < unTracks; iTracks++) { // select analysis candidates
- AliESDtrack* track = fESD->GetTrack(iTracks);
- if (!PhiTrack(track)) continue;
- Bool_t charge = kFALSE;
- if (track->Charge() > 0) {
- charge = kTRUE;
- if(fQA) {fEventStats->Fill(1); fPtP->Fill(track->Pt());}
- }
- if (track->Charge() < 0) {
- if(fQA) fEventStats->Fill(2);
- fPtN->Fill(track->Pt());
- }
- if (IsKaon(track)) {
- if (charge) {
- up[unp] = track;
- unp++;
- if(fQA) {fEventStats->Fill(3);fPtKP->Fill(track->Pt());}
- }
- if (!charge) {
- un[unn] = track;
- unn++;
- if(fQA) {fEventStats->Fill(4); fPtKN->Fill(track->Pt());}
- }
- }
- }
- for (Int_t pTracks = 0; pTracks < unp ; pTracks++) { // perform analysis
- for (Int_t nTracks = 0; nTracks < unn ; nTracks++) {
- if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(up[pTracks], un[nTracks]))) continue;
- if (fCandidateEtaPtCut && (!CheckCandidateEtaPtCut(up[pTracks], un[nTracks]))) continue;
- PtSelector(0, up[pTracks], un[nTracks]);
- Double_t pt = PhiPt(up[pTracks], un[nTracks]);
- Double_t mass = InvariantMass(up[pTracks], un[nTracks]);
- TVector3 a(up[pTracks]->Px(), up[pTracks]->Py(), up[pTracks]->Pz());
- TVector3 b(un[nTracks]->Px(), un[nTracks]->Py(), up[pTracks]->Pz());
- TVector3 c = a + b;
- Double_t phi = c.Phi();
- Double_t eta = c.Eta();
- Int_t nIDs[2];
- nIDs[0] = up[pTracks]->GetID();
- nIDs[1] = un[nTracks]->GetID();
- MakeTrack(mass, pt, phi, eta, 2, nIDs);
- }
- }
- if (fDebug > 0) printf("I received %d candidates\n", fCandidates->GetEntriesFast()); // insert candidates in flow event
- for (int iCand = 0; iCand != fCandidates->GetEntriesFast(); ++iCand) {
- AliFlowCandidateTrack *cand = dynamic_cast<AliFlowCandidateTrack*>(fCandidates->At(iCand));
- if (!cand) continue;
- if (fDebug > 1) printf(" >Checking at candidate %d with %d daughters: mass %f\n", iCand, cand->GetNDaughters(), cand->Mass());
- for (int iDau = 0; iDau != cand->GetNDaughters(); ++iDau) {
- if (fDebug > 1) printf(" >Daughter %d with fID %d", iDau, cand->GetIDDaughter(iDau));
- for (int iRPs = 0; iRPs != fFlowEvent->NumberOfTracks(); ++iRPs) {
- AliFlowTrack *iRP = dynamic_cast<AliFlowTrack*>(fFlowEvent->GetTrack(iRPs));
- if (!iRP) continue;
- if (!iRP->InRPSelection()) continue;
- if (cand->GetIDDaughter(iDau) == iRP->GetID()) {
- if (fDebug>1) printf(" was in RP set");
- iRP->SetForRPSelection(kFALSE);
- fFlowEvent->SetNumberOfRPs(fFlowEvent->GetNumberOfRPs() - 1);
- }
- }
- if (fDebug>1) printf("\n");
- }
- cand->SetForPOISelection(kTRUE);
- fFlowEvent->InsertTrack(((AliFlowTrack*) cand));
- }
- if (fDebug > 0) printf("TPCevent %d\n", fFlowEvent->NumberOfTracks());
-
- for (Int_t pTracks = 0; pTracks < unp ; pTracks++) {
- for (Int_t nTracks = pTracks + 1; nTracks < unp ; nTracks++) {
- if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(up[pTracks], up[nTracks]))) continue;
- if (fCandidateMinEta && (!CheckCandidateEtaPtCut(up[pTracks], up[nTracks]))) continue;
- PtSelector(1, up[pTracks], up[nTracks]);
- }
- }
- for (Int_t nTracks = 0; nTracks < unn ; nTracks++) {
- for (Int_t pTracks = nTracks + 1; pTracks < unn ; pTracks++) {
- if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(un[nTracks], un[pTracks]))) continue;
- if (fCandidateEtaPtCut && (!CheckCandidateEtaPtCut(un[nTracks], un[pTracks]))) continue;
- PtSelector(2, un[nTracks], un[pTracks]);
- }
- }
- PostData(1, fOutputList);
- PostData(2, fFlowEvent);
- }
+}
+//_____________________________________________________________________________
+void AliAnalysisTaskPhiFlow::Exec(Option_t* c)
+{
+ // skip the event selection for SE task (e.g. for MC productions)
+ if(fSkipEventSelection) AliAnalysisTaskPhiFlow::UserExec(c);
+ // exec of task se will do event selection and call UserExec
+ else AliAnalysisTaskSE::Exec(c);
}
//_____________________________________________________________________________
void AliAnalysisTaskPhiFlow::ReconstructionWithEventMixing(TObjArray* MixingCandidates) const
if(!fTypeMixing) { // mix with like-sign kaons
for(Int_t pBuf = 0; pBuf < buffer_unp; pBuf++) {
if(fDebug > 1 ) cout << Form(" buffered k+ track %d", pBuf) << endl;
- if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(mix_up[pMix], buffer_up[pBuf]))) continue;
+// if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(mix_up[pMix], buffer_up[pBuf]))) continue;
if (fCandidateMinEta && (!CheckCandidateEtaPtCut(mix_up[pMix], buffer_up[pBuf]))) continue;
PtSelector(1, mix_up[pMix], buffer_up[pBuf]);
}
else if(fTypeMixing) { // mix with unlike-sign kaons
for(Int_t nBuf = 0; nBuf < buffer_unn; nBuf++) {
if(fDebug > 1 ) cout << Form(" buffered k- track %d", nBuf) << endl;
- if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(mix_up[pMix], buffer_un[nBuf]))) continue;
+// if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(mix_up[pMix], buffer_un[nBuf]))) continue;
if (fCandidateMinEta && (!CheckCandidateEtaPtCut(mix_up[pMix], buffer_un[nBuf]))) continue;
PtSelector(2, mix_up[pMix], buffer_un[nBuf]);
}
if(!fTypeMixing) { // mix with like-sign kaons
for(Int_t nBuf = 0; nBuf < buffer_unn; nBuf++) {
if(fDebug > 1 ) cout << Form(" buffered k- track %d", nBuf) << endl;
- if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(mix_un[nMix], buffer_un[nBuf]))) continue;
+// if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(mix_un[nMix], buffer_un[nBuf]))) continue;
if (fCandidateMinEta && (!CheckCandidateEtaPtCut(mix_un[nMix], buffer_un[nBuf]))) continue;
PtSelector(2, mix_un[nMix], buffer_un[nBuf]);
}
else if(fTypeMixing) { // mix with unlike-sign kaons
for(Int_t pBuf = 0; pBuf < buffer_unp; pBuf++) {
if(fDebug > 1 ) cout << Form(" buffered k+ track %d", pBuf) << endl;
- if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(mix_un[nMix], buffer_up[pBuf]))) continue;
+// if (fApplyDeltaDipCut && (!CheckDeltaDipAngle(mix_un[nMix], buffer_up[pBuf]))) continue;
if (fCandidateMinEta && (!CheckCandidateEtaPtCut(mix_un[nMix], buffer_up[pBuf]))) continue;
PtSelector(1, mix_un[nMix], buffer_up[pBuf]);
}
if (fDebug>1) cout << " Received MC kaon " << endl;
Double_t b[2] = { -99., -99.};
Double_t bCov[3] = { -99., -99., -99.};
- track->PropagateToDCA(fAOD->GetPrimaryVertex(), fAOD->GetMagneticField(), 100., b, bCov);
+ AliAODTrack copy(*track);
+ if(!copy.PropagateToDCA(fAOD->GetPrimaryVertex(), fAOD->GetMagneticField(), 100., b, bCov)) return;
// find corresponding mc particle
AliAODMCParticle *partMC = (AliAODMCParticle*) arrayMC->At(TMath::Abs(track->GetLabel()));
if (!partMC) {