#include "TRandom3.h"\r
#include "TProfile.h"\r
#include "TProfile2D.h"\r
+#include "TLorentzVector.h"\r
+#include "TParticle.h"\r
\r
+#include "AliVEventHandler.h"\r
#include "AliAnalysisTaskSE.h"\r
+#include "AliAnalysisManager.h"\r
\r
#include "AliVEvent.h"\r
#include "AliESDInputHandler.h"\r
#include "AliMCEvent.h"\r
#include "AliESD.h"\r
#include "AliESDEvent.h"\r
+#include "AliPID.h"\r
#include "AliPIDResponse.h"\r
#include "AliESDVZERO.h"\r
#include "AliESDUtils.h"\r
#include "AliMCParticle.h"\r
+#include "AliAODMCParticle.h"\r
+#include "AliAODEvent.h"\r
+#include "AliAODVertex.h"\r
+#include "AliAODTrack.h"\r
#include "AliVTrack.h"\r
+#include "AliESDtrack.h"\r
+#include "AliESDtrackCuts.h"\r
+#include "AliAODTrack.h"\r
+#include "AliStack.h"\r
+#include "AliMCEvent.h"\r
\r
#include "AliFlowCandidateTrack.h"\r
#include "AliFlowEvent.h"\r
#include "AliFlowTrackCuts.h"\r
#include "AliFlowVector.h"\r
#include "AliFlowCommonConstants.h"\r
+#include "AliKFParticle.h"\r
+#include "AliKFVertex.h"\r
\r
#include "AliHFEcuts.h"\r
#include "AliHFEpid.h"\r
AliAnalysisTaskSE(),\r
fListHist(0x0), \r
fAODAnalysis(kFALSE),\r
+ fUseFlagAOD(kFALSE),\r
+ fApplyCut(kTRUE),\r
+ fFlags(1<<4),\r
fVZEROEventPlane(kFALSE),\r
fVZEROEventPlaneA(kFALSE),\r
fVZEROEventPlaneC(kFALSE),\r
fUseMCReactionPlane(kFALSE),\r
fMCPID(kFALSE),\r
fNoPID(kFALSE),\r
+ fChi2OverNDFCut(3.0),\r
+ fMaxdca(3.0),\r
+ fMaxopeningtheta(0.02),\r
+ fMaxopeningphi(0.1),\r
+ fMaxopening3D(0.1),\r
+ fMaxInvmass(0.1),\r
+ fSetMassConstraint(kFALSE),\r
fDebugLevel(0),\r
fcutsRP(0),\r
fcutsPOI(0),\r
fHFECuts(0),\r
fPID(0),\r
+ fPIDTOFOnly(0),\r
fPIDqa(0),\r
- fflowEvent(0x0),\r
+ fflowEvent(NULL),\r
+ fHFEBackgroundCuts(0),\r
+ fPIDBackground(0),\r
+ fPIDBackgroundqa(0),\r
+ fAlgorithmMA(kTRUE),\r
+ fArraytrack(NULL),\r
+ fCounterPoolBackground(0),\r
fHFEVZEROEventPlane(0x0),\r
fHistEV(0),\r
fEventPlane(0x0),\r
fCosRes(0x0),\r
fSinRes(0x0),\r
fProfileCosRes(0x0),\r
+ fTrackingCuts(0x0),\r
+ fDeltaPhiMapsBeforePID(0x0),\r
+ fCosPhiMapsBeforePID(0x0),\r
fDeltaPhiMaps(0x0),\r
+ fDeltaPhiMapsContamination(0x0),\r
fCosPhiMaps(0x0),\r
- fProfileCosPhiMaps(0x0)\r
+ fProfileCosPhiMaps(0x0),\r
+ fDeltaPhiMapsTaggedPhotonic(0x0),\r
+ fCosPhiMapsTaggedPhotonic(0x0),\r
+ fDeltaPhiMapsTaggedNonPhotonic(0x0),\r
+ fCosPhiMapsTaggedNonPhotonic(0x0),\r
+ fDeltaPhiMapsTaggedPhotonicLS(0x0),\r
+ fCosPhiMapsTaggedPhotonicLS(0x0),\r
+ fMCSourceDeltaPhiMaps(0x0),\r
+ fOppSignDeltaPhiMaps(0x0),\r
+ fSameSignDeltaPhiMaps(0x0),\r
+ fOppSignAngle(0x0),\r
+ fSameSignAngle(0x0)\r
{\r
// Constructor\r
\r
AliAnalysisTaskHFEFlow:: AliAnalysisTaskHFEFlow(const char *name) :\r
AliAnalysisTaskSE(name),\r
fListHist(0x0),\r
- fAODAnalysis(kFALSE), \r
+ fAODAnalysis(kFALSE),\r
+ fUseFlagAOD(kFALSE),\r
+ fApplyCut(kTRUE),\r
+ fFlags(1<<4), \r
fVZEROEventPlane(kFALSE),\r
fVZEROEventPlaneA(kFALSE),\r
fVZEROEventPlaneC(kFALSE),\r
fUseMCReactionPlane(kFALSE),\r
fMCPID(kFALSE),\r
fNoPID(kFALSE),\r
+ fChi2OverNDFCut(3.0),\r
+ fMaxdca(3.0),\r
+ fMaxopeningtheta(0.02),\r
+ fMaxopeningphi(0.1),\r
+ fMaxopening3D(0.1),\r
+ fMaxInvmass(0.1),\r
+ fSetMassConstraint(kFALSE),\r
fDebugLevel(0),\r
fcutsRP(0),\r
fcutsPOI(0),\r
fHFECuts(0),\r
fPID(0),\r
+ fPIDTOFOnly(0),\r
fPIDqa(0),\r
- fflowEvent(0x0),\r
+ fflowEvent(NULL),\r
+ fHFEBackgroundCuts(0),\r
+ fPIDBackground(0),\r
+ fPIDBackgroundqa(0),\r
+ fAlgorithmMA(kTRUE), \r
+ fArraytrack(NULL),\r
+ fCounterPoolBackground(0),\r
fHFEVZEROEventPlane(0x0),\r
fHistEV(0),\r
fEventPlane(0x0),\r
fCosRes(0x0),\r
fSinRes(0x0),\r
fProfileCosRes(0x0),\r
+ fTrackingCuts(0x0),\r
+ fDeltaPhiMapsBeforePID(0x0),\r
+ fCosPhiMapsBeforePID(0x0),\r
fDeltaPhiMaps(0x0),\r
+ fDeltaPhiMapsContamination(0x0),\r
fCosPhiMaps(0x0),\r
- fProfileCosPhiMaps(0x0)\r
+ fProfileCosPhiMaps(0x0),\r
+ fDeltaPhiMapsTaggedPhotonic(0x0),\r
+ fCosPhiMapsTaggedPhotonic(0x0),\r
+ fDeltaPhiMapsTaggedNonPhotonic(0x0),\r
+ fCosPhiMapsTaggedNonPhotonic(0x0),\r
+ fDeltaPhiMapsTaggedPhotonicLS(0x0),\r
+ fCosPhiMapsTaggedPhotonicLS(0x0),\r
+ fMCSourceDeltaPhiMaps(0x0),\r
+ fOppSignDeltaPhiMaps(0x0),\r
+ fSameSignDeltaPhiMaps(0x0),\r
+ fOppSignAngle(0x0),\r
+ fSameSignAngle(0x0)\r
{\r
//\r
// named ctor\r
fPID = new AliHFEpid("hfePid");\r
fPIDqa = new AliHFEpidQAmanager;\r
\r
+ fPIDBackground = new AliHFEpid("hfePidBackground");\r
+ fPIDBackgroundqa = new AliHFEpidQAmanager;\r
+\r
+ fPIDTOFOnly = new AliHFEpid("hfePidTOFOnly");\r
+\r
DefineInput(0,TChain::Class());\r
DefineOutput(1, TList::Class());\r
for(Int_t bincless = 0; bincless < fNbBinsCentralityQCumulant; bincless++) {\r
//____________________________________________________________\r
AliAnalysisTaskHFEFlow::AliAnalysisTaskHFEFlow(const AliAnalysisTaskHFEFlow &ref):\r
AliAnalysisTaskSE(ref),\r
- fListHist(0x0),\r
+ fListHist(NULL),\r
fAODAnalysis(ref.fAODAnalysis), \r
+ fUseFlagAOD(ref.fUseFlagAOD),\r
+ fApplyCut(ref.fApplyCut),\r
+ fFlags(ref.fFlags),\r
fVZEROEventPlane(ref.fVZEROEventPlane),\r
fVZEROEventPlaneA(ref.fVZEROEventPlaneA),\r
fVZEROEventPlaneC(ref.fVZEROEventPlaneC),\r
fUseMCReactionPlane(ref.fUseMCReactionPlane),\r
fMCPID(ref.fMCPID),\r
fNoPID(ref.fNoPID),\r
+ fChi2OverNDFCut(ref.fChi2OverNDFCut),\r
+ fMaxdca(ref.fMaxdca),\r
+ fMaxopeningtheta(ref.fMaxopeningtheta),\r
+ fMaxopeningphi(ref.fMaxopeningphi),\r
+ fMaxopening3D(ref.fMaxopening3D),\r
+ fMaxInvmass(ref.fMaxInvmass),\r
+ fSetMassConstraint(ref.fSetMassConstraint),\r
fDebugLevel(ref.fDebugLevel),\r
- fcutsRP(0),\r
- fcutsPOI(0),\r
- fHFECuts(0),\r
- fPID(0),\r
- fPIDqa(0),\r
- fflowEvent(0x0),\r
- fHFEVZEROEventPlane(0x0),\r
- fHistEV(0),\r
- fEventPlane(0x0),\r
- fEventPlaneaftersubtraction(0x0),\r
- fCosSin2phiep(0x0),\r
- fCos2phie(0x0),\r
- fSin2phie(0x0),\r
- fCos2phiep(0x0),\r
- fSin2phiep(0x0),\r
- fSin2phiephiep(0x0),\r
- fCosResabc(0x0),\r
- fSinResabc(0x0),\r
- fProfileCosResab(0x0),\r
- fProfileCosResac(0x0),\r
- fProfileCosResbc(0x0),\r
- fCosRes(0x0),\r
- fSinRes(0x0),\r
- fProfileCosRes(0x0),\r
- fDeltaPhiMaps(0x0),\r
- fCosPhiMaps(0x0),\r
- fProfileCosPhiMaps(0x0)\r
+ fcutsRP(NULL),\r
+ fcutsPOI(NULL),\r
+ fHFECuts(NULL),\r
+ fPID(NULL),\r
+ fPIDTOFOnly(NULL),\r
+ fPIDqa(NULL),\r
+ fflowEvent(NULL),\r
+ fHFEBackgroundCuts(NULL),\r
+ fPIDBackground(NULL),\r
+ fPIDBackgroundqa(NULL),\r
+ fAlgorithmMA(ref.fAlgorithmMA),\r
+ fArraytrack(NULL),\r
+ fCounterPoolBackground(ref.fCounterPoolBackground),\r
+ fHFEVZEROEventPlane(NULL),\r
+ fHistEV(NULL),\r
+ fEventPlane(NULL),\r
+ fEventPlaneaftersubtraction(NULL),\r
+ fCosSin2phiep(NULL),\r
+ fCos2phie(NULL),\r
+ fSin2phie(NULL),\r
+ fCos2phiep(NULL),\r
+ fSin2phiep(NULL),\r
+ fSin2phiephiep(NULL),\r
+ fCosResabc(NULL),\r
+ fSinResabc(NULL),\r
+ fProfileCosResab(NULL),\r
+ fProfileCosResac(NULL),\r
+ fProfileCosResbc(NULL),\r
+ fCosRes(NULL),\r
+ fSinRes(NULL),\r
+ fProfileCosRes(NULL),\r
+ fTrackingCuts(NULL),\r
+ fDeltaPhiMapsBeforePID(NULL),\r
+ fCosPhiMapsBeforePID(NULL),\r
+ fDeltaPhiMaps(NULL),\r
+ fDeltaPhiMapsContamination(NULL),\r
+ fCosPhiMaps(NULL),\r
+ fProfileCosPhiMaps(NULL),\r
+ fDeltaPhiMapsTaggedPhotonic(NULL),\r
+ fCosPhiMapsTaggedPhotonic(NULL),\r
+ fDeltaPhiMapsTaggedNonPhotonic(NULL),\r
+ fCosPhiMapsTaggedNonPhotonic(NULL),\r
+ fDeltaPhiMapsTaggedPhotonicLS(NULL),\r
+ fCosPhiMapsTaggedPhotonicLS(NULL),\r
+ fMCSourceDeltaPhiMaps(NULL),\r
+ fOppSignDeltaPhiMaps(NULL),\r
+ fSameSignDeltaPhiMaps(NULL),\r
+ fOppSignAngle(NULL),\r
+ fSameSignAngle(NULL)\r
{\r
//\r
// Copy Constructor\r
//\r
AliAnalysisTaskHFEFlow &target = dynamic_cast<AliAnalysisTaskHFEFlow &>(o);\r
target.fAODAnalysis = fAODAnalysis;\r
+ target.fUseFlagAOD = fUseFlagAOD;\r
+ target.fApplyCut = fApplyCut;\r
+ target.fFlags = fFlags;\r
target.fVZEROEventPlane = fVZEROEventPlane;\r
target.fVZEROEventPlaneA = fVZEROEventPlaneA;\r
target.fVZEROEventPlaneC = fVZEROEventPlaneC;\r
target.fUseMCReactionPlane = fUseMCReactionPlane;\r
target.fMCPID = fMCPID;\r
target.fNoPID = fNoPID;\r
+ target.fChi2OverNDFCut = fChi2OverNDFCut;\r
+ target.fMaxdca = fMaxdca;\r
+ target.fMaxopeningtheta = fMaxopeningtheta;\r
+ target.fMaxopeningphi = fMaxopeningphi;\r
+ target.fMaxopening3D = fMaxopening3D;\r
+ target.fMaxInvmass = fMaxInvmass;\r
+ target.fSetMassConstraint = fSetMassConstraint;\r
+ target.fAlgorithmMA = fAlgorithmMA;\r
+ target.fCounterPoolBackground = fCounterPoolBackground;\r
target.fDebugLevel = fDebugLevel;\r
target.fcutsRP = fcutsRP;\r
target.fcutsPOI = fcutsPOI;\r
//\r
// Destructor\r
//\r
+ if(fArraytrack) delete fArraytrack;\r
if(fListHist) delete fListHist;\r
if(fcutsRP) delete fcutsRP;\r
if(fcutsPOI) delete fcutsPOI;\r
if(fHFECuts) delete fHFECuts;\r
if(fPID) delete fPID;\r
+ if(fPIDTOFOnly) delete fPIDTOFOnly;\r
if(fPIDqa) delete fPIDqa;\r
if(fflowEvent) delete fflowEvent;\r
- if(fHFEVZEROEventPlane) delete fHFEVZEROEventPlane;\r
- if(fHistEV) delete fHistEV;\r
- if(fEventPlane) delete fEventPlane;\r
- if(fEventPlaneaftersubtraction) delete fEventPlaneaftersubtraction;\r
- if(fCosSin2phiep) delete fCosSin2phiep;\r
- if(fCos2phie) delete fCos2phie;\r
- if(fSin2phie) delete fSin2phie;\r
- if(fCos2phiep) delete fCos2phiep;\r
- if(fSin2phiep) delete fSin2phiep;\r
- if(fSin2phiephiep) delete fSin2phiephiep;\r
- if(fCosResabc) delete fCosResabc;\r
- if(fSinResabc) delete fSinResabc;\r
- if(fProfileCosResab) delete fProfileCosResab;\r
- if(fProfileCosResac) delete fProfileCosResac;\r
- if(fProfileCosResbc) delete fProfileCosResbc;\r
- if(fCosRes) delete fCosRes;\r
- if(fSinRes) delete fSinRes;\r
- if(fProfileCosRes) delete fProfileCosRes;\r
- if(fDeltaPhiMaps) delete fDeltaPhiMaps;\r
- if(fCosPhiMaps) delete fCosPhiMaps;\r
- if(fProfileCosPhiMaps) delete fProfileCosPhiMaps;\r
+ if(fHFEBackgroundCuts) delete fHFEBackgroundCuts;\r
+ if(fPIDBackground) delete fPIDBackground;\r
+ if(fPIDBackgroundqa) delete fPIDBackgroundqa;\r
+ //if(fHFEVZEROEventPlane) delete fHFEVZEROEventPlane;\r
+ \r
\r
}\r
//________________________________________________________________________\r
//kPure - no mixing, kTrackWithMCkine, kTrackWithMCPID, kTrackWithMCpt\r
//AliFlowTrackCuts::trackParameterMix rpmix = AliFlowTrackCuts::kPure;\r
//AliFlowTrackCuts::trackParameterMix poimix = AliFlowTrackCuts::kPure;\r
+ \r
+ // AOD or ESD\r
+ AliVEventHandler *inputHandler = dynamic_cast<AliVEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());\r
+ if(!TString(inputHandler->IsA()->GetName()).CompareTo("AliAODInputHandler")){\r
+ SetAODAnalysis(kTRUE);\r
+ //printf("Put AOD analysis on\n");\r
+ } else {\r
+ SetAODAnalysis(kFALSE);\r
+ }\r
+\r
\r
// RP TRACK CUTS:\r
fcutsRP = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010();\r
fPID->InitializePID();\r
fPIDqa->Initialize(fPID);\r
fPID->SortDetectors();\r
+\r
+ if(!fPIDTOFOnly->GetNumberOfPIDdetectors()) fPIDTOFOnly->AddDetector("TPC", 0);\r
+ fPIDTOFOnly->InitializePID();\r
+ fPIDTOFOnly->SortDetectors();\r
+\r
+ // HFE Background cuts\r
+\r
+ if(!fHFEBackgroundCuts){\r
+ fHFEBackgroundCuts = new AliESDtrackCuts();\r
+ fHFEBackgroundCuts->SetName("nackgroundcuts");\r
+ //Configure Default Track Cuts\r
+ fHFEBackgroundCuts->SetAcceptKinkDaughters(kFALSE);\r
+ fHFEBackgroundCuts->SetRequireTPCRefit(kTRUE);\r
+ fHFEBackgroundCuts->SetEtaRange(-0.9,0.9);\r
+ fHFEBackgroundCuts->SetRequireSigmaToVertex(kTRUE);\r
+ fHFEBackgroundCuts->SetMaxChi2PerClusterTPC(4.0);\r
+ fHFEBackgroundCuts->SetMinNClustersTPC(50);\r
+ fHFEBackgroundCuts->SetPtRange(0.3,1e10);\r
+ }\r
+ \r
+ // PID background HFE\r
+ if(!fPIDBackground->GetNumberOfPIDdetectors()) fPIDBackground->AddDetector("TPC", 0);\r
+ fPIDBackground->InitializePID();\r
+ fPIDBackgroundqa->Initialize(fPIDBackground);\r
+ fPIDBackground->SortDetectors();\r
\r
+\r
+\r
//**************************\r
// Bins for the THnSparse\r
//**************************\r
\r
+ /*\r
Int_t nBinsPt = 44;\r
Double_t minPt = 0.1;\r
Double_t maxPt = 20.0;\r
Double_t binLimPt[nBinsPt+1];\r
for(Int_t i=0; i<=nBinsPt; i++) binLimLogPt[i]=(Double_t)TMath::Log10(minPt) + (TMath::Log10(maxPt)-TMath::Log10(minPt))/nBinsPt*(Double_t)i ;\r
for(Int_t i=0; i<=nBinsPt; i++) binLimPt[i]=(Double_t)TMath::Power(10,binLimLogPt[i]);\r
+ */\r
+\r
+ Int_t nBinsPt = 35;\r
+ Double_t binLimPt[36] = {0., 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1., 1.1, 1.2,\r
+ 1.3, 1.4, 1.5, 1.75, 2., 2.25, 2.5, 2.75, 3., 3.5, 4., 4.5, 5.,\r
+ 5.5, 6., 7., 8., 10., 12., 14., 16., 18., 20.};\r
+\r
\r
Int_t nBinsPtPlus = fNbBinsPtQCumulant;\r
Double_t minPtPlus = fMinPtQCumulant;\r
Double_t binLimEta[nBinsEta+1];\r
for(Int_t i=0; i<=nBinsEta; i++) binLimEta[i]=(Double_t)minEta + (maxEta-minEta)/nBinsEta*(Double_t)i ;\r
\r
+ Int_t nBinsStep = 7;\r
+ Double_t minStep = 0.;\r
+ Double_t maxStep = 7.;\r
+ Double_t binLimStep[nBinsStep+1];\r
+ for(Int_t i=0; i<=nBinsStep; i++) binLimStep[i]=(Double_t)minStep + (maxStep-minStep)/nBinsStep*(Double_t)i ;\r
+\r
Int_t nBinsEtaLess = 2;\r
Double_t binLimEtaLess[nBinsEtaLess+1];\r
for(Int_t i=0; i<=nBinsEtaLess; i++) binLimEtaLess[i]=(Double_t)minEta + (maxEta-minEta)/nBinsEtaLess*(Double_t)i ;\r
Double_t binLimCMore[nBinsCMore+1];\r
for(Int_t i=0; i<=nBinsCMore; i++) binLimCMore[i]=(Double_t)minCMore + (maxCMore-minCMore)/nBinsCMore*(Double_t)i ;\r
\r
- Int_t nBinsPhi = 25;\r
+ Int_t nBinsPhi = 12;\r
Double_t minPhi = 0.0;\r
Double_t maxPhi = TMath::Pi();\r
Double_t binLimPhi[nBinsPhi+1];\r
//printf("bin phi is %f for %d\n",binLimPhi[i],i);\r
}\r
\r
+ Int_t nBinsPhiLess = 2.0;\r
+ Double_t minPhiLess = 0.0;\r
+ Double_t maxPhiLess = 2.0;\r
+ Double_t binLimPhiLess[nBinsPhiLess+1];\r
+ for(Int_t i=0; i<=nBinsPhiLess; i++) {\r
+ binLimPhiLess[i]=(Double_t)minPhiLess + (maxPhiLess-minPhiLess)/nBinsPhiLess*(Double_t)i ;\r
+ }\r
+\r
+ Int_t nBinsTPCdEdx = 140;\r
+ Double_t minTPCdEdx = -12.0;\r
+ Double_t maxTPCdEdx = 12.0;\r
+ Double_t binLimTPCdEdx[nBinsTPCdEdx+1];\r
+ for(Int_t i=0; i<=nBinsTPCdEdx; i++) {\r
+ binLimTPCdEdx[i]=(Double_t)minTPCdEdx + (maxTPCdEdx-minTPCdEdx)/nBinsTPCdEdx*(Double_t)i ;\r
+ }\r
+\r
+ Int_t nBinsAngle = 40;\r
+ Double_t minAngle = 0.0;\r
+ Double_t maxAngle = 1.0;\r
+ Double_t binLimAngle[nBinsAngle+1];\r
+ for(Int_t i=0; i<=nBinsAngle; i++) {\r
+ binLimAngle[i]=(Double_t)minAngle + (maxAngle-minAngle)/nBinsAngle*(Double_t)i ;\r
+ //printf("bin phi is %f for %d\n",binLimPhi[i],i);\r
+ }\r
+\r
Int_t nBinsCharge = 2;\r
Double_t minCharge = -1.0;\r
Double_t maxCharge = 1.0;\r
Double_t binLimCharge[nBinsCharge+1];\r
for(Int_t i=0; i<=nBinsCharge; i++) binLimCharge[i]=(Double_t)minCharge + (maxCharge-minCharge)/nBinsCharge*(Double_t)i ;\r
+\r
+ Int_t nBinsSource = 10;\r
+ Double_t minSource = 0.;\r
+ Double_t maxSource = 10.;\r
+ Double_t binLimSource[nBinsSource+1];\r
+ for(Int_t i=0; i<=nBinsSource; i++) binLimSource[i]=(Double_t)minSource + (maxSource-minSource)/nBinsSource*(Double_t)i ;\r
+\r
+ Int_t nBinsInvMass = 50;\r
+ Double_t minInvMass = 0.;\r
+ Double_t maxInvMass = 0.3;\r
+ Double_t binLimInvMass[nBinsInvMass+1];\r
+ for(Int_t i=0; i<=nBinsInvMass; i++) binLimInvMass[i]=(Double_t)minInvMass + (maxInvMass-minInvMass)/nBinsInvMass*(Double_t)i ;\r
\r
//******************\r
// Histograms\r
//******************\r
\r
fListHist = new TList();\r
+ fListHist->SetOwner();\r
\r
// Histos\r
fHistEV = new TH2D("fHistEV", "events", 3, 0, 3, 3, 0,3);\r
// Profile cosres centrality\r
fProfileCosRes = new TProfile("ProfileCosRes","ProfileCosRes",nBinsCMore,binLimCMore);\r
fProfileCosRes->Sumw2();\r
+\r
+ // Debugging tracking steps\r
+ const Int_t nDimTrStep=2;\r
+ Int_t nBinTrStep[nDimTrStep] = {nBinsPt,nBinsStep};\r
+ fTrackingCuts = new THnSparseF("TrackingCuts","TrackingCuts",nDimTrStep,nBinTrStep);\r
+ fTrackingCuts->SetBinEdges(0,binLimPt);\r
+ fTrackingCuts->SetBinEdges(1,binLimStep);\r
+ fTrackingCuts->Sumw2();\r
\r
// Maps delta phi\r
const Int_t nDimg=5;\r
fDeltaPhiMaps->SetBinEdges(4,binLimEtaLess);\r
fDeltaPhiMaps->Sumw2(); \r
\r
+ // Maps delta phi contamination\r
+ const Int_t nDimgcont=4;\r
+ Int_t nBingcont[nDimgcont] = {nBinsPhiLess,nBinsC,nBinsPt, nBinsTPCdEdx};\r
+ fDeltaPhiMapsContamination = new THnSparseF("DeltaPhiMapsContamination","DeltaPhiMapsContamination",nDimgcont,nBingcont);\r
+ fDeltaPhiMapsContamination->SetBinEdges(0,binLimPhiLess);\r
+ fDeltaPhiMapsContamination->SetBinEdges(1,binLimC);\r
+ fDeltaPhiMapsContamination->SetBinEdges(2,binLimPt);\r
+ fDeltaPhiMapsContamination->SetBinEdges(3,binLimTPCdEdx);\r
+ fDeltaPhiMapsContamination->Sumw2(); \r
+\r
+\r
+ //\r
+ const Int_t nDimgb=3;\r
+ Int_t nBingb[nDimgb] = {nBinsPhi,nBinsC,nBinsPt};\r
+\r
+ fDeltaPhiMapsBeforePID = new THnSparseF("DeltaPhiMapsBeforePID","DeltaPhiMapsBeforePID",nDimgb,nBingb);\r
+ fDeltaPhiMapsBeforePID->SetBinEdges(0,binLimPhi);\r
+ fDeltaPhiMapsBeforePID->SetBinEdges(1,binLimC);\r
+ fDeltaPhiMapsBeforePID->SetBinEdges(2,binLimPt);\r
+ fDeltaPhiMapsBeforePID->Sumw2(); \r
+\r
+ \r
+ fDeltaPhiMapsTaggedPhotonic = new THnSparseF("DeltaPhiMapsTaggedPhotonic","DeltaPhiMapsTaggedPhotonic",nDimgb,nBingb);\r
+ fDeltaPhiMapsTaggedPhotonic->SetBinEdges(0,binLimPhi);\r
+ fDeltaPhiMapsTaggedPhotonic->SetBinEdges(1,binLimC);\r
+ fDeltaPhiMapsTaggedPhotonic->SetBinEdges(2,binLimPt);\r
+ fDeltaPhiMapsTaggedPhotonic->Sumw2(); \r
+\r
+ fDeltaPhiMapsTaggedNonPhotonic = new THnSparseF("DeltaPhiMapsTaggedNonPhotonic","DeltaPhiMapsTaggedNonPhotonic",nDimgb,nBingb);\r
+ fDeltaPhiMapsTaggedNonPhotonic->SetBinEdges(0,binLimPhi);\r
+ fDeltaPhiMapsTaggedNonPhotonic->SetBinEdges(1,binLimC);\r
+ fDeltaPhiMapsTaggedNonPhotonic->SetBinEdges(2,binLimPt);\r
+ fDeltaPhiMapsTaggedNonPhotonic->Sumw2(); \r
+\r
+ fDeltaPhiMapsTaggedPhotonicLS = new THnSparseF("DeltaPhiMapsTaggedPhotonicLS","DeltaPhiMapsTaggedPhotonicLS",nDimgb,nBingb);\r
+ fDeltaPhiMapsTaggedPhotonicLS->SetBinEdges(0,binLimPhi);\r
+ fDeltaPhiMapsTaggedPhotonicLS->SetBinEdges(1,binLimC);\r
+ fDeltaPhiMapsTaggedPhotonicLS->SetBinEdges(2,binLimPt);\r
+ fDeltaPhiMapsTaggedPhotonicLS->Sumw2(); \r
+\r
// Maps cos phi\r
const Int_t nDimh=5;\r
Int_t nBinh[nDimh] = {nBinsCos,nBinsC,nBinsPt,nBinsCharge,nBinsEtaLess};\r
fCosPhiMaps->SetBinEdges(4,binLimEtaLess);\r
fCosPhiMaps->Sumw2();\r
\r
+ const Int_t nDimhb=3;\r
+ Int_t nBinhb[nDimhb] = {nBinsCos,nBinsC,nBinsPt};\r
+\r
+ fCosPhiMapsBeforePID = new THnSparseF("CosPhiMapsBeforePID","CosPhiMapsBeforePID",nDimhb,nBinhb);\r
+ fCosPhiMapsBeforePID->SetBinEdges(0,binLimCos);\r
+ fCosPhiMapsBeforePID->SetBinEdges(1,binLimC);\r
+ fCosPhiMapsBeforePID->SetBinEdges(2,binLimPt);\r
+ fCosPhiMapsBeforePID->Sumw2();\r
+\r
+ fCosPhiMapsTaggedPhotonic = new THnSparseF("CosPhiMapsTaggedPhotonic","CosPhiMapsTaggedPhotonic",nDimhb,nBinhb);\r
+ fCosPhiMapsTaggedPhotonic->SetBinEdges(0,binLimCos);\r
+ fCosPhiMapsTaggedPhotonic->SetBinEdges(1,binLimC);\r
+ fCosPhiMapsTaggedPhotonic->SetBinEdges(2,binLimPt);\r
+ fCosPhiMapsTaggedPhotonic->Sumw2();\r
+\r
+ fCosPhiMapsTaggedNonPhotonic = new THnSparseF("CosPhiMapsTaggedNonPhotonic","CosPhiMapsTaggedNonPhotonic",nDimhb,nBinhb);\r
+ fCosPhiMapsTaggedNonPhotonic->SetBinEdges(0,binLimCos);\r
+ fCosPhiMapsTaggedNonPhotonic->SetBinEdges(1,binLimC);\r
+ fCosPhiMapsTaggedNonPhotonic->SetBinEdges(2,binLimPt);\r
+ fCosPhiMapsTaggedNonPhotonic->Sumw2();\r
+ \r
+ fCosPhiMapsTaggedPhotonicLS = new THnSparseF("CosPhiMapsTaggedPhotonicLS","CosPhiMapsTaggedPhotonicLS",nDimhb,nBinhb);\r
+ fCosPhiMapsTaggedPhotonicLS->SetBinEdges(0,binLimCos);\r
+ fCosPhiMapsTaggedPhotonicLS->SetBinEdges(1,binLimC);\r
+ fCosPhiMapsTaggedPhotonicLS->SetBinEdges(2,binLimPt);\r
+ fCosPhiMapsTaggedPhotonicLS->Sumw2();\r
+\r
// Profile Maps cos phi\r
fProfileCosPhiMaps = new TProfile2D("ProfileCosPhiMaps","ProfileCosPhiMaps",nBinsC,binLimC,nBinsPt,binLimPt);\r
fProfileCosPhiMaps->Sumw2();\r
\r
+ // Background study\r
+ const Int_t nDimMCSource=3;\r
+ Int_t nBinMCSource[nDimMCSource] = {nBinsC,nBinsPt,nBinsSource};\r
+ fMCSourceDeltaPhiMaps = new THnSparseF("MCSourceDeltaPhiMaps","MCSourceDeltaPhiMaps",nDimMCSource,nBinMCSource);\r
+ fMCSourceDeltaPhiMaps->SetBinEdges(0,binLimC);\r
+ fMCSourceDeltaPhiMaps->SetBinEdges(1,binLimPt);\r
+ fMCSourceDeltaPhiMaps->SetBinEdges(2,binLimSource);\r
+ fMCSourceDeltaPhiMaps->Sumw2();\r
+\r
+ // Maps invmass opposite\r
+ const Int_t nDimOppSign=5;\r
+ Int_t nBinOppSign[nDimOppSign] = {nBinsPhi,nBinsC,nBinsPt,nBinsInvMass,nBinsSource};\r
+ fOppSignDeltaPhiMaps = new THnSparseF("OppSignDeltaPhiMaps","OppSignDeltaPhiMaps",nDimOppSign,nBinOppSign);\r
+ fOppSignDeltaPhiMaps->SetBinEdges(0,binLimPhi);\r
+ fOppSignDeltaPhiMaps->SetBinEdges(1,binLimC);\r
+ fOppSignDeltaPhiMaps->SetBinEdges(2,binLimPt);\r
+ fOppSignDeltaPhiMaps->SetBinEdges(3,binLimInvMass);\r
+ fOppSignDeltaPhiMaps->SetBinEdges(4,binLimSource);\r
+ fOppSignDeltaPhiMaps->Sumw2();\r
+\r
+ // Maps invmass same sign\r
+ const Int_t nDimSameSign=5;\r
+ Int_t nBinSameSign[nDimSameSign] = {nBinsPhi,nBinsC,nBinsPt,nBinsInvMass,nBinsSource};\r
+ fSameSignDeltaPhiMaps = new THnSparseF("SameSignDeltaPhiMaps","SameSignDeltaPhiMaps",nDimSameSign,nBinSameSign);\r
+ fSameSignDeltaPhiMaps->SetBinEdges(0,binLimPhi);\r
+ fSameSignDeltaPhiMaps->SetBinEdges(1,binLimC);\r
+ fSameSignDeltaPhiMaps->SetBinEdges(2,binLimPt);\r
+ fSameSignDeltaPhiMaps->SetBinEdges(3,binLimInvMass);\r
+ fSameSignDeltaPhiMaps->SetBinEdges(4,binLimSource);\r
+ fSameSignDeltaPhiMaps->Sumw2();\r
+\r
+ // Maps angle same sign\r
+ const Int_t nDimAngleSameSign=3;\r
+ Int_t nBinAngleSameSign[nDimAngleSameSign] = {nBinsAngle,nBinsC,nBinsSource};\r
+ fSameSignAngle = new THnSparseF("SameSignAngleMaps","SameSignAngleMaps",nDimAngleSameSign,nBinAngleSameSign);\r
+ fSameSignAngle->SetBinEdges(0,binLimAngle);\r
+ fSameSignAngle->SetBinEdges(1,binLimC);\r
+ fSameSignAngle->SetBinEdges(2,binLimSource);\r
+ fSameSignAngle->Sumw2();\r
+\r
+ // Maps angle opp sign\r
+ const Int_t nDimAngleOppSign=3;\r
+ Int_t nBinAngleOppSign[nDimAngleOppSign] = {nBinsAngle,nBinsC,nBinsSource};\r
+ fOppSignAngle = new THnSparseF("OppSignAngleMaps","OppSignAngleMaps",nDimAngleOppSign,nBinAngleOppSign);\r
+ fOppSignAngle->SetBinEdges(0,binLimAngle);\r
+ fOppSignAngle->SetBinEdges(1,binLimC);\r
+ fOppSignAngle->SetBinEdges(2,binLimSource);\r
+ fOppSignAngle->Sumw2();\r
+\r
\r
//**************************\r
// Add to the list\r
\r
//fListHist->Add(qaCutsRP);\r
fListHist->Add(fPIDqa->MakeList("HFEpidQA"));\r
+ fListHist->Add(fPIDBackgroundqa->MakeList("HFEpidBackgroundQA"));\r
fListHist->Add(fHistEV);\r
fListHist->Add(fProfileCosRes);\r
fListHist->Add(fProfileCosResab);\r
fListHist->Add(fCosResabc);\r
fListHist->Add(fSinRes);\r
fListHist->Add(fSinResabc);\r
+ fListHist->Add(fTrackingCuts);\r
+ fListHist->Add(fDeltaPhiMapsBeforePID);\r
+ fListHist->Add(fCosPhiMapsBeforePID);\r
fListHist->Add(fDeltaPhiMaps);\r
+ fListHist->Add(fDeltaPhiMapsContamination);\r
fListHist->Add(fCosPhiMaps);\r
fListHist->Add(fProfileCosPhiMaps);\r
+ fListHist->Add(fDeltaPhiMapsTaggedPhotonic);\r
+ fListHist->Add(fCosPhiMapsTaggedPhotonic);\r
+ fListHist->Add(fDeltaPhiMapsTaggedNonPhotonic);\r
+ fListHist->Add(fCosPhiMapsTaggedNonPhotonic);\r
+ fListHist->Add(fDeltaPhiMapsTaggedPhotonicLS);\r
+ fListHist->Add(fCosPhiMapsTaggedPhotonicLS);\r
+ fListHist->Add(fMCSourceDeltaPhiMaps);\r
+ fListHist->Add(fOppSignDeltaPhiMaps);\r
+ fListHist->Add(fSameSignDeltaPhiMaps);\r
+ fListHist->Add(fSameSignAngle);\r
+ fListHist->Add(fOppSignAngle);\r
+\r
\r
if(fHFEVZEROEventPlane && (fDebugLevel > 2)) fListHist->Add(fHFEVZEROEventPlane->GetOutputList());\r
\r
Double_t valuensparseg[5];\r
Double_t valuensparseh[5];\r
Double_t valuensparsehprofile[3];\r
-\r
+ Double_t valuensparseMCSourceDeltaPhiMaps[3];\r
+ Double_t valuetrackingcuts[2];\r
+ Double_t valuedeltaphicontamination[4];\r
+ \r
AliMCEvent *mcEvent = MCEvent();\r
AliMCParticle *mctrack = NULL;\r
- \r
+ \r
/////////////////\r
// centrality\r
/////////////////\r
//AliESDEvent *esd = dynamic_cast<AliESDEvent*>(InputEvent());\r
//if(!esd) return;\r
AliCentrality *centrality = fInputEvent->GetCentrality();\r
+ //printf("Got the centrality\n");\r
if(!centrality) return;\r
cntr = centrality->GetCentralityPercentile("V0M");\r
if((0.0< cntr) && (cntr<5.0)) binct = 0.5;\r
valuensparseh[1] = binct; \r
valuensparsehprofile[1] = binct; \r
valuecossinephiep[2] = binctMore;\r
-\r
+ valuensparseMCSourceDeltaPhiMaps[0] = binct;\r
+ valuedeltaphicontamination[1] = binct;\r
+ \r
//////////////////////\r
// run number\r
//////////////////////\r
\r
Int_t runnumber = fInputEvent->GetRunNumber();\r
+ //printf("Run number %d\n",runnumber);\r
\r
if(!fPID->IsInitialized()){\r
// Initialize PID with the given run number\r
fPID->InitializePID(runnumber);\r
}\r
+ if(!fPIDTOFOnly->IsInitialized()){\r
+ // Initialize PID with the given run number\r
+ fPIDTOFOnly->InitializePID(runnumber);\r
+ }\r
+\r
+ //\r
+ if(!fPIDBackground->IsInitialized()){\r
+ // Initialize PID with the given run number\r
+ fPIDBackground->InitializePID(runnumber);\r
+ }\r
+\r
+ fHFECuts->SetRecEvent(fInputEvent);\r
+ if(mcEvent) fHFECuts->SetMCEvent(mcEvent);\r
\r
\r
//////////\r
return;\r
}\r
fPID->SetPIDResponse(pidResponse);\r
+ fPIDTOFOnly->SetPIDResponse(pidResponse);\r
+ fPIDBackground->SetPIDResponse(pidResponse);\r
\r
fHistEV->Fill(binctt,0.0);\r
\r
// Event cut\r
//////////////////\r
if(!fHFECuts->CheckEventCuts("fEvRecCuts", fInputEvent)) {\r
+ //printf("Do not pass the event cut\n");\r
PostData(1, fListHist);\r
return;\r
}\r
else {\r
\r
eventPlaneV0 = TVector2::Phi_0_2pi(vEPa->GetEventplane("V0", fInputEvent,2));\r
+ //printf("eventPlaneV0 %f\n",eventPlaneV0);\r
if(eventPlaneV0 > TMath::Pi()) eventPlaneV0 = eventPlaneV0 - TMath::Pi();\r
+ //printf("eventPlaneV0 %f\n",eventPlaneV0);\r
eventPlaneV0A = TVector2::Phi_0_2pi(vEPa->GetEventplane("V0A", fInputEvent,2));\r
if(eventPlaneV0A > TMath::Pi()) eventPlaneV0A = eventPlaneV0A - TMath::Pi();\r
eventPlaneV0C = TVector2::Phi_0_2pi(vEPa->GetEventplane("V0C", fInputEvent,2));\r
//}\r
//}\r
\r
- if((!standardQ) || (!qsub1a) || (!qsub2a)) return;\r
+ if((!standardQ) || (!qsub1a) || (!qsub2a)) {\r
+ //printf("No event plane\n");\r
+ return;\r
+ }\r
\r
///////////////////////\r
// AliFlowEvent \r
fEventPlane->Fill(&valuensparsea[0]);\r
\r
// Fill\r
- if(fDebugLevel > 3) fCosSin2phiep->Fill(&valuecossinephiep[0]);\r
+ if(fDebugLevel > 5) fCosSin2phiep->Fill(&valuecossinephiep[0]);\r
\r
if(!fVZEROEventPlane) {\r
valuensparsef[0] = diffsub1sub2a;\r
fCosRes->Fill(&valuensparsef[0]);\r
valuensparsefsin[0] = diffsub1sub2asin;\r
- fSinRes->Fill(&valuensparsefsin[0]);\r
- if(fDebugLevel > 1) {\r
+ if(fDebugLevel > 5) fSinRes->Fill(&valuensparsefsin[0]);\r
+ if(fDebugLevel > 5) {\r
fProfileCosRes->Fill(valuensparsef[1],valuensparsef[0]);\r
}\r
}\r
valuensparsefbissin[0] = diffsubasubbsin;\r
valuensparsefbissin[1] = diffsubbsubcsin;\r
valuensparsefbissin[2] = diffsubasubcsin;\r
- fSinResabc->Fill(&valuensparsefbissin[0]);\r
- if(fDebugLevel > 1) {\r
+ if(fDebugLevel > 5) fSinResabc->Fill(&valuensparsefbissin[0]);\r
+ if(fDebugLevel > 5) {\r
fProfileCosResab->Fill(valuensparsefbis[3],valuensparsefbis[0]);\r
fProfileCosResac->Fill(valuensparsefbis[3],valuensparsefbis[1]);\r
fProfileCosResbc->Fill(valuensparsefbis[3],valuensparsefbis[2]);\r
\r
}\r
\r
- //////////////////////////\r
- // Loop over ESD track\r
- //////////////////////////\r
- \r
-\r
+ ////////////////////////////////////////\r
+ // Loop to determine pool background\r
+ /////////////////////////////////////////\r
+ if( fArraytrack ){ \r
+ fArraytrack->~TArrayI();\r
+ new(fArraytrack) TArrayI(nbtracks);\r
+ }\r
+ else { \r
+ fArraytrack = new TArrayI(nbtracks);\r
+ }\r
+ fCounterPoolBackground = 0;\r
for(Int_t k = 0; k < nbtracks; k++){\r
\r
AliVTrack *track = (AliVTrack *) fInputEvent->GetTrack(k);\r
if(!track) continue;\r
+ \r
+ // Track cuts\r
+ Bool_t survivedbackground = kTRUE;\r
+ if(fAODAnalysis) {\r
+ AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);\r
+ if(aodtrack) {\r
+ AliESDtrack esdTrack(aodtrack);\r
+ // set the TPC cluster info\r
+ esdTrack.SetTPCClusterMap(aodtrack->GetTPCClusterMap());\r
+ esdTrack.SetTPCSharedMap(aodtrack->GetTPCSharedMap());\r
+ esdTrack.SetTPCPointsF(aodtrack->GetTPCNclsF());\r
+ // needed to calculate the impact parameters\r
+ AliAODEvent *aodeventu = dynamic_cast<AliAODEvent *>(fInputEvent);\r
+ if(aodeventu) {\r
+ AliAODVertex *vAOD = aodeventu->GetPrimaryVertex();\r
+ Double_t bfield = aodeventu->GetMagneticField();\r
+ Double_t pos[3],cov[6];\r
+ vAOD->GetXYZ(pos);\r
+ vAOD->GetCovarianceMatrix(cov);\r
+ const AliESDVertex vESD(pos,cov,100.,100);\r
+ esdTrack.RelateToVertex(&vESD,bfield,3.);\r
+ } \r
+ if(!fHFEBackgroundCuts->IsSelected(&esdTrack)) {\r
+ survivedbackground = kFALSE;\r
+ }\r
+ }\r
+ }\r
+ else {\r
+ AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);\r
+ if(esdtrack) {\r
+ if(!fHFEBackgroundCuts->IsSelected(esdtrack)) survivedbackground = kFALSE;\r
+ }\r
+ }\r
+ // PID\r
+ if(survivedbackground) {\r
+ // PID track cuts\r
+ AliHFEpidObject hfetrack2;\r
+ if(!fAODAnalysis) hfetrack2.SetAnalysisType(AliHFEpidObject::kESDanalysis);\r
+ else hfetrack2.SetAnalysisType(AliHFEpidObject::kAODanalysis);\r
+ hfetrack2.SetRecTrack(track);\r
+ hfetrack2.SetCentrality((Int_t)binct);\r
+ //printf("centrality %f and %d\n",binct,hfetrack.GetCentrality());\r
+ hfetrack2.SetPbPb();\r
+ if(fPIDBackground->IsSelected(&hfetrack2,0x0,"recTrackCont",fPIDBackgroundqa)) {\r
+ fArraytrack->AddAt(k,fCounterPoolBackground);\r
+ fCounterPoolBackground++;\r
+ //printf("fCounterPoolBackground %d, track %d\n",fCounterPoolBackground,k);\r
+ }\r
+ }\r
+ }\r
\r
- Bool_t survived = kTRUE;\r
- for(Int_t icut = AliHFEcuts::kStepRecKineITSTPC; icut <= AliHFEcuts::kStepHFEcutsTRD; icut++){\r
- if(!fHFECuts->CheckParticleCuts(icut + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)){\r
- survived = kFALSE;\r
- break;\r
+ // Look at kink mother in case of AOD\r
+ Int_t numberofvertices = 1;\r
+ AliAODEvent *aodevent = NULL;\r
+ Int_t numberofmotherkink = 0; \r
+ if(fAODAnalysis) {\r
+ aodevent = dynamic_cast<AliAODEvent *>(fInputEvent);\r
+ if(aodevent) {\r
+ numberofvertices = aodevent->GetNumberOfVertices();\r
+ }\r
+ }\r
+ Double_t listofmotherkink[numberofvertices];\r
+ if(fAODAnalysis && aodevent) {\r
+ for(Int_t ivertex=0; ivertex < numberofvertices; ivertex++) {\r
+ AliAODVertex *aodvertex = aodevent->GetVertex(ivertex);\r
+ if(!aodvertex) continue;\r
+ if(aodvertex->GetType()==AliAODVertex::kKink) {\r
+ AliAODTrack *mother = (AliAODTrack *) aodvertex->GetParent();\r
+ if(!mother) continue;\r
+ Int_t idmother = mother->GetID();\r
+ listofmotherkink[numberofmotherkink] = idmother;\r
+ //printf("ID %d\n",idmother);\r
+ numberofmotherkink++;\r
}\r
- //printf("Pass the cut %d\n",icut);\r
}\r
+ }\r
+ \r
+ //////////////////////////\r
+ // Loop over track\r
+ //////////////////////////\r
+ \r
+ for(Int_t k = 0; k < nbtracks; k++){\r
+ \r
+ AliVTrack *track = (AliVTrack *) fInputEvent->GetTrack(k);\r
+ if(!track) continue;\r
\r
- if(!survived) continue;\r
-\r
- // Apply PID\r
- if(!fNoPID) {\r
- // Apply PID for Data\r
- if(!fMCPID) {\r
- AliHFEpidObject hfetrack;\r
- if(!fAODAnalysis) hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);\r
- else hfetrack.SetAnalysisType(AliHFEpidObject::kAODanalysis);\r
- hfetrack.SetRecTrack(track);\r
- hfetrack.SetCentrality((Int_t)binct);\r
- //printf("centrality %f and %d\n",binct,hfetrack.GetCentrality());\r
- hfetrack.SetPbPb();\r
- if(!fPID->IsSelected(&hfetrack,0x0,"recTrackCont",fPIDqa)) {\r
- continue;\r
+ if(fAODAnalysis) {\r
+ AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);\r
+ if(!aodtrack){\r
+ AliError("AOD track is not there");\r
+ continue;\r
+ } \r
+ //printf("Find AOD track on\n");\r
+ if(fUseFlagAOD){\r
+ if(aodtrack->GetFlags() != fFlags) continue; // Only process AOD tracks where the HFE is set\r
+ }\r
+ }\r
+ \r
+ if(fApplyCut) {\r
+\r
+ valuetrackingcuts[0] = track->Pt(); \r
+ valuetrackingcuts[1] = 0;\r
+\r
+ // RecKine: ITSTPC cuts \r
+ if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepRecKineITSTPC + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
+ if(fDebugLevel > 3) fTrackingCuts->Fill(&valuetrackingcuts[0]);\r
+ \r
+ // Reject kink mother\r
+ if(fAODAnalysis) {\r
+ Bool_t kinkmotherpass = kTRUE;\r
+ for(Int_t kinkmother = 0; kinkmother < numberofmotherkink; kinkmother++) {\r
+ if(track->GetID() == listofmotherkink[kinkmother]) {\r
+ kinkmotherpass = kFALSE;\r
+ continue;\r
+ }\r
}\r
+ if(!kinkmotherpass) continue;\r
}\r
else {\r
- if(!mcEvent) continue;\r
- if(!(mctrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(track->GetLabel()))))) continue;\r
- //printf("PdgCode %d\n",TMath::Abs(mctrack->Particle()->GetPdgCode()));\r
- if(TMath::Abs(mctrack->Particle()->GetPdgCode())!=11) continue;\r
+ AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);\r
+ if(esdtrack){ \r
+ if(esdtrack->GetKinkIndex(0) != 0) continue; \r
+ } // Quick and dirty fix to reject both kink mothers and daughters\r
}\r
+ \r
+ valuetrackingcuts[1] = 1; \r
+ if(fDebugLevel > 3) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
+ // RecPrim\r
+ if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepRecPrim + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
+ valuetrackingcuts[1] = 2; \r
+ if(fDebugLevel > 3) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
+\r
+ // HFEcuts: ITS layers cuts\r
+ if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsITS + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
+ valuetrackingcuts[1] = 3; \r
+ if(fDebugLevel > 3) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
+\r
+ // HFE cuts: TOF PID and mismatch flag\r
+ if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTOF + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
+ valuetrackingcuts[1] = 4; \r
+ if(fDebugLevel > 3) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
+ \r
+ // HFE cuts: TPC PID cleanup\r
+ if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTPC + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
+ valuetrackingcuts[1] = 5; \r
+ if(fDebugLevel > 3) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
+ \r
+ // HFEcuts: Nb of tracklets TRD0\r
+ if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
+ valuetrackingcuts[1] = 6; \r
+ if(fDebugLevel > 3) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
+ \r
}\r
\r
+ //printf("Survived\n");\r
\r
- /////////////////////////////////////////////////////////////////////////////\r
- // Add candidate to AliFlowEvent for POI and subtract from RP if needed\r
- ////////////////////////////////////////////////////////////////////////////\r
- Int_t idtrack = static_cast<AliVTrack*>(track)->GetID();\r
- Bool_t found = kFALSE;\r
- Int_t numberoffound = 0;\r
- //printf("A: Number of tracks %d\n",fflowEvent->NumberOfTracks());\r
- for(Int_t iRPs=0; iRPs< fflowEvent->NumberOfTracks(); iRPs++) {\r
- AliFlowTrack *iRP = (AliFlowTrack*) (fflowEvent->GetTrack(iRPs));\r
- //if(!iRP->InRPSelection()) continue;\r
- if( TMath::Abs(idtrack) == TMath::Abs(iRP->GetID()) ) {\r
- iRP->SetForPOISelection(kTRUE);\r
- found = kTRUE;\r
- numberoffound ++;\r
- }\r
- }\r
- //printf("Found %d mal\n",numberoffound);\r
- if(!found) {\r
- AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*) MakeTrack(massElectron,track->Pt(),track->Phi(), track->Eta());\r
- sTrack->SetID(idtrack);\r
- fflowEvent->AddTrack(sTrack);\r
- //printf("Add the track\n");\r
- }\r
- //printf("B: Number of tracks %d\n",fflowEvent->NumberOfTracks());\r
- \r
- \r
/////////////////////////////////////////////////////////\r
- // Subtract electron candidate from TPC event plane\r
+ // Subtract candidate from TPC event plane\r
////////////////////////////////////////////////////////\r
Float_t eventplanesubtracted = 0.0; \r
\r
}\r
else eventplanesubtracted = eventPlanea;\r
\r
- ////////////////////////////////////////\r
- // Fill pt and eta for the THnSparseF\r
- ///////////////////////////////////////\r
-\r
- valuensparsee[2] = track->Pt();\r
- valuensparsee[3] = track->Eta(); \r
- valuensparseg[2] = track->Pt();\r
- valuensparseh[2] = track->Pt();\r
- valuensparsehprofile[2] = track->Pt();\r
- if(track->Charge() > 0.0) {\r
- valuensparseg[3] = 0.2;\r
- valuensparseh[3] = 0.2;\r
- }\r
- else {\r
- valuensparseg[3] = -0.2;\r
- valuensparseh[3] = -0.2;\r
- }\r
- valuensparseh[4] = track->Eta();\r
- valuensparseg[4] = track->Eta();\r
-\r
- //printf("charge %d\n",(Int_t)track->Charge());\r
-\r
+ ///////////////////////////////////////////\r
+ // Event plane\r
+ //////////////////////////////////////////\r
Bool_t fillEventPlane = kTRUE;\r
if(!fVZEROEventPlane){\r
//if((!qsub1a) || (!qsub2a) || (!eventplanedefined)) fillEventPlane = kFALSE;\r
}\r
}\r
\r
- \r
///////////////\r
// MC\r
//////////////\r
// Suppose phi track is between 0.0 and phi\r
Double_t deltaphi = TVector2::Phi_0_2pi(phitrack - eventplanesubtracted);\r
if(deltaphi > TMath::Pi()) deltaphi = deltaphi - TMath::Pi();\r
- \r
+\r
+ ////////////////////////////////\r
+ // Determine the deltaphi bin\r
+ ///////////////////////////////\r
+\r
+ // in-plane\r
+ if(((deltaphi<(TMath::Pi()/4.)) && (deltaphi>0.0)) || ((deltaphi>(3*TMath::Pi()/4.)) && (deltaphi<TMath::Pi()))) valuedeltaphicontamination[0] = 0.5;\r
+ // out-of-plane\r
+ if((deltaphi>(TMath::Pi()/4.)) && (deltaphi<(3*TMath::Pi()/4.))) valuedeltaphicontamination[0] = 1.5;\r
+\r
+ ////////////////////////////////////////\r
+ // Define variables\r
+ ///////////////////////////////////////\r
+\r
+\r
+ valuedeltaphicontamination[2] = track->Pt();\r
+ valuensparsee[2] = track->Pt();\r
+ valuensparsee[3] = track->Eta(); \r
+ valuensparseg[2] = track->Pt();\r
+ valuensparseh[2] = track->Pt();\r
+ valuensparsehprofile[2] = track->Pt();\r
+ valuensparseMCSourceDeltaPhiMaps[1] = track->Pt();\r
+ if(track->Charge() > 0.0) {\r
+ valuensparseg[3] = 0.2;\r
+ valuensparseh[3] = 0.2;\r
+ }\r
+ else {\r
+ valuensparseg[3] = -0.2;\r
+ valuensparseh[3] = -0.2;\r
+ }\r
+ valuensparseh[4] = track->Eta();\r
+ valuensparseg[4] = track->Eta();\r
+\r
+ //printf("charge %d\n",(Int_t)track->Charge());\r
+\r
+ ////////////////////////\r
+ // Fill before PID\r
+ ///////////////////////\r
+ \r
+ if(fDebugLevel > 5) { \r
+ \r
+ valuensparseg[0] = deltaphi;\r
+ if(fillEventPlane) fDeltaPhiMapsBeforePID->Fill(&valuensparseg[0]);\r
+ \r
+ //\r
+ valuensparseh[0] = TMath::Cos(2*TVector2::Phi_mpi_pi(phitrack-eventplanesubtracted));\r
+ if(fillEventPlane) {\r
+ fCosPhiMapsBeforePID->Fill(&valuensparseh[0]);\r
+ }\r
+ }\r
+ \r
+ ////////////////////////\r
+ // Apply PID\r
+ ////////////////////////\r
+ if(!fNoPID) {\r
+ // Apply PID for Data\r
+ if(!fMCPID) {\r
+ // pid object\r
+ AliHFEpidObject hfetrack;\r
+ if(!fAODAnalysis) hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);\r
+ else hfetrack.SetAnalysisType(AliHFEpidObject::kAODanalysis);\r
+ hfetrack.SetRecTrack(track);\r
+ hfetrack.SetCentrality((Int_t)binct);\r
+ //printf("centrality %f and %d\n",binct,hfetrack.GetCentrality());\r
+ hfetrack.SetPbPb();\r
+\r
+ // Only TOF PID\r
+ if(fPIDTOFOnly->IsSelected(&hfetrack,0x0,"recTrackCont",0x0)) {\r
+ Float_t nsigma = pidResponse->NumberOfSigmasTPC(track, AliPID::kElectron);\r
+ valuedeltaphicontamination[3] = nsigma;\r
+ fDeltaPhiMapsContamination->Fill(&valuedeltaphicontamination[0]);\r
+ }\r
+\r
+ // Complete PID TOF+TPC\r
+ if(!fPID->IsSelected(&hfetrack,0x0,"recTrackCont",fPIDqa)) {\r
+ continue;\r
+ }\r
+ }\r
+ else {\r
+ if(!mcEvent) continue;\r
+ if(!(mctrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(track->GetLabel()))))) continue;\r
+ //printf("PdgCode %d\n",TMath::Abs(mctrack->Particle()->GetPdgCode()));\r
+ if(TMath::Abs(mctrack->Particle()->GetPdgCode())!=11) continue;\r
+ }\r
+ }\r
+\r
+\r
+ /////////////////////////////////////////////////////////////////////////////\r
+ // Add candidate to AliFlowEvent for POI and subtract from RP if needed\r
+ ////////////////////////////////////////////////////////////////////////////\r
+ Int_t idtrack = static_cast<AliVTrack*>(track)->GetID();\r
+ Bool_t found = kFALSE;\r
+ Int_t numberoffound = 0;\r
+ //printf("A: Number of tracks %d\n",fflowEvent->NumberOfTracks());\r
+ for(Int_t iRPs=0; iRPs< fflowEvent->NumberOfTracks(); iRPs++) {\r
+ AliFlowTrack *iRP = (AliFlowTrack*) (fflowEvent->GetTrack(iRPs));\r
+ //if(!iRP->InRPSelection()) continue;\r
+ if( TMath::Abs(idtrack) == TMath::Abs(iRP->GetID()) ) {\r
+ iRP->SetForPOISelection(kTRUE);\r
+ found = kTRUE;\r
+ numberoffound ++;\r
+ }\r
+ }\r
+ //printf("Found %d mal\n",numberoffound);\r
+ if(!found) {\r
+ AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*) MakeTrack(massElectron,track->Pt(),track->Phi(), track->Eta());\r
+ sTrack->SetID(idtrack);\r
+ fflowEvent->AddTrack(sTrack);\r
+ //printf("Add the track\n");\r
+ }\r
+ //printf("B: Number of tracks %d\n",fflowEvent->NumberOfTracks());\r
+ \r
+ \r
+ \r
/////////////////////\r
// Fill THnSparseF\r
/////////////////////\r
\r
//\r
valuensparseabis[0] = eventplanesubtracted;\r
- if(fillEventPlane) fEventPlaneaftersubtraction->Fill(&valuensparseabis[0]);\r
+ if((fillEventPlane) && (fDebugLevel > 5)) fEventPlaneaftersubtraction->Fill(&valuensparseabis[0]);\r
\r
\r
- if(fDebugLevel > 3) \r
+ if(fDebugLevel > 5) \r
{\r
//\r
valuensparsee[0] = TMath::Cos(2*phitrack);\r
}\r
}\r
\r
- \r
+ if(fDebugLevel > 1) {\r
+ // background\r
+ Int_t source = 0;\r
+ Int_t indexmother = -1;\r
+ source = FindMother(TMath::Abs(track->GetLabel()),mcEvent, indexmother);\r
+ valuensparseMCSourceDeltaPhiMaps[2] = source;\r
+ if(mcEvent) fMCSourceDeltaPhiMaps->Fill(&valuensparseMCSourceDeltaPhiMaps[0]);\r
+ Int_t taggedvalue = LookAtNonHFE(k,track,fInputEvent,mcEvent,binct,deltaphi,source,indexmother);\r
+ if(fillEventPlane) {\r
+ // No opposite charge partner found in the invariant mass choosen\r
+ if((taggedvalue!=2) && (taggedvalue!=6)) {\r
+ fDeltaPhiMapsTaggedNonPhotonic->Fill(&valuensparseg[0]);\r
+ fCosPhiMapsTaggedNonPhotonic->Fill(&valuensparseh[0]);\r
+ }\r
+ // One opposite charge partner found in the invariant mass choosen\r
+ if((taggedvalue==2) || (taggedvalue==6)) {\r
+ fDeltaPhiMapsTaggedPhotonic->Fill(&valuensparseg[0]);\r
+ fCosPhiMapsTaggedPhotonic->Fill(&valuensparseh[0]);\r
+ }\r
+ // One same charge partner found in the invariant mass choosen\r
+ if((taggedvalue==4) || (taggedvalue==6)) {\r
+ fDeltaPhiMapsTaggedPhotonicLS->Fill(&valuensparseg[0]);\r
+ fCosPhiMapsTaggedPhotonicLS->Fill(&valuensparseh[0]);\r
+ }\r
+ }\r
+ }\r
\r
}\r
\r
for(Int_t bincless = 0; bincless < fNbBinsCentralityQCumulant; bincless++) {\r
if((fBinCentralityLess[bincless]< cntr) && (cntr < fBinCentralityLess[bincless+1])) PostData(bincless+2,fflowEvent);\r
}\r
+\r
+ if(fArraytrack) {\r
+ delete fArraytrack;\r
+ fArraytrack = NULL;\r
+ }\r
\r
PostData(1, fListHist);\r
\r
-\r
\r
}\r
//______________________________________________________________________________\r
}\r
return phiend;\r
}\r
+//_____________________________________________________________________________________________\r
+Int_t AliAnalysisTaskHFEFlow::LookAtNonHFE(Int_t iTrack1, AliVTrack *track1, AliVEvent *vEvent, AliMCEvent *mcEvent,Int_t binct,Double_t deltaphi,Int_t source,Int_t indexmother)\r
+{ \r
+ //\r
+ // Look At Non HFE\r
+ //\r
+\r
+ // return -1 if nothing\r
+ // return 2 if opposite charge within the mass range found\r
+ // return 4 if like charge within the mass range found\r
+ // return 6 if opposite charge and like charge within the mass range found\r
+ //\r
+\r
+ Int_t taggedphotonic = -1;\r
+\r
+ Bool_t oppositetaggedphotonic = kFALSE;\r
+ Bool_t sametaggedphotonic = kFALSE;\r
+\r
+ //printf("fCounterPoolBackground %d in LookAtNonHFE!!!\n",fCounterPoolBackground);\r
+ if(!fArraytrack) return taggedphotonic;\r
+ //printf("process track %d\n",iTrack1);\r
+ \r
+ TVector3 v3D1;\r
+ TVector3 v3D2;\r
+\r
+ Double_t valuensparseDeltaPhiMaps[5];\r
+ Double_t valueangle[3];\r
+\r
+ valuensparseDeltaPhiMaps[1] = binct;\r
+ valuensparseDeltaPhiMaps[2] = track1->Pt();\r
+ valuensparseDeltaPhiMaps[0] = deltaphi;\r
+ valuensparseDeltaPhiMaps[4] = source;\r
+ \r
+ valueangle[2] = source;\r
+ valueangle[1] = binct;\r
+\r
+ // Pdg code\r
+ Int_t pdg1 = CheckPdg(TMath::Abs(track1->GetLabel()),mcEvent);\r
+ Int_t numberfound = 0;\r
+\r
+ //Magnetic Field\r
+ Double_t bfield = vEvent->GetMagneticField();\r
+\r
+ // Get Primary vertex\r
+ const AliVVertex *pVtx = vEvent->GetPrimaryVertex();\r
+ \r
+ for(Int_t idex = 0; idex < fCounterPoolBackground; idex++) \r
+ {\r
+\r
+ Int_t iTrack2 = fArraytrack->At(idex);\r
+ //printf("track %d\n",iTrack2);\r
+ AliVTrack* track2 = (AliVTrack *) vEvent->GetTrack(iTrack2);\r
+ if (!track2) \r
+ {\r
+ printf("ERROR: Could not receive track %d\n", iTrack2);\r
+ continue;\r
+ }\r
+ if(iTrack2==iTrack1) continue;\r
+ //printf("Different\n");\r
+\r
+ // Reset the MC info\r
+ valueangle[2] = source;\r
+ valuensparseDeltaPhiMaps[4] = source;\r
+\r
+ // track cuts and PID already done\r
+\r
+ // if MC look\r
+ Int_t pdg2 = -100;\r
+ if(mcEvent) {\r
+ Int_t source2 = 0;\r
+ Int_t indexmother2 = -1;\r
+ source2 = FindMother(TMath::Abs(track2->GetLabel()),mcEvent, indexmother2);\r
+ pdg2 = CheckPdg(TMath::Abs(track2->GetLabel()),mcEvent);\r
+ if(source2 >=0 ) {\r
+ if((indexmother2 == indexmother) && (source == source2) && ((pdg1*pdg2)<0.0)) {\r
+ if(source == kElectronfromconversion) {\r
+ valueangle[2] = kElectronfromconversionboth;\r
+ valuensparseDeltaPhiMaps[4] = kElectronfromconversionboth;\r
+ numberfound++;\r
+ }\r
+ if(source == kElectronfrompi0) {\r
+ valueangle[2] = kElectronfrompi0both;\r
+ valuensparseDeltaPhiMaps[4] = kElectronfrompi0both;\r
+ }\r
+ if(source == kElectronfrometa) {\r
+ valueangle[2] = kElectronfrometaboth;\r
+ valuensparseDeltaPhiMaps[4] = kElectronfrometaboth;\r
+ }\r
+ }\r
+ }\r
+ }\r
+ \r
+ if(fAlgorithmMA && (!fAODAnalysis))\r
+ {\r
+ // tracks\r
+ AliESDtrack *esdtrack2 = dynamic_cast<AliESDtrack *>(track2); \r
+ AliESDtrack *esdtrack1 = dynamic_cast<AliESDtrack *>(track1); \r
+ if((!esdtrack2) || (!esdtrack1)) continue;\r
+\r
+ //Variables\r
+ Double_t p1[3];\r
+ Double_t p2[3];\r
+ Double_t xt1; //radial position track 1 at the DCA point\r
+ Double_t xt2; //radial position track 2 at the DCA point\r
+ //DCA track1-track2\r
+ Double_t dca12 = esdtrack2->GetDCA(esdtrack1,bfield,xt2,xt1);\r
+\r
+ // Cut dca\r
+ if(dca12 > fMaxdca) continue;\r
+ \r
+ //Momento of the track extrapolated to DCA track-track \r
+ //Track1\r
+ Bool_t hasdcaT1 = esdtrack1->GetPxPyPzAt(xt1,bfield,p1);\r
+ //Track2\r
+ Bool_t hasdcaT2 = esdtrack2->GetPxPyPzAt(xt2,bfield,p2);\r
+ \r
+ if(!hasdcaT1 || !hasdcaT2) AliWarning("It could be a problem in the extrapolation");\r
+ \r
+ //track1-track2 Invariant Mass\r
+ Double_t eMass = 0.000510998910; //Electron mass in GeV\r
+ Double_t pP1 = sqrt(p1[0]*p1[0]+p1[1]*p1[1]+p1[2]*p1[2]); //Track 1 momentum\r
+ Double_t pP2 = sqrt(p2[0]*p2[0]+p2[1]*p2[1]+p2[2]*p2[2]); //Track 2 momentum\r
+ Double_t eE1 = TMath::Sqrt(pP1*pP1+eMass*eMass);\r
+ Double_t eE2 = TMath::Sqrt(pP2*pP2+eMass*eMass);\r
+ \r
+ //TLorentzVector v1(p1[0],p1[1],p1[2],sqrt(eMass*eMass+pP1*pP1));\r
+ //TLorentzVector v2(p2[0],p2[1],p2[2],sqrt(eMass*eMass+pP2*pP2));\r
+ //Double_t imass = (v1+v2).M(); //Invariant Mass\r
+ //Double_t angle3D = v1.Angle(v2.Vect()); //Opening Angle (Total Angle)\r
+ \r
+ // daughter\r
+ v3D1.SetXYZ(p1[0],p1[1],p1[2]);\r
+ v3D2.SetXYZ(p2[0],p2[1],p2[2]);\r
+ Double_t openingangle = TVector2::Phi_0_2pi(v3D2.Angle(v3D1));\r
+ \r
+ // mother\r
+ TVector3 motherrec = v3D1 + v3D2;\r
+ Double_t invmass = TMath::Sqrt((eE1+eE2)*(eE1+eE2)-(motherrec.Px()*motherrec.Px()+motherrec.Py()*motherrec.Py()+motherrec.Pz()*motherrec.Pz()));\r
+ \r
+ // xy\r
+ //TVector3 vectordiff = v3D1 - v3D2;\r
+ //Double_t diffphi = TVector2::Phi_0_2pi(vectordiff.Phi());\r
+ //Double_t massxy = TMath::Sqrt((eE1+eE2)*(eE1+eE2)-(pP1*pP1+pP2*pP2+2*pP1*pP2*TMath::Cos(diffphi)));\r
+\r
+ // rz\r
+ //Double_t difftheta = TVector2::Phi_0_2pi(vectordiff.Eta());\r
+ //Double_t massrz = TMath::Sqrt((eE1+eE2)*(eE1+eE2)-(pP1*pP1+pP2*pP2+2*pP1*pP2*TMath::Cos(difftheta)));\r
+ \r
+\r
+ Float_t fCharge1 = track1->Charge();\r
+ Float_t fCharge2 = track2->Charge();\r
+\r
+ // Fill Histo\r
+ //valueangle[0] = diffphi;\r
+ //valueangle[1] = difftheta;\r
+ valueangle[0] = openingangle;\r
+ if((fCharge1*fCharge2)>0.0) fSameSignAngle->Fill(&valueangle[0]);\r
+ else fOppSignAngle->Fill(&valueangle[0]);\r
+\r
+ // Cut\r
+ if(openingangle > fMaxopening3D) continue;\r
+ //if(difftheta > fMaxopeningtheta) continue;\r
+ //if(diffphi > fMaxopeningphi) continue;\r
+\r
+ // Invmass\r
+ valuensparseDeltaPhiMaps[3] = invmass;\r
+ if((fCharge1*fCharge2)>0.0) fSameSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);\r
+ else fOppSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);\r
+ \r
+ // Cut\r
+ if(invmass < fMaxInvmass) {\r
+ if((fCharge1*fCharge2)<0.0) oppositetaggedphotonic=kTRUE;\r
+ if((fCharge1*fCharge2)>0.0) sametaggedphotonic=kTRUE;\r
+ }\r
+\r
+\r
+ }\r
+ else \r
+ {\r
+ Int_t fPDGtrack1 = 11; \r
+ Int_t fPDGtrack2 = 11;\r
+ \r
+ Float_t fCharge1 = track1->Charge();\r
+ Float_t fCharge2 = track2->Charge();\r
+ \r
+ if(fCharge1>0) fPDGtrack1 = -11;\r
+ if(fCharge2>0) fPDGtrack2 = -11;\r
+ \r
+ AliKFParticle ktrack1(*track1, fPDGtrack1);\r
+ AliKFParticle ktrack2(*track2, fPDGtrack2);\r
+ AliKFParticle recoGamma(ktrack1, ktrack2);\r
+ \r
+ //Reconstruction Cuts\r
+ if(recoGamma.GetNDF()<1) continue;\r
+ Double_t chi2OverNDF = recoGamma.GetChi2()/recoGamma.GetNDF();\r
+ if(TMath::Sqrt(TMath::Abs(chi2OverNDF))>fChi2OverNDFCut) continue;\r
+\r
+ // DCA\r
+ //Double_t dca12 = ktrack1.GetDistanceFromParticle(ktrack2);\r
+ //if(dca12 > fMaxdca) continue; \r
+\r
+ // if set mass constraint\r
+ if(fSetMassConstraint && pVtx) {\r
+ AliKFVertex primV(*pVtx);\r
+ primV += recoGamma;\r
+ primV -= ktrack1;\r
+ primV -= ktrack2;\r
+ recoGamma.SetProductionVertex(primV);\r
+ recoGamma.SetMassConstraint(0,0.0001);\r
+ } \r
+\r
+ //Invariant Mass\r
+ Double_t imass; \r
+ Double_t width;\r
+ recoGamma.GetMass(imass,width);\r
+ \r
+ //Opening Angle (Total Angle)\r
+ Double_t angle = ktrack1.GetAngle(ktrack2);\r
+ valueangle[0] = angle;\r
+ if((fCharge1*fCharge2)>0.0) fSameSignAngle->Fill(&valueangle[0]);\r
+ else fOppSignAngle->Fill(&valueangle[0]);\r
+\r
+ // Cut\r
+ if(angle > fMaxopening3D) continue; \r
+\r
+ // Invmass\r
+ valuensparseDeltaPhiMaps[3] = imass;\r
+ if((fCharge1*fCharge2)>0.0) fSameSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);\r
+ else {\r
+ fOppSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);\r
+ /*\r
+ if(valueangle[2] == kElectronfromconversionboth) {\r
+ printf("Reconstructed charge1 %f, charge 2 %f and invmass %f\n",fCharge1,fCharge2,imass);\r
+ printf("MC charge1 %d, charge 2 %d\n",pdg1,pdg2);\r
+ printf("DCA %f\n",dca12);\r
+ printf("Number of found %d\n",numberfound);\r
+ }\r
+ */\r
+ }\r
+ \r
+ // Cut\r
+ if(imass < fMaxInvmass) {\r
+ if((fCharge1*fCharge2)<0.0) oppositetaggedphotonic=kTRUE;\r
+ if((fCharge1*fCharge2)>0.0) sametaggedphotonic=kTRUE;\r
+ }\r
+ }\r
+ }\r
+ \r
+ if(oppositetaggedphotonic && sametaggedphotonic){\r
+ taggedphotonic = 6;\r
+ }\r
+\r
+ if(!oppositetaggedphotonic && sametaggedphotonic){\r
+ taggedphotonic = 4;\r
+ }\r
+\r
+ if(oppositetaggedphotonic && !sametaggedphotonic){\r
+ taggedphotonic = 2;\r
+ }\r
+\r
+ \r
+ return taggedphotonic;\r
+}\r
+//_________________________________________________________________________\r
+Int_t AliAnalysisTaskHFEFlow::FindMother(Int_t tr, AliMCEvent *mcEvent, Int_t &indexmother){\r
+ //\r
+ // Find the mother if MC\r
+ //\r
+\r
+ if(!mcEvent) return 0;\r
+\r
+ Int_t pdg = CheckPdg(tr,mcEvent);\r
+ if(TMath::Abs(pdg)!= 11) {\r
+ indexmother = -1;\r
+ return kNoElectron;\r
+ }\r
+ \r
+ indexmother = IsMotherGamma(tr,mcEvent);\r
+ if(indexmother > 0) return kElectronfromconversion;\r
+ indexmother = IsMotherPi0(tr,mcEvent);\r
+ if(indexmother > 0) return kElectronfrompi0;\r
+ indexmother = IsMotherC(tr,mcEvent);\r
+ if(indexmother > 0) return kElectronfromC;\r
+ indexmother = IsMotherB(tr,mcEvent);\r
+ if(indexmother > 0) return kElectronfromB;\r
+ indexmother = IsMotherEta(tr,mcEvent);\r
+ if(indexmother > 0) return kElectronfrometa;\r
+ \r
+ return kElectronfromother;\r
+\r
+\r
+}\r
+//____________________________________________________________________________________________________________\r
+Int_t AliAnalysisTaskHFEFlow::CheckPdg(Int_t tr, AliMCEvent* mcEvent) {\r
+\r
+ //\r
+ // Return the pdg of the particle\r
+ //\r
+\r
+\r
+ Int_t pdgcode = -1;\r
+ if(tr < 0) return pdgcode;\r
+\r
+ if(!mcEvent) return pdgcode;\r
+\r
+ AliVParticle *mctrack = mcEvent->GetTrack(tr);\r
+ \r
+ \r
+ if(mctrack->IsA() == AliMCParticle::Class()) {\r
+ AliMCParticle *mctrackesd = NULL;\r
+ if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return pdgcode;\r
+ pdgcode = mctrackesd->PdgCode();\r
+ }\r
+\r
+ if(mctrack->IsA() == AliAODMCParticle::Class()) {\r
+ AliAODMCParticle *mctrackaod = NULL;\r
+ if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return pdgcode;\r
+ pdgcode = mctrackaod->GetPdgCode();\r
+ }\r
+ \r
+ return pdgcode;\r
+\r
+ \r
+}\r
+//____________________________________________________________________________________________________________\r
+Int_t AliAnalysisTaskHFEFlow::IsMotherGamma(Int_t tr, AliMCEvent* mcEvent) {\r
+\r
+ //\r
+ // Return the lab of gamma mother or -1 if not gamma\r
+ //\r
+\r
+ if(tr < 0) return -1;\r
+ AliVParticle *mctrack = mcEvent->GetTrack(tr);\r
+ \r
+ if(mctrack->IsA() == AliMCParticle::Class()) {\r
+ AliMCParticle *mctrackesd = NULL;\r
+ if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
+ TParticle *particle = 0x0;\r
+ particle = mctrackesd->Particle();\r
+ // Take mother\r
+ if(!particle) return -1;\r
+ Int_t imother = particle->GetFirstMother(); \r
+ if(imother < 0) return -1; \r
+ AliMCParticle *mothertrack = NULL;\r
+ if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
+ TParticle * mother = mothertrack->Particle();\r
+ if(!mother) return -1;\r
+ // Check gamma \r
+ Int_t pdg = mother->GetPdgCode();\r
+ if(TMath::Abs(pdg) == 22) return imother;\r
+ if(TMath::Abs(pdg) == 11) {\r
+ return IsMotherGamma(imother,mcEvent);\r
+ }\r
+ return -1;\r
+ }\r
+\r
+ if(mctrack->IsA() == AliAODMCParticle::Class()) {\r
+ AliAODMCParticle *mctrackaod = NULL;\r
+ if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
+ // Take mother\r
+ Int_t imother = mctrackaod->GetMother();\r
+ if(imother < 0) return -1; \r
+ AliAODMCParticle *mothertrack = NULL;\r
+ if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
+ // Check gamma \r
+ Int_t pdg = mothertrack->GetPdgCode();\r
+ if(TMath::Abs(pdg) == 22) return imother;\r
+ if(TMath::Abs(pdg) == 11) {\r
+ return IsMotherGamma(imother,mcEvent);\r
+ }\r
+ return -1;\r
+\r
+ }\r
+ \r
+ return -1;\r
+\r
+ \r
+}\r
+//\r
+//____________________________________________________________________________________________________________\r
+Int_t AliAnalysisTaskHFEFlow::IsMotherPi0(Int_t tr, AliMCEvent* mcEvent) {\r
+\r
+ //\r
+ // Return the lab of pi0 mother or -1 if not pi0\r
+ //\r
+\r
+ if(tr < 0) return -1;\r
+ AliVParticle *mctrack = mcEvent->GetTrack(tr);\r
+ \r
+ if(mctrack->IsA() == AliMCParticle::Class()) {\r
+ AliMCParticle *mctrackesd = NULL;\r
+ if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
+ TParticle *particle = 0x0;\r
+ particle = mctrackesd->Particle();\r
+ // Take mother\r
+ if(!particle) return -1;\r
+ Int_t imother = particle->GetFirstMother(); \r
+ if(imother < 0) return -1; \r
+ AliMCParticle *mothertrack = NULL;\r
+ if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
+ TParticle * mother = mothertrack->Particle();\r
+ if(!mother) return -1;\r
+ // Check gamma \r
+ Int_t pdg = mother->GetPdgCode();\r
+ if(TMath::Abs(pdg) == 111) return imother;\r
+ if(TMath::Abs(pdg) == 11) {\r
+ return IsMotherPi0(imother,mcEvent);\r
+ }\r
+ return -1;\r
+ }\r
+\r
+ if(mctrack->IsA() == AliAODMCParticle::Class()) {\r
+ AliAODMCParticle *mctrackaod = NULL;\r
+ if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
+ // Take mother\r
+ Int_t imother = mctrackaod->GetMother();\r
+ if(imother < 0) return -1; \r
+ AliAODMCParticle *mothertrack = NULL;\r
+ if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
+ // Check gamma \r
+ Int_t pdg = mothertrack->GetPdgCode();\r
+ if(TMath::Abs(pdg) == 111) return imother;\r
+ if(TMath::Abs(pdg) == 11) {\r
+ return IsMotherPi0(imother,mcEvent);\r
+ }\r
+ return -1;\r
+ }\r
+\r
+ return -1;\r
+ \r
+}\r
+//____________________________________________________________________________________________________________\r
+Int_t AliAnalysisTaskHFEFlow::IsMotherC(Int_t tr, AliMCEvent* mcEvent) {\r
+\r
+ //\r
+ // Return the lab of signal mother or -1 if not signal\r
+ //\r
+\r
+ if(tr < 0) return -1;\r
+ AliVParticle *mctrack = mcEvent->GetTrack(tr);\r
+ \r
+ if(mctrack->IsA() == AliMCParticle::Class()) {\r
+ AliMCParticle *mctrackesd = NULL;\r
+ if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
+ TParticle *particle = 0x0;\r
+ particle = mctrackesd->Particle();\r
+ // Take mother\r
+ if(!particle) return -1;\r
+ Int_t imother = particle->GetFirstMother(); \r
+ if(imother < 0) return -1; \r
+ AliMCParticle *mothertrack = NULL;\r
+ if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
+ TParticle * mother = mothertrack->Particle();\r
+ if(!mother) return -1;\r
+ // Check gamma \r
+ Int_t pdg = mother->GetPdgCode();\r
+ if((TMath::Abs(pdg)==411) || (TMath::Abs(pdg)==421) || (TMath::Abs(pdg)==431) || (TMath::Abs(pdg)==4122) || (TMath::Abs(pdg)==4132) || (TMath::Abs(pdg)==4232) || (TMath::Abs(pdg)==43320)) return imother;\r
+ if(TMath::Abs(pdg) == 11) {\r
+ return IsMotherC(imother,mcEvent);\r
+ }\r
+ return -1;\r
+ }\r
+\r
+ if(mctrack->IsA() == AliAODMCParticle::Class()) {\r
+ AliAODMCParticle *mctrackaod = NULL;\r
+ if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
+ // Take mother\r
+ Int_t imother = mctrackaod->GetMother();\r
+ if(imother < 0) return -1; \r
+ AliAODMCParticle *mothertrack = NULL;\r
+ if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
+ // Check gamma \r
+ Int_t pdg = mothertrack->GetPdgCode();\r
+ if((TMath::Abs(pdg)==411) || (TMath::Abs(pdg)==421) || (TMath::Abs(pdg)==431) || (TMath::Abs(pdg)==4122) || (TMath::Abs(pdg)==4132) || (TMath::Abs(pdg)==4232) || (TMath::Abs(pdg)==43320)) return imother;\r
+ if(TMath::Abs(pdg) == 11) {\r
+ return IsMotherC(imother,mcEvent);\r
+ }\r
+ return -1;\r
+ }\r
+\r
+ return -1;\r
+\r
+}\r
+//____________________________________________________________________________________________________________\r
+Int_t AliAnalysisTaskHFEFlow::IsMotherB(Int_t tr, AliMCEvent* mcEvent) {\r
+\r
+ //\r
+ // Return the lab of signal mother or -1 if not signal\r
+ //\r
+\r
+ if(tr < 0) return -1;\r
+ AliVParticle *mctrack = mcEvent->GetTrack(tr);\r
+ \r
+ if(mctrack->IsA() == AliMCParticle::Class()) {\r
+ AliMCParticle *mctrackesd = NULL;\r
+ if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
+ TParticle *particle = 0x0;\r
+ particle = mctrackesd->Particle();\r
+ // Take mother\r
+ if(!particle) return -1;\r
+ Int_t imother = particle->GetFirstMother(); \r
+ if(imother < 0) return -1; \r
+ AliMCParticle *mothertrack = NULL;\r
+ if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
+ TParticle * mother = mothertrack->Particle();\r
+ if(!mother) return -1;\r
+ // Check gamma \r
+ Int_t pdg = mother->GetPdgCode();\r
+ if((TMath::Abs(pdg)==511) || (TMath::Abs(pdg)==521) || (TMath::Abs(pdg)==531) || (TMath::Abs(pdg)==5122) || (TMath::Abs(pdg)==5132) || (TMath::Abs(pdg)==5232) || (TMath::Abs(pdg)==53320)) return imother; \r
+ if(TMath::Abs(pdg) == 11) {\r
+ return IsMotherB(imother,mcEvent);\r
+ }\r
+ return -1;\r
+ }\r
+\r
+ if(mctrack->IsA() == AliAODMCParticle::Class()) {\r
+ AliAODMCParticle *mctrackaod = NULL;\r
+ if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
+ // Take mother\r
+ Int_t imother = mctrackaod->GetMother();\r
+ if(imother < 0) return -1; \r
+ AliAODMCParticle *mothertrack = NULL;\r
+ if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
+ // Check gamma \r
+ Int_t pdg = mothertrack->GetPdgCode();\r
+ if((TMath::Abs(pdg)==511) || (TMath::Abs(pdg)==521) || (TMath::Abs(pdg)==531) || (TMath::Abs(pdg)==5122) || (TMath::Abs(pdg)==5132) || (TMath::Abs(pdg)==5232) || (TMath::Abs(pdg)==53320)) return imother;\r
+ if(TMath::Abs(pdg) == 11) {\r
+ return IsMotherB(imother,mcEvent);\r
+ }\r
+ return -1;\r
+ }\r
+\r
+ return -1;\r
+\r
+}\r
+//____________________________________________________________________________________________________________\r
+Int_t AliAnalysisTaskHFEFlow::IsMotherEta(Int_t tr, AliMCEvent* mcEvent) {\r
+\r
+ //\r
+ // Return the lab of pi0 mother or -1 if not pi0\r
+ //\r
+\r
+ if(tr < 0) return -1;\r
+ AliVParticle *mctrack = mcEvent->GetTrack(tr);\r
+ \r
+ if(mctrack->IsA() == AliMCParticle::Class()) {\r
+ AliMCParticle *mctrackesd = NULL;\r
+ if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
+ TParticle *particle = 0x0;\r
+ particle = mctrackesd->Particle();\r
+ // Take mother\r
+ if(!particle) return -1;\r
+ Int_t imother = particle->GetFirstMother(); \r
+ if(imother < 0) return -1; \r
+ AliMCParticle *mothertrack = NULL;\r
+ if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
+ TParticle * mother = mothertrack->Particle();\r
+ if(!mother) return -1;\r
+ // Check gamma \r
+ Int_t pdg = mother->GetPdgCode();\r
+ if(TMath::Abs(pdg) == 221) return imother;\r
+ if(TMath::Abs(pdg) == 11) {\r
+ return IsMotherEta(imother,mcEvent);\r
+ }\r
+ return -1;\r
+ }\r
+\r
+ if(mctrack->IsA() == AliAODMCParticle::Class()) {\r
+ AliAODMCParticle *mctrackaod = NULL;\r
+ if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
+ // Take mother\r
+ Int_t imother = mctrackaod->GetMother();\r
+ if(imother < 0) return -1; \r
+ AliAODMCParticle *mothertrack = NULL;\r
+ if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
+ // Check gamma \r
+ Int_t pdg = mothertrack->GetPdgCode();\r
+ if(TMath::Abs(pdg) == 221) return imother;\r
+ if(TMath::Abs(pdg) == 11) {\r
+ return IsMotherEta(imother,mcEvent);\r
+ }\r
+ return -1;\r
+ }\r
+\r
+ return -1;\r
+ \r
+}\r
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff