memset(container, 0, sizeof(Double_t) * 10);
// container for the output THnSparse
Double_t dataE[6]; // [pT, eta, Phi, type, 'C' or 'B']
- Double_t dataDca[8]; // [source, pT, dca, dcaSig, centrality]
- Double_t eDca[6]; // [source, pT, dca, dcaSig, centrality]
+ Double_t dataDca[6]; // [source, pT, dca, centrality]
Int_t nElectronCandidates = 0;
AliESDtrack *track = NULL, *htrack = NULL;
AliMCParticle *mctrack = NULL;
}
}
- // check if it is the proton from the lambda decay, and yes fill the dca info
- Double_t hfeimpactR4p=0., hfeimpactnsigmaR4p=0.;
+ Double_t hfeimpactR4all=0., hfeimpactnsigmaR4all=0.;
+ Int_t sourceDca =-1;
if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) == 211)){
if(track->Pt()>4.){
- fExtraCuts->GetHFEImpactParameters(track, hfeimpactR4p, hfeimpactnsigmaR4p);
- fQACollection->Fill("pionDcaSig",track->Pt(),hfeimpactnsigmaR4p);
- }
- }
- else if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) == 2212)){
- fExtraCuts->GetHFEImpactParameters(track, hfeimpactR4p, hfeimpactnsigmaR4p);
- fQACollection->Fill("protonDcaSig",track->Pt(),hfeimpactnsigmaR4p);
- Int_t glabel=TMath::Abs(mctrack->GetMother());
- if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
- if(TMath::Abs(mctrackmother->Particle()->GetPdgCode())==3122){
- fQACollection->Fill("secondaryprotonDcaSig",track->Pt(),hfeimpactnsigmaR4p);
- }
+ fExtraCuts->GetHFEImpactParameters(track, hfeimpactR4all, hfeimpactnsigmaR4all);
+ dataDca[0]=0; //pion
+ dataDca[1]=track->Pt();
+ dataDca[2]=hfeimpactR4all;
+ dataDca[3]=fCentralityF;
+ dataDca[4] = v0pid;
+ dataDca[5] = double(track->Charge());
+ fQACollection->Fill("Dca", dataDca);
}
}
else if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) == 11)){ // to increas statistics for Martin
if(signal){
- fExtraCuts->GetHFEImpactParameters(track, hfeimpactR4p, hfeimpactnsigmaR4p);
- Int_t sourceDca =-1;
+ fExtraCuts->GetHFEImpactParameters(track, hfeimpactR4all, hfeimpactnsigmaR4all);
if(fSignalCuts->IsCharmElectron(track)){
sourceDca=1;
}
else if(fSignalCuts->IsJpsiElectron(track)){
sourceDca=5;
}
- else if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) != 11)){
- sourceDca=6;
- fQACollection->Fill("hadronsBeforeIPcut",track->Pt());
- fQACollection->Fill("hadronsBeforeIPcutMC",mctrack->Pt());
- }
else {
- sourceDca=7;
+ sourceDca=6;
}
- eDca[0]=sourceDca;
- eDca[1]=track->Pt();
- eDca[2]=hfeimpactR4p;
- eDca[3]=hfeimpactnsigmaR4p;
- eDca[4]=fCentralityF;
- eDca[5] = track->Charge()/3;
-
- fQACollection->Fill("eDca", eDca);
+ dataDca[0]=sourceDca;
+ dataDca[1]=track->Pt();
+ dataDca[2]=hfeimpactR4all;
+ dataDca[3]=fCentralityF;
+ dataDca[4] = v0pid;
+ dataDca[5] = double(track->Charge());
+ if(signal) fQACollection->Fill("Dca", dataDca);
}
}
}
}
} // end of electron background analysis
-
- Int_t sourceDca =-1;
if (GetPlugin(kDEstep)) {
Double_t weightElecBgV0[kBgLevels] = {0.,0.,0.,};
Int_t elecSource = 0;
- // minjung for IP QA(temporary ~ 2weeks)
Double_t hfeimpactR=0., hfeimpactnsigmaR=0.;
fExtraCuts->GetHFEImpactParameters(track, hfeimpactR, hfeimpactnsigmaR);
- sourceDca=0;
if(HasMCData())
{
- if(fSignalCuts->IsCharmElectron(track)){
- sourceDca=1;
- }
- else if(fSignalCuts->IsBeautyElectron(track)){
- sourceDca=2;
- }
- else if(fSignalCuts->IsGammaElectron(track)){
- sourceDca=3;
- }
- else if(fSignalCuts->IsNonHFElectron(track)){
- sourceDca=4;
- }
- else if(fSignalCuts->IsJpsiElectron(track)){
- sourceDca=5;
- }
- else if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) != 11)){
- sourceDca=6;
- fQACollection->Fill("hadronsBeforeIPcut",track->Pt());
- fQACollection->Fill("hadronsBeforeIPcutMC",mctrack->Pt());
- }
- else {
- sourceDca=7;
- }
-
+ if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) != 11)){
+ fQACollection->Fill("hadronsBeforeIPcut",track->Pt());
+ }
if(fMCQA && signal) {
fMCQA->SetContainerStep(0);
}
} // end of MC
- dataDca[0]=sourceDca;
+ dataDca[0]=-1; //for data, don't know the origin
dataDca[1]=track->Pt();
dataDca[2]=hfeimpactR;
- dataDca[3]=hfeimpactnsigmaR;
- dataDca[4]=fCentralityF;
- dataDca[5] = 49;
- Double_t xr[3]={49,49,49};
- if(HasMCData()) {
- mctrack->XvYvZv(xr);
- dataDca[5] = TMath::Sqrt(xr[0]*xr[0]+xr[1]*xr[1]);
- }
- dataDca[6] = v0pid;
- dataDca[7] = track->Charge()/3;
-
- // printf("Entries dca: [%.3f|%.3f|%.3f|%f|%f]\n", dataDca[0],dataDca[1],dataDca[2],dataDca[3],dataDca[4]);
+ dataDca[3]=fCentralityF;
+ dataDca[4] = v0pid;
+ dataDca[5] = double(track->Charge());
if (!HasMCData()) fQACollection->Fill("Dca", dataDca);
- else if(signal) fQACollection->Fill("Dca", dataDca);
-
// Fill Containers for impact parameter analysis
if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepHFEcutsDca + AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack,track)) continue;
if(HasMCData()){
if(mctrack && (TMath::Abs(mctrack->Particle()->GetPdgCode()) != 11)){
fQACollection->Fill("hadronsAfterIPcut",track->Pt());
- fQACollection->Fill("hadronsAfterIPcutMC",mctrack->Pt());
}
}
}
//printf("pass\n");
fContainer->NewEvent();
-
+
+ // Look for kink mother
+ Int_t numberofvertices = fAOD->GetNumberOfVertices();
+ Double_t listofmotherkink[numberofvertices];
+ Int_t numberofmotherkink = 0;
+ for(Int_t ivertex=0; ivertex < numberofvertices; ivertex++) {
+ AliAODVertex *aodvertex = fAOD->GetVertex(ivertex);
+ if(!aodvertex) continue;
+ if(aodvertex->GetType()==AliAODVertex::kKink) {
+ AliAODTrack *mother = (AliAODTrack *) aodvertex->GetParent();
+ if(!mother) continue;
+ Int_t idmother = mother->GetID();
+ listofmotherkink[numberofmotherkink] = idmother;
+ //printf("ID %d\n",idmother);
+ numberofmotherkink++;
+ }
+ }
+ //printf("Number of kink mother in the events %d\n",numberofmotherkink);
+
+
+ // Loop over tracks
AliAODTrack *track = NULL;
AliAODMCParticle *mctrack = NULL;
Double_t dataE[6]; // [pT, eta, Phi, Charge, type, 'C' or 'B']
if(fApplyCutAOD) {
// RecKine: ITSTPC cuts
if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
+
+ // Reject kink mother
+ Bool_t kinkmotherpass = kTRUE;
+ for(Int_t kinkmother = 0; kinkmother < numberofmotherkink; kinkmother++) {
+ if(track->GetID() == listofmotherkink[kinkmother]) {
+ kinkmotherpass = kFALSE;
+ continue;
+ }
+ }
+ if(!kinkmotherpass) continue;
// RecPrim
if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
fQACollection->CreateTH1Farray("hadronsBeforeIPcut", "Hadrons before IP cut", nBinPt, kPtRange);
fQACollection->CreateTH1Farray("hadronsAfterIPcut", "Hadrons after IP cut", nBinPt, kPtRange);
- fQACollection->CreateTH1Farray("hadronsBeforeIPcutMC", "Hadrons before IP cut; MC p_{t}", nBinPt, kPtRange);
- fQACollection->CreateTH1Farray("hadronsAfterIPcutMC", "Hadrons after IP cut; MC p_{t} ", nBinPt, kPtRange);
fQACollection->CreateTH2Farray("Ke3Kecorr", "Ke3 decay e and K correlation; Ke3K p_{t}; Ke3e p_{t}; ", nBinPt, kPtRange, 20,0.,20.);
fQACollection->CreateTH2Farray("Ke3K0Lecorr", "Ke3 decay e and K0L correlation; Ke3K0L p_{t}; Ke3e p_{t}; ", nBinPt, kPtRange, 20,0.,20.);
fQACollection->CreateTH1Farray("Kptspectra", "Charged Kaons: MC p_{t} ", nBinPt, kPtRange);
fQACollection->CreateTH1Farray("K0Lptspectra", "K0L: MC p_{t} ", nBinPt, kPtRange);
- const Double_t kDCAbound[2] = {-5., 5.};
- const Double_t kDCAsigbound[2] = {-50., 50.};
+ const Double_t kDCAbound[2] = {-0.2, 0.2};
- const Int_t nDimDca=8;
- const Int_t nBinDca[nDimDca] = { 9, nBinPt, 2000, 2000, 12, 500, 6, 2};
- const Int_t nDimeDca=6;
- const Int_t nBineDca[nDimeDca] = { 9, nBinPt, 2000, 2000, 12, 2};
- Double_t minimaDca[nDimDca] = { -1., 0., kDCAbound[0], kDCAsigbound[0], -1., 0, -1, -1.1};
- Double_t maximaDca[nDimDca] = { 8., 20., kDCAbound[1], kDCAsigbound[1], 11., 50, 5, 1.1};
+ const Int_t nDimDca=6;
+ const Int_t nBinDca[nDimDca] = { 8, nBinPt, 800, 12, 6, 2};
+ Double_t minimaDca[nDimDca] = { -1., 0., kDCAbound[0], -1., -1, -1.1};
+ Double_t maximaDca[nDimDca] = { 7., 20., kDCAbound[1], 11., 5, 1.1};
Double_t *sourceBins = AliHFEtools::MakeLinearBinning(nBinDca[0], minimaDca[0], maximaDca[0]);
Double_t *dcaBins = AliHFEtools::MakeLinearBinning(nBinDca[2], minimaDca[2], maximaDca[2]);
- Double_t *dcaSigBins = AliHFEtools::MakeLinearBinning(nBinDca[3], minimaDca[3], maximaDca[3]);
- Double_t *centralityBins = AliHFEtools::MakeLinearBinning(nBinDca[4], minimaDca[4], maximaDca[4]);
- Double_t *eProdRBins = AliHFEtools::MakeLinearBinning(nBinDca[5], minimaDca[5], maximaDca[5]);
- Double_t *v0PIDBins = AliHFEtools::MakeLinearBinning(nBinDca[6], minimaDca[6], maximaDca[6]);
- Double_t *chargeBins = AliHFEtools::MakeLinearBinning(nBinDca[7], minimaDca[7], maximaDca[7]);
+ Double_t *centralityBins = AliHFEtools::MakeLinearBinning(nBinDca[3], minimaDca[3], maximaDca[3]);
+ Double_t *v0PIDBins = AliHFEtools::MakeLinearBinning(nBinDca[4], minimaDca[4], maximaDca[4]);
+ Double_t *chargeBins = AliHFEtools::MakeLinearBinning(nBinDca[5], minimaDca[5], maximaDca[5]);
- fQACollection->CreateTHnSparseNoLimits("Dca", "Dca; source (0-all, 1-charm,etc); pT [GeV/c]; dca; dcasig; centrality bin; eProdR; v0pid; charge", nDimDca, nBinDca);
+ fQACollection->CreateTHnSparseNoLimits("Dca", "Dca; source (0-all, 1-charm,etc); pT [GeV/c]; dca; centrality bin; v0pid; charge", nDimDca, nBinDca);
((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(0, sourceBins);
((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(1, kPtRange);
((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(2, dcaBins);
- ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(3, dcaSigBins);
- ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(4, centralityBins);
- ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(5, eProdRBins);
- ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(6, v0PIDBins);
- ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(7, chargeBins);
-
- fQACollection->CreateTHnSparseNoLimits("eDca", "eDca; source (0-all, 1-charm,etc); pT [GeV/c]; dca; dcasig; centrality bin; charge", nDimeDca, nBineDca);
- ((THnSparse*)(fQACollection->Get("eDca")))->SetBinEdges(0, sourceBins);
- ((THnSparse*)(fQACollection->Get("eDca")))->SetBinEdges(1, kPtRange);
- ((THnSparse*)(fQACollection->Get("eDca")))->SetBinEdges(2, dcaBins);
- ((THnSparse*)(fQACollection->Get("eDca")))->SetBinEdges(3, dcaSigBins);
- ((THnSparse*)(fQACollection->Get("eDca")))->SetBinEdges(4, centralityBins);
- ((THnSparse*)(fQACollection->Get("eDca")))->SetBinEdges(5, chargeBins);
-
- fQACollection->CreateTH2Farray("secondaryprotonDcaSig", "secondary proton dca significance: dca sig ", nBinPt, kPtRange, 2000, kDCAsigbound[0], kDCAsigbound[1]);
- fQACollection->CreateTH2Farray("protonDcaSig", "proton dca significance: dca sig ", nBinPt, kPtRange, 2000, kDCAsigbound[0], kDCAsigbound[1]);
- fQACollection->CreateTH2Farray("pionDcaSig", "pion dca significance: dca sig ", nBinPt, kPtRange, 2000, kDCAsigbound[0], kDCAsigbound[1]);
+ ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(3, centralityBins);
+ ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(4, v0PIDBins);
+ ((THnSparse*)(fQACollection->Get("Dca")))->SetBinEdges(5, chargeBins);
break;
}
#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 "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
\r
#include "AliHFEcuts.h"\r
#include "AliHFEpid.h"\r
fUseMCReactionPlane(kFALSE),\r
fMCPID(kFALSE),\r
fNoPID(kFALSE),\r
+ fChi2OverNDFCut(999),\r
+ fMaxdca(3.0),\r
+ fMaxopeningtheta(0.02),\r
+ fMaxopeningphi(0.1),\r
+ fMaxopening3D(0.1),\r
+ fMaxInvmass(0.1),\r
fDebugLevel(0),\r
fcutsRP(0),\r
fcutsPOI(0),\r
fHFECuts(0),\r
fPID(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
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
fUseMCReactionPlane(kFALSE),\r
fMCPID(kFALSE),\r
fNoPID(kFALSE),\r
+ fChi2OverNDFCut(999),\r
+ fMaxdca(3.0),\r
+ fMaxopeningtheta(0.02),\r
+ fMaxopeningphi(0.1),\r
+ fMaxopening3D(0.1),\r
+ fMaxInvmass(0.1),\r
fDebugLevel(0),\r
fcutsRP(0),\r
fcutsPOI(0),\r
fHFECuts(0),\r
fPID(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
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
DefineInput(0,TChain::Class());\r
DefineOutput(1, TList::Class());\r
for(Int_t bincless = 0; bincless < fNbBinsCentralityQCumulant; bincless++) {\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
fDebugLevel(ref.fDebugLevel),\r
fcutsRP(0),\r
fcutsPOI(0),\r
fHFECuts(0),\r
fPID(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
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
// Copy Constructor\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.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(fPID) delete fPID;\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
fPID->InitializePID();\r
fPIDqa->Initialize(fPID);\r
fPID->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
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 = 6;\r
+ Double_t minStep = 0.;\r
+ Double_t maxStep = 6.;\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 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 = 8;\r
+ Double_t minSource = 0.;\r
+ Double_t maxSource = 8.;\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
+ 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(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
+ \r
AliMCEvent *mcEvent = MCEvent();\r
AliMCParticle *mctrack = NULL;\r
- \r
+ \r
/////////////////\r
// centrality\r
/////////////////\r
valuensparseh[1] = binct; \r
valuensparsehprofile[1] = binct; \r
valuecossinephiep[2] = binctMore;\r
-\r
+ valuensparseMCSourceDeltaPhiMaps[0] = binct;\r
+ \r
//////////////////////\r
// run number\r
//////////////////////\r
fPID->InitializePID(runnumber);\r
}\r
\r
+ if(!fPIDBackground->IsInitialized()){\r
+ // Initialize PID with the given run number\r
+ fPIDBackground->InitializePID(runnumber);\r
+ }\r
+\r
\r
//////////\r
// PID\r
return;\r
}\r
fPID->SetPIDResponse(pidResponse);\r
+ fPIDBackground->SetPIDResponse(pidResponse);\r
\r
fHistEV->Fill(binctt,0.0);\r
\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 pos[3],cov[6];\r
+ vAOD->GetXYZ(pos);\r
+ vAOD->GetCovarianceMatrix(cov);\r
+ const AliESDVertex vESD(pos,cov,100.,100);\r
+ esdTrack.RelateToVertex(&vESD,0.,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
+ // 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
+ }\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(fAODAnalysis) {\r
AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);\r
if(!aodtrack){\r
AliError("AOD track is not there");\r
- return;\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
- //printf("Check flag on\n");\r
}\r
}\r
\r
if(fApplyCut) {\r
- //printf("Apply cut\n");\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
- //printf("Cut %d\n",icut);\r
- break;\r
- }\r
- //printf("Pass the cut %d\n",icut);\r
- }\r
- if(!survived) continue;\r
- }\r
\r
- //printf("Survived\n");\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
+ 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] = 0; \r
+ if(fDebugLevel > 3) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
+ // RecPrim\r
+ if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepRecPrim + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
+ valuetrackingcuts[1] = 1; \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] = 2; \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] = 3; \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] = 4; \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] = 5; \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
+ // Define variables\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
+ 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 > 4) { \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
+ 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
+ }\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
+ //Magnetic Field\r
+ Double_t bfield = vEvent->GetMagneticField();\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
+ // track cuts and PID already done\r
+\r
+ // if MC look\r
+ if(mcEvent) {\r
+ Int_t source2 = 0;\r
+ Int_t indexmother2 = -1;\r
+ source2 = FindMother(TMath::Abs(track2->GetLabel()),mcEvent, indexmother2);\r
+ if((indexmother2 == indexmother) && (source == source2)) {\r
+ if(source == kElectronfromconversion) {\r
+ valueangle[2] = kElectronfromconversionboth;\r
+ valuensparseDeltaPhiMaps[4] = kElectronfromconversionboth;\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
+ 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 fKFtrack1(*track1, fPDGtrack1);\r
+ AliKFParticle fKFtrack2(*track2, fPDGtrack2);\r
+ AliKFParticle fRecoGamma(fKFtrack1, fKFtrack2);\r
+ \r
+ //Reconstruction Cuts\r
+ if(fRecoGamma.GetNDF()<1) continue;\r
+ Double_t chi2OverNDF = fRecoGamma.GetChi2()/fRecoGamma.GetNDF();\r
+ if(TMath::Sqrt(TMath::Abs(chi2OverNDF))>fChi2OverNDFCut) continue;\r
+ \r
+ //Invariant Mass\r
+ Double_t imass; \r
+ Double_t width;\r
+ fRecoGamma.GetMass(imass,width);\r
+ \r
+ //Opening Angle (Total Angle)\r
+ Double_t angle = fKFtrack1.GetAngle(fKFtrack2);\r
+ valueangle[0] = angle;\r
+ if((fCharge1*fCharge2)>0.0) fSameSignAngle->Fill(&valueangle[0]);\r
+ else fOppSignAngle->Fill(&valueangle[0]); \r
+\r
+ // Invmass\r
+ valuensparseDeltaPhiMaps[3] = imass;\r
+ if((fCharge1*fCharge2)>0.0) fSameSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);\r
+ else fOppSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);\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
+\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
+//_________________________________________________________________________\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
+ 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::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
#include <AliAnalysisTaskSE.h>\r
\r
class TList;\r
+class AliVTrack;\r
+class AliVEvent;\r
+class AliESDtrack;\r
+class AliESDEvent;\r
+class AliMCEvent;\r
class AliFlowTrackCuts;\r
class AliFlowCandidateTrack;\r
class AliHFEcuts;\r
class THnSparse;\r
class AliHFEpidQAmanager;\r
class AliFlowEvent;\r
+class AliESDtrackCuts;\r
class AliHFEVZEROEventPlane;\r
+class TArrayI;\r
\r
class AliAnalysisTaskHFEFlow: public AliAnalysisTaskSE {\r
\r
public:\r
+\r
+ typedef enum{\r
+ kElectronfromconversion = 0,\r
+ kElectronfromconversionboth = 1,\r
+ kElectronfrompi0 = 2,\r
+ kElectronfrompi0both = 3,\r
+ kElectronfromC = 4,\r
+ kElectronfromB = 5,\r
+ kElectronfrometa = 6,\r
+ kElectronfrometaboth = 7,\r
+ kElectronfromother = 8\r
+ } FlowSource_t;\r
+ \r
+ typedef enum{\r
+ kS = 0,\r
+ kOp = 1\r
+ } FlowSign_t;\r
+\r
+\r
+\r
+\r
AliAnalysisTaskHFEFlow();\r
AliAnalysisTaskHFEFlow(const char *name);\r
AliAnalysisTaskHFEFlow(const AliAnalysisTaskHFEFlow &ref);\r
\r
AliHFEpid *GetPID() const { return fPID; }\r
AliHFEpidQAmanager *GetPIDQAManager() const { return fPIDqa; }\r
+ AliHFEpid *GetPIDBackground() const { return fPIDBackground; }\r
+ AliHFEpidQAmanager *GetPIDBackgroundQAManager() const { return fPIDBackgroundqa; }\r
+\r
\r
void SetHFECuts(AliHFEcuts * const cuts) { fHFECuts = cuts; };\r
+ void SetHFEBackgroundCuts(AliESDtrackCuts * const cuts) { fHFEBackgroundCuts = cuts; };\r
void SetSubEtaGapTPC(Bool_t subEtaGapTPC) { fSubEtaGapTPC = subEtaGapTPC; };\r
void SetEtaGap(Double_t etaGap) { fEtaGap = etaGap; };\r
void SetVZEROEventPlane(Bool_t vzeroEventPlane) { fVZEROEventPlane = vzeroEventPlane; };\r
\r
AliFlowCandidateTrack *MakeTrack( Double_t mass, Double_t pt, Double_t phi, Double_t eta) ;\r
Double_t GetPhiAfterAddV2(Double_t phi,Double_t reactionPlaneAngle) const;\r
+\r
+ void SetMaxInvmass(Double_t maxInvmass) { fMaxInvmass = maxInvmass; };\r
+ void SetMaxopening3D(Double_t maxOpening3D) { fMaxopening3D = maxOpening3D; };\r
+ void SetMaxopeningtheta(Double_t maxOpeningtheta) { fMaxopeningtheta = maxOpeningtheta; };\r
+ void SetMaxopeningphi(Double_t maxOpeningphi) { fMaxopeningphi = maxOpeningphi; };\r
+\r
+ Int_t LookAtNonHFE(Int_t iTrack1, AliVTrack *track1, AliVEvent *fESD, AliMCEvent *mcEvent,Int_t binct,Double_t deltaphi,Int_t source,Int_t indexmother);\r
\r
private:\r
TList *fListHist; //! TH list\r
Bool_t fMCPID; // MC PID for electrons\r
Bool_t fNoPID; // No PID for checks\r
\r
+ Double_t fChi2OverNDFCut; // Limit chi2\r
+ Double_t fMaxdca; // Limit dca\r
+ Double_t fMaxopeningtheta; // Limit opening angle in theta\r
+ Double_t fMaxopeningphi; // Limit opening angle in phi\r
+ Double_t fMaxopening3D; // Limit opening 3D\r
+ Double_t fMaxInvmass; // Limit invariant mass\r
+ \r
+\r
Int_t fDebugLevel; // Debug Level \r
\r
// Cuts for FLOW PWG2\r
AliHFEcuts *fHFECuts; // HFE cuts\r
AliHFEpid *fPID; // PID cuts \r
AliHFEpidQAmanager *fPIDqa; // QA Manager\r
- AliFlowEvent *fflowEvent; //! Flow event \r
+ AliFlowEvent *fflowEvent; //! Flow event \r
+\r
+ // Cuts for background study\r
+ AliESDtrackCuts *fHFEBackgroundCuts; // HFE background cuts\r
+ AliHFEpid *fPIDBackground; // PID background cuts \r
+ AliHFEpidQAmanager *fPIDBackgroundqa; // QA Manager Background \r
+ Bool_t fAlgorithmMA; // algorithm MA\r
+\r
+ // List of tracks\r
+ TArrayI *fArraytrack; // list of tracks\r
+ Int_t fCounterPoolBackground; // number of tracks\r
\r
// VZERO Event plane after calibration 2010\r
AliHFEVZEROEventPlane *fHFEVZEROEventPlane; // VZERO event plane calibrated\r
THnSparseF *fSinRes; //! Res\r
TProfile *fProfileCosRes; //! Profile Res\r
\r
+ // Debuging Cuts step by step all centrality together: pt, step (6)\r
+ THnSparseF *fTrackingCuts; //! Tracking Cuts\r
+\r
+ // Before PID cut\r
+ // G Maps delta phi as function of deltaphi, centrality, pt\r
+ THnSparseF *fDeltaPhiMapsBeforePID; //! Delta phi\r
+ // H Maps cos phi : cos, centrality, pt\r
+ THnSparseF *fCosPhiMapsBeforePID; //! Cos\r
+\r
// G Maps delta phi as function of deltaphi, centrality, pt\r
THnSparseF *fDeltaPhiMaps; //! Delta phi\r
- \r
// H Maps cos phi : cos, centrality, pt\r
THnSparseF *fCosPhiMaps; //! Cos\r
TProfile2D *fProfileCosPhiMaps; //! Profile Cos\r
- \r
+\r
+ // Background study: not statistic but tagged \r
+ THnSparseF *fDeltaPhiMapsTaggedPhotonic; //! Delta phi\r
+ THnSparseF *fCosPhiMapsTaggedPhotonic; //! Cos\r
+ THnSparseF *fDeltaPhiMapsTaggedNonPhotonic; //! Delta phi\r
+ THnSparseF *fCosPhiMapsTaggedNonPhotonic; //! Cos\r
+ THnSparseF *fDeltaPhiMapsTaggedPhotonicLS; //! Delta phi\r
+ THnSparseF *fCosPhiMapsTaggedPhotonicLS; //! Cos\r
+\r
+ // Background study: centrality, pt, source\r
+ THnSparseF *fMCSourceDeltaPhiMaps; //! Source MC\r
+ // Background study: deltaphi, centrality, pt, minv, source\r
+ THnSparseF *fOppSignDeltaPhiMaps; //! Delta phi\r
+ THnSparseF *fSameSignDeltaPhiMaps; //! Delta phi\r
+ // Background study: angle, centrality, source\r
+ THnSparseF *fOppSignAngle; // ! Opening Angles\r
+ THnSparseF *fSameSignAngle; // ! Opening Angles\r
+\r
+ Int_t FindMother(Int_t tr, AliMCEvent *mcEvent, Int_t &indexmother);\r
+ Int_t IsMotherGamma(Int_t tr, AliMCEvent* mcEvent);\r
+ Int_t IsMotherPi0(Int_t tr, AliMCEvent* mcEvent);\r
+ Int_t IsMotherC(Int_t tr, AliMCEvent* mcEvent);\r
+ Int_t IsMotherB(Int_t tr, AliMCEvent* mcEvent);\r
+ Int_t IsMotherEta(Int_t tr, AliMCEvent* mcEvent);\r
+ \r
\r
ClassDef(AliAnalysisTaskHFEFlow, 1); // analysisclass\r
};\r
// Cut event
if(fInputEvent->IsPileupFromSPD(3, 0.8, 3., 2., 5)){
- AliDebug(1, "Event flagged as pileup\n");
+ AliDebug(1, "Event flagged as pileup\n");
return;
}
if(!fTrackCuts->CheckEventCuts("fCutsEvRec", fInputEvent)){
// Get Production Vertex in radial direction
Double_t vx = mcpart->Particle()->Vx(),
- vy = mcpart->Particle()->Vy();
+ vy = mcpart->Particle()->Vy();
Double_t productionVertex = TMath::Sqrt(vx*vx+vy*vy);
// Get Mother PDG code of the particle
- (*fDebugTree) << "MCDebug"
+ (*fDebugTree) << "MCDebug"
<< "centrality=" << centrality
<< "MBtrigger=" << isMBTrigger
<< "CentralTrigger=" << isCentralTrigger
<< "pdg=" << pdg
<< "ProductionVertex=" << productionVertex
<< "motherPdg=" << motherPdg
- << "source=" << source
+ << "source=" << source
<< "\n";
}
}
// Cut track (Only basic track cuts)
if(!fTrackCuts->CheckParticleCuts(AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kStepRecKineITSTPC, track)) continue;
// Debug streaming of PID-related quantities
+ copyTrack.~AliESDtrack();
new(©Track) AliESDtrack(*track);
if(fTPCpid->HasEtaCorrection()) fTPCpid->ApplyEtaCorrection(©Track, AliHFEpidObject::kESDanalysis); // Apply Eta Correction on copy track
Double_t nSigmaTOF = pid->NumberOfSigmasTOF(track, AliPID::kElectron);
mcptTPC = ref->Pt();
}
-
+
TParticle *mctrack1 = mctrack->Particle();
mesonID=GetElecSourceMC(mctrack1);
if(mesonID==AliHFEmcQA::kGammaPi0 || mesonID==AliHFEmcQA::kPi0) mArr=0; //pion
else if(mesonID==AliHFEmcQA::kGammaRho0 || mesonID==AliHFEmcQA::kRho0) mArr=5; //rho
mctrack->XvYvZv(xr);
-
+
eR= TMath::Sqrt(xr[0]*xr[0]+xr[1]*xr[1]);
eZ = xr[2];
TParticle *mctrackt = mctrack->Particle();
AliMCParticle *mctrackmother = NULL;
if(!(mArr<0)){
- if(mesonID>=AliHFEmcQA::kGammaPi0) { // conversion electron, be careful with the enum odering
- Int_t glabel=TMath::Abs(mctrack->GetMother()); // gamma label
-
if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
- glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's label
-
if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
- mesonPt = mctrackmother->Pt(); //meson pt
- bgcategory = 1.;
- mctrackmother->XvYvZv(xr);
- mesonR = TMath::Sqrt(xr[0]*xr[0]+xr[1]*xr[1]);
- mesonZ = xr[2];
-
- mctrackt = mctrackmother->Particle();
- if(mctrackt){
- mesonunique = mctrackt->GetUniqueID();
- }
- if(glabel>fMCEvent->GetNumberOfPrimaries()) {
- bgcategory = 2.;
- glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's mother
-
if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
- mesonMomPdg=mctrackmother->PdgCode();
- mesonMomPt=mctrackmother->Pt();
- if(TMath::Abs(mctrackmother->PdgCode())==310){
- bgcategory = 3.;
- glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's mother's mother
-
if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
- mesonGMomPdg=mctrackmother->PdgCode();
- glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's mother
-
if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
- mesonGGMomPdg=mctrackmother->PdgCode();
- }
- }
- }
- }
- }
- }
- }
- }
- else{ // nonHFE except for the conversion electron
- Int_t glabel=TMath::Abs(mctrack->GetMother());
-
if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
- mesonPt = mctrackmother->Pt(); //meson pt
- bgcategory = -1.;
- mctrackmother->XvYvZv(xr);
- mesonR = TMath::Sqrt(xr[0]*xr[0]+xr[1]*xr[1]);
- mesonZ = xr[2];
-
- mctrackt = mctrackmother->Particle();
- if(mctrackt){
- mesonunique = mctrackt->GetUniqueID();
- }
- if(glabel>fMCEvent->GetNumberOfPrimaries()) {
- bgcategory = -2.;
- glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's mother
-
if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
- mesonMomPdg=mctrackmother->PdgCode();
- mesonMomPt=mctrackmother->Pt();
- if(TMath::Abs(mctrackmother->PdgCode())==310){
- bgcategory = -3.;
- glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's mother's mother
-
if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
+ if(mesonID>=AliHFEmcQA::kGammaPi0) { // conversion electron, be careful with the enum odering
+ Int_t glabel=TMath::Abs(mctrack->GetMother()); // gamma label
+ if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
+ glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's label
+ if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
+ mesonPt = mctrackmother->Pt(); //meson pt
+ bgcategory = 1.;
+ mctrackmother->XvYvZv(xr);
+ mesonR = TMath::Sqrt(xr[0]*xr[0]+xr[1]*xr[1]);
+ mesonZ = xr[2];
+
+ mctrackt = mctrackmother->Particle();
+ if(mctrackt){
+ mesonunique = mctrackt->GetUniqueID();
+ }
+ if(glabel>fMCEvent->GetNumberOfPrimaries()) {
+ bgcategory = 2.;
+ glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's mother
+ if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
+ mesonMomPdg=mctrackmother->PdgCode();
+ mesonMomPt=mctrackmother->Pt();
+ if(TMath::Abs(mctrackmother->PdgCode())==310){
+ bgcategory = 3.;
+ glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's mother's mother
+ if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
+ mesonGMomPdg=mctrackmother->PdgCode();
+ glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's mother
+ if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
+ mesonGGMomPdg=mctrackmother->PdgCode();
+ }
+ }
+ }
+ }
+ }
+ }
+ }
+ }
+ else{ // nonHFE except for the conversion electron
+ Int_t glabel=TMath::Abs(mctrack->GetMother());
+ if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
+ mesonPt = mctrackmother->Pt(); //meson pt
+ bgcategory = -1.;
+ mctrackmother->XvYvZv(xr);
+ mesonR = TMath::Sqrt(xr[0]*xr[0]+xr[1]*xr[1]);
+ mesonZ = xr[2];
+
+ mctrackt = mctrackmother->Particle();
+ if(mctrackt){
+ mesonunique = mctrackt->GetUniqueID();
+ }
+ if(glabel>fMCEvent->GetNumberOfPrimaries()) {
+ bgcategory = -2.;
+ glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's mother
+ if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
+ mesonMomPdg=mctrackmother->PdgCode();
+ mesonMomPt=mctrackmother->Pt();
+ if(TMath::Abs(mctrackmother->PdgCode())==310){
+ bgcategory = -3.;
+ glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's mother's mother
+ if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
mesonGMomPdg=mctrackmother->PdgCode();
- glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's mother
-
if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
- mesonGGMomPdg=mctrackmother->PdgCode();
- }
- }
- }
- }
- }
- }
- }
+ glabel=TMath::Abs(mctrackmother->GetMother()); // gamma's mother's mother
+ if((mctrackmother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(glabel)))){
+ mesonGGMomPdg=mctrackmother->PdgCode();
+ }
+ }
+ }
+ }
+ }
+ }
+ }
}
Double_t trddEdxSum[6];
for(Int_t a=0;a<6;a++) { trddEdxSum[a]= 0.;}
for(Int_t itl = 0; itl < 6; itl++){
- Int_t nSliceNonZero = 0;
+ Int_t nSliceNonZero = 0;
trddEdxSum[itl] = track->GetTRDslice(itl, 0); // in new reconstruction slice 0 contains the total charge
for(Int_t islice = 0; islice < 8; islice++){
if(track->GetTRDslice(itl, islice) > 0.001) nSliceNonZero++;
Double_t hfebCov[3] = {-999.,-999.,-999.};
fExtraCuts->GetHFEImpactParameters(track, hfeb, hfebCov);
+ Double_t tofdx= -999.0;
+ Double_t tofdz= -999.0;
+ tofdx=track->GetTOFsignalDx();
+ tofdz=track->GetTOFsignalDz();
+
// Fill Tree
(*fDebugTree) << "PIDdebug"
<< "centrality=" << centrality
<< "nclustersITS=" << nclustersITS
<< "nclusters=" << nclustersTRD
<< "chi2matching=" << chi2matching
- << "chi2PerClusterITS=" << chi2PerClusterITS
+ << "chi2PerClusterITS=" << chi2PerClusterITS
<< "its0=" << hasClusterITS[0]
<< "its1=" << hasClusterITS[1]
<< "its2=" << hasClusterITS[2]
<< "trd3=" << hasTrackletTRD[3]
<< "trd4=" << hasTrackletTRD[4]
<< "trd5=" << hasTrackletTRD[5]
- << "TRDdEdxl0=" << trddEdxSum[0]
- << "TRDdEdxl1=" << trddEdxSum[1]
- << "TRDdEdxl2=" << trddEdxSum[2]
- << "TRDdEdxl3=" << trddEdxSum[3]
- << "TRDdEdxl4=" << trddEdxSum[4]
- << "TRDdEdxl5=" << trddEdxSum[5]
+ << "TRDdEdxl0=" << trddEdxSum[0]
+ << "TRDdEdxl1=" << trddEdxSum[1]
+ << "TRDdEdxl2=" << trddEdxSum[2]
+ << "TRDdEdxl3=" << trddEdxSum[3]
+ << "TRDdEdxl4=" << trddEdxSum[4]
+ << "TRDdEdxl5=" << trddEdxSum[5]
<< "TOFsigmaEl=" << nSigmaTOF
<< "TPCsigmaEl=" << nSigmaTPC
<< "TPCdEdx=" << tPCdEdx
<< "TRDlikeEl=" << likeEleTRD
<< "TRDlikeEln=" << likeEleTRDn
- << "trdtruncmean1=" << trdtruncmean1
+ << "trdtruncmean1=" << trdtruncmean1
<< "trdtruncmean2=" << trdtruncmean2
<< "dcaR=" << b[0]
<< "dcaZ=" << b[1]
<< "hfedcacovZ=" << hfebCov[2]
<< "vx=" << vtx[0]
<< "vy=" << vtx[1]
- << "vz=" << vtx[2]
- << "ncontrib=" << ncontrib
- << "mesonID=" << mesonID
- << "eR=" << eR
- << "mesonR=" << mesonR
- << "eZ=" << eZ
- << "mesonZ=" << mesonZ
- << "unique=" << unique
- << "mesonunique=" << mesonunique
- << "bgcategory=" << bgcategory
- << "mesonpt=" << mesonPt
+ << "vz=" << vtx[2]
+ << "tofdx=" << tofdx
+ << "tofdz=" << tofdz
+ << "ncontrib=" << ncontrib
+ << "mesonID=" << mesonID
+ << "eR=" << eR
+ << "mesonR=" << mesonR
+ << "eZ=" << eZ
+ << "mesonZ=" << mesonZ
+ << "unique=" << unique
+ << "mesonunique=" << mesonunique
+ << "bgcategory=" << bgcategory
+ << "mesonpt=" << mesonPt
<< "mesonMomPdg=" << mesonMomPdg
<< "mesonGMomPdg=" << mesonGMomPdg
<< "mesonGGMomPdg=" << mesonGGMomPdg
TParticle *partMotherCopy = mctrack->Particle();
Int_t maPdgcode = partMother->GetPdgCode();
- // if the mother is charmed hadron
- if ( (int(abs(maPdgcode)/100.)%10) == AliHFEmcQA::kCharm || (int(abs(maPdgcode)/1000.)%10) == AliHFEmcQA::kCharm ) {
+ // if the mother is charmed hadron
+ if ( (int(abs(maPdgcode)/100.)%10) == AliHFEmcQA::kCharm || (int(abs(maPdgcode)/1000.)%10) == AliHFEmcQA::kCharm ) {
- for (Int_t i=0; i<fNparents; i++){
+ for (Int_t i=0; i<fNparents; i++){
if (abs(maPdgcode)==fParentSelect[0][i]){
isFinalOpenCharm = kTRUE;
}
- }
- if (!isFinalOpenCharm) return -1;
+ }
+ if (!isFinalOpenCharm) return -1;
- // iterate until you find B hadron as a mother or become top ancester
- for (Int_t i=1; i<fgkMaxIter; i++){
+ // iterate until you find B hadron as a mother or become top ancester
+ for (Int_t i=1; i<fgkMaxIter; i++){
Int_t jLabel = partMother->GetFirstMother();
if (jLabel == -1){
Int_t grandMaPDG = grandMa->GetPdgCode();
for (Int_t j=0; j<fNparents; j++){
- if (abs(grandMaPDG)==fParentSelect[1][j]){
- origin = AliHFEmcQA::kBeautyCharm;
- return origin;
- }
+ if (abs(grandMaPDG)==fParentSelect[1][j]){
+ origin = AliHFEmcQA::kBeautyCharm;
+ return origin;
+ }
}
partMother = grandMa;
- } // end of iteration
- } // end of if
- else if ( (int(abs(maPdgcode)/100.)%10) == AliHFEmcQA::kBeauty || (int(abs(maPdgcode)/1000.)%10) == AliHFEmcQA::kBeauty ) {
- for (Int_t i=0; i<fNparents; i++){
+ } // end of iteration
+ } // end of if
+ else if ( (int(abs(maPdgcode)/100.)%10) == AliHFEmcQA::kBeauty || (int(abs(maPdgcode)/1000.)%10) == AliHFEmcQA::kBeauty ) {
+ for (Int_t i=0; i<fNparents; i++){
if (abs(maPdgcode)==fParentSelect[1][i]){
origin = AliHFEmcQA::kDirectBeauty;
return origin;
}
- }
- } // end of if
- else if ( abs(maPdgcode) == 22 ) { //conversion
-
- tmpMomLabel = partMotherCopy->GetFirstMother();
-
if(!(mctrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tmpMomLabel))))) return -1;
- partMother = mctrack->Particle();
- maPdgcode = partMother->GetPdgCode();
- if ( abs(maPdgcode) == 111 ) {
- origin = AliHFEmcQA::kGammaPi0;
- return origin;
- }
- else if ( abs(maPdgcode) == 221 ) {
- origin = AliHFEmcQA::kGammaEta;
- return origin;
- }
- else if ( abs(maPdgcode) == 223 ) {
- origin = AliHFEmcQA::kGammaOmega;
- return origin;
- }
- else if ( abs(maPdgcode) == 333 ) {
- origin = AliHFEmcQA::kGammaPhi;
- return origin;
- }
- else if ( abs(maPdgcode) == 331 ) {
- origin = AliHFEmcQA::kGammaEtaPrime;
- return origin;
- }
- else if ( abs(maPdgcode) == 113 ) {
- origin = AliHFEmcQA::kGammaRho0;
- return origin;
- }
- else origin = AliHFEmcQA::kElse;
- //origin = kGamma; // finer category above
- return origin;
-
- } // end of if
- else if ( abs(maPdgcode) == 111 ) {
- origin = AliHFEmcQA::kPi0;
- return origin;
- } // end of if
- else if ( abs(maPdgcode) == 221 ) {
- origin = AliHFEmcQA::kEta;
- return origin;
- } // end of if
- else if ( abs(maPdgcode) == 223 ) {
- origin = AliHFEmcQA::kOmega;
- return origin;
- } // end of if
- else if ( abs(maPdgcode) == 333 ) {
- origin = AliHFEmcQA::kPhi;
- return origin;
- } // end of if
- else if ( abs(maPdgcode) == 331 ) {
- origin = AliHFEmcQA::kEtaPrime;
- return origin;
- } // end of if
- else if ( abs(maPdgcode) == 113 ) {
- origin = AliHFEmcQA::kRho0;
- return origin;
- } // end of if
- else{
+ }
+ } // end of if
+ else if ( abs(maPdgcode) == 22 ) { //conversion
+
+ tmpMomLabel = partMotherCopy->GetFirstMother();
+ if(!(mctrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tmpMomLabel))))) return -1;
+ partMother = mctrack->Particle();
+ maPdgcode = partMother->GetPdgCode();
+ if ( abs(maPdgcode) == 111 ) {
+ origin = AliHFEmcQA::kGammaPi0;
+ return origin;
+ }
+ else if ( abs(maPdgcode) == 221 ) {
+ origin = AliHFEmcQA::kGammaEta;
+ return origin;
+ }
+ else if ( abs(maPdgcode) == 223 ) {
+ origin = AliHFEmcQA::kGammaOmega;
+ return origin;
+ }
+ else if ( abs(maPdgcode) == 333 ) {
+ origin = AliHFEmcQA::kGammaPhi;
+ return origin;
+ }
+ else if ( abs(maPdgcode) == 331 ) {
+ origin = AliHFEmcQA::kGammaEtaPrime;
+ return origin;
+ }
+ else if ( abs(maPdgcode) == 113 ) {
+ origin = AliHFEmcQA::kGammaRho0;
+ return origin;
+ }
+ else origin = AliHFEmcQA::kElse;
+ //origin = kGamma; // finer category above
+ return origin;
+
+ } // end of if
+ else if ( abs(maPdgcode) == 111 ) {
+ origin = AliHFEmcQA::kPi0;
+ return origin;
+ } // end of if
+ else if ( abs(maPdgcode) == 221 ) {
+ origin = AliHFEmcQA::kEta;
+ return origin;
+ } // end of if
+ else if ( abs(maPdgcode) == 223 ) {
+ origin = AliHFEmcQA::kOmega;
+ return origin;
+ } // end of if
+ else if ( abs(maPdgcode) == 333 ) {
+ origin = AliHFEmcQA::kPhi;
+ return origin;
+ } // end of if
+ else if ( abs(maPdgcode) == 331 ) {
+ origin = AliHFEmcQA::kEtaPrime;
+ return origin;
+ } // end of if
+ else if ( abs(maPdgcode) == 113 ) {
+ origin = AliHFEmcQA::kRho0;
+ return origin;
+ } // end of if
+ else{
origin = AliHFEmcQA::kElse;
- }
- return origin;
+ }
+ return origin;
}
//______________________________________________________
Bool_t AliHFEextraCuts::IsKinkMother(AliVTrack *track){
//
- // Is kink Mother
+ // Is kink Mother: only for ESD since need to loop over vertices for AOD
//
+ //
+
TClass *type = track->IsA();
if(type == AliESDtrack::Class()){
AliESDtrack *esdtrack = static_cast<AliESDtrack *>(track);
if(esdtrack->GetKinkIndex(0)!=0) return kTRUE;
else return kFALSE;
}
- else if(type == AliAODTrack::Class()){
- AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
- if(aodtrack){
- //printf("find a track\n");
- AliAODVertex *aodvertex = aodtrack->GetProdVertex();
- if(!aodvertex) return kFALSE;
- Int_t n = aodvertex->GetNDaughters();
- //printf("Number of daughters %d\n",n);
- for(Int_t k=0; k < n; k++) {
- AliAODTrack *aodtrackdaughter = (AliAODTrack *) aodvertex->GetDaughter(k);
- if(aodtrackdaughter){
- AliAODVertex *aodvertexdaughter = aodtrackdaughter->GetProdVertex();
- if(!aodvertexdaughter) continue;
- //if(aodvertexdaughter->GetType()==AliAODVertex::kKink) printf("Daughter of type %d\n",(Int_t)aodvertexdaughter->GetType());
- if(aodvertexdaughter->GetType()==AliAODVertex::kKink) return kTRUE;
- }
- }
- return kFALSE;
- }
- else return kFALSE;
- }
return kFALSE;
fNumberOfIterations(5),
fChargeChoosen(kAllCharge),
fNCentralityBinAtTheEnd(0),
+ fTestCentralityLow(-1),
+ fTestCentralityHigh(-1),
fHadronEffbyIPcut(NULL),
fConversionEffbgc(NULL),
fNonHFEEffbgc(NULL),
fEfficiencyesdTOFPIDD[k] = 0;
fEfficiencyIPCharmD[k] = 0;
fEfficiencyIPBeautyD[k] = 0;
+ fEfficiencyIPBeautyesdD[k] = 0;
fEfficiencyIPConversionD[k] = 0;
fEfficiencyIPNonhfeD[k] = 0;
fBeautyEff[k] = 0;
fEfficiencyCharmSigD[k] = 0;
fEfficiencyBeautySigD[k] = 0;
+ fEfficiencyBeautySigesdD[k] = 0;
}
}
memset(fEtaRange, 0, sizeof(Double_t) * 2);
AliCFContainer *contaminationcontainer = datahfecontainer->GetCFContainer("hadronicBackground");
if((!datacontainer) || (!contaminationcontainer)) return kFALSE;
- AliCFContainer *datacontainerD = GetSlicedContainer(datacontainer, fNbDimensions, dims, -1, fChargeChoosen);
- AliCFContainer *contaminationcontainerD = GetSlicedContainer(contaminationcontainer, fNbDimensions, dims, -1, fChargeChoosen);
+ AliCFContainer *datacontainerD = GetSlicedContainer(datacontainer, fNbDimensions, dims, -1, fChargeChoosen,fTestCentralityLow,fTestCentralityHigh);
+ AliCFContainer *contaminationcontainerD = GetSlicedContainer(contaminationcontainer, fNbDimensions, dims, -1, fChargeChoosen,fTestCentralityLow,fTestCentralityHigh);
if((!datacontainerD) || (!contaminationcontainerD)) return kFALSE;
SetContainer(datacontainerD,AliHFEspectrum::kDataContainer);
if(fBeamType==1)
{
- fConvSourceContainer[iSource][iLevel][icentr] = GetSlicedContainer(convtempContainer, fNbDimensions, dims, -1, fChargeChoosen,icentr+1);
- fNonHFESourceContainer[iSource][iLevel][icentr] = GetSlicedContainer(nonHFEtempContainer, fNbDimensions, dims, -1, fChargeChoosen,icentr+1);
+ fConvSourceContainer[iSource][iLevel][icentr] = GetSlicedContainer(convtempContainer, fNbDimensions, dims, -1, fChargeChoosen,icentr+1,icentr+1);
+ fNonHFESourceContainer[iSource][iLevel][icentr] = GetSlicedContainer(nonHFEtempContainer, fNbDimensions, dims, -1, fChargeChoosen,icentr+1,icentr+1);
}
// if((!fConvSourceContainer[iSource][iLevel][icentr])||(!fNonHFESourceContainer[iSource][iLevel])) return kFALSE;
}
Int_t source = -1;
if(!fInclusiveSpectrum) source = 1; //beauty
- AliCFContainer *mccontaineresdD = GetSlicedContainer(mccontaineresd, fNbDimensions, dims, source, fChargeChoosen);
- AliCFContainer *mccontainermcD = GetSlicedContainer(mccontainermc, fNbDimensions, dims, source, fChargeChoosen);
+ AliCFContainer *mccontaineresdD = GetSlicedContainer(mccontaineresd, fNbDimensions, dims, source, fChargeChoosen,fTestCentralityLow,fTestCentralityHigh);
+ AliCFContainer *mccontainermcD = GetSlicedContainer(mccontainermc, fNbDimensions, dims, source, fChargeChoosen,fTestCentralityLow,fTestCentralityHigh);
if((!mccontaineresdD) || (!mccontainermcD)) return kFALSE;
SetContainer(mccontainermcD,AliHFEspectrum::kMCContainerMC);
SetContainer(mccontaineresdD,AliHFEspectrum::kMCContainerESD);
mccontainermcD = GetSlicedContainer(mccontainermcbg, fNbDimensions, dims, source, fChargeChoosen);
SetContainer(mccontainermcD,AliHFEspectrum::kMCContainerCharmMC);
- SetParameterizedEff(mccontainermc, mccontainermcbg, mccontaineresd, mccontaineresdbg, dims);
-
if(!fNonHFEsyst){
AliCFContainer *nonHFEweightedContainerD = GetSlicedContainer(nonHFEweightedContainer, fNbDimensions, dims, -1, fChargeChoosen);
SetContainer(nonHFEweightedContainerD,AliHFEspectrum::kMCWeightedContainerNonHFEESD);
AliCFContainer *convweightedContainerD = GetSlicedContainer(convweightedContainer, fNbDimensions, dims, -1, fChargeChoosen);
SetContainer(convweightedContainerD,AliHFEspectrum::kMCWeightedContainerConversionESD);
+ }
- if(fBeamType==0){
- AliCFEffGrid *nonHFEEffGrid = (AliCFEffGrid*) GetEfficiency(GetContainer(kMCWeightedContainerNonHFEESD),1,0);
- fNonHFEEffbgc = (TH1D *) nonHFEEffGrid->Project(0);
-
- AliCFEffGrid *conversionEffGrid = (AliCFEffGrid*) GetEfficiency(GetContainer(kMCWeightedContainerConversionESD),1,0);
- fConversionEffbgc = (TH1D *) conversionEffGrid->Project(0);
- }
+ SetParameterizedEff(mccontainermc, mccontainermcbg, mccontaineresd, mccontaineresdbg, dims);
- }
}
// MC container: correlation matrix
THnSparseF *mccorrelation = 0x0;
else mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepafterPID"); // we confirmed that we get same result by using it instead of correlationstepafterDE
//else mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepafterDE");
if(!mccorrelation) return kFALSE;
- THnSparseF *mccorrelationD = GetSlicedCorrelation(mccorrelation, fNbDimensions, dims);
+ THnSparseF *mccorrelationD = GetSlicedCorrelation(mccorrelation, fNbDimensions, dims,fTestCentralityLow,fTestCentralityHigh);
if(!mccorrelationD) {
printf("No correlation\n");
return kFALSE;
if(v0hfecontainer) {
AliCFContainer *containerV0 = v0hfecontainer->GetCFContainer("taggedTrackContainerReco");
if(!containerV0) return kFALSE;
- AliCFContainer *containerV0Electron = GetSlicedContainer(containerV0, fNbDimensions, dims, AliPID::kElectron);
+ AliCFContainer *containerV0Electron = GetSlicedContainer(containerV0, fNbDimensions, dims, AliPID::kElectron,fChargeChoosen,fTestCentralityLow,fTestCentralityHigh);
if(!containerV0Electron) return kFALSE;
SetContainer(containerV0Electron,AliHFEspectrum::kDataContainerV0);
}
if(fWriteToFile) ccontaminationspectrum->SaveAs("contaminationspectrum.eps");
}
+ TFile *file = TFile::Open("tests.root","recreate");
+ datacontainerD->Write("data");
+ mccontainermcD->Write("mcefficiency");
+ mccorrelationD->Write("correlationmatrix");
+ file->Close();
return kTRUE;
}
{
correctedspectrum->GetAxis(0)->SetRange(i+1,i+1);
correctedspectrumDc[i] = Normalize(correctedspectrum,i);
- correctedspectrumDc[i]->SetTitle(Form("UnfoldingCorrectedSpectrum_%i",i));
- correctedspectrumDc[i]->SetName(Form("UnfoldingCorrectedSpectrum_%i",i));
+ if(correctedspectrumDc[i]){
+ correctedspectrumDc[i]->SetTitle(Form("UnfoldingCorrectedSpectrum_%i",i));
+ correctedspectrumDc[i]->SetName(Form("UnfoldingCorrectedSpectrum_%i",i));
+ correctedspectrumDc[i]->Write();
+ }
alltogetherCorrection->GetAxis(0)->SetRange(i+1,i+1);
alltogetherspectrumDc[i] = Normalize(alltogetherCorrection,i);
- alltogetherspectrumDc[i]->SetTitle(Form("AlltogetherSpectrum_%i",i));
- alltogetherspectrumDc[i]->SetName(Form("AlltogetherSpectrum_%i",i));
- correctedspectrumDc[i]->Write();
- alltogetherspectrumDc[i]->Write();
-
-
+ if(alltogetherspectrumDc[i]){
+ alltogetherspectrumDc[i]->SetTitle(Form("AlltogetherSpectrum_%i",i));
+ alltogetherspectrumDc[i]->SetName(Form("AlltogetherSpectrum_%i",i));
+ alltogetherspectrumDc[i]->Write();
+ }
+
TH1D *centrcrosscheck = correctedspectrum->Projection(0);
centrcrosscheck->SetTitle(Form("centrality_%i",i));
centrcrosscheck->SetName(Form("centrality_%i",i));
AliCFDataGrid *dataspectrumbeforesubstraction = (AliCFDataGrid *) ((AliCFDataGrid *)GetSpectrum(GetContainer(kDataContainer),fStepData))->Clone();
dataspectrumbeforesubstraction->SetName("dataspectrumbeforesubstraction");
+
// Background Estimate
AliCFContainer *backgroundContainer = GetContainer(kBackgroundData);
if(!backgroundContainer){
if(!fInclusiveSpectrum){
//Background subtraction for IP analysis
+
+ TH1D *incElecCent[kCentrality-1];
+ TH1D *charmCent[kCentrality-1];
+ TH1D *convCent[kCentrality-1];
+ TH1D *nonHFECent[kCentrality-1];
+ TH1D *subtractedCent[kCentrality-1];
TH1D *measuredTH1Draw = (TH1D *) dataspectrumbeforesubstraction->Project(ptpr);
CorrectFromTheWidth(measuredTH1Draw);
- TCanvas *rawspectra = new TCanvas("rawspectra","rawspectra",600,500);
+ if(fBeamType==1){
+ THnSparseF* sparseIncElec = (THnSparseF *) dataspectrumbeforesubstraction->GetGrid();
+ for(Int_t icent = 1; icent < kCentrality-1; icent++){
+ sparseIncElec->GetAxis(0)->SetRange(icent,icent);
+ incElecCent[icent-1] = (TH1D *) sparseIncElec->Projection(ptpr);
+ CorrectFromTheWidth(incElecCent[icent-1]);
+ }
+ }
+ TCanvas *rawspectra = new TCanvas("rawspectra","rawspectra",500,400);
rawspectra->cd();
rawspectra->SetLogy();
- TLegend *lRaw = new TLegend(0.65,0.65,0.95,0.95);
+ gStyle->SetOptStat(0);
+ TLegend *lRaw = new TLegend(0.55,0.55,0.85,0.85);
measuredTH1Draw->SetMarkerStyle(20);
measuredTH1Draw->Draw();
+ measuredTH1Draw->GetXaxis()->SetRangeUser(0.0,7.9);
lRaw->AddEntry(measuredTH1Draw,"measured raw spectrum");
TH1D* htemp;
Int_t* bins=new Int_t[2];
lRaw->AddEntry(charmbg,"charm elecs");
// subtract charm background
spectrumSubtracted->Add(charmbgContainer,-1.0);
+ if(fBeamType==1){
+ THnSparseF* sparseCharmElec = (THnSparseF *) charmbgContainer->GetGrid();
+ for(Int_t icent = 1; icent < kCentrality-1; icent++){
+ sparseCharmElec->GetAxis(0)->SetRange(icent,icent);
+ charmCent[icent-1] = (TH1D *) sparseCharmElec->Projection(ptpr);
+ CorrectFromTheWidth(charmCent[icent-1]);
+ }
+ }
}
if(fIPanaConversionBgSubtract){
// Conversion background
lRaw->AddEntry(conversionbg,"conversion elecs");
// subtract conversion background
spectrumSubtracted->Add(conversionbgContainer,-1.0);
- printf("Conversion background subtraction is preliminary!\n");
+ if(fBeamType==1){
+ THnSparseF* sparseconvElec = (THnSparseF *) conversionbgContainer->GetGrid();
+ for(Int_t icent = 1; icent < kCentrality-1; icent++){
+ sparseconvElec->GetAxis(0)->SetRange(icent,icent);
+ convCent[icent-1] = (TH1D *) sparseconvElec->Projection(ptpr);
+ CorrectFromTheWidth(convCent[icent-1]);
+ }
+ }
}
if(fIPanaNonHFEBgSubtract){
// NonHFE background
lRaw->AddEntry(nonhfebg,"non-HF elecs");
// subtract Dalitz/dielectron background
spectrumSubtracted->Add(nonHFEbgContainer,-1.0);
- printf("Non HFE background subtraction is preliminary!\n");
+ if(fBeamType==1){
+ THnSparseF* sparseNonHFEElec = (THnSparseF *) nonHFEbgContainer->GetGrid();
+ for(Int_t icent = 1; icent < kCentrality-1; icent++){
+ sparseNonHFEElec->GetAxis(0)->SetRange(icent,icent);
+ nonHFECent[icent-1] = (TH1D *) sparseNonHFEElec->Projection(ptpr);
+ CorrectFromTheWidth(nonHFECent[icent-1]);
+ }
+ }
}
+
TH1D *rawbgsubtracted = (TH1D *) spectrumSubtracted->Project(ptpr);
CorrectFromTheWidth(rawbgsubtracted);
rawbgsubtracted->SetMarkerStyle(24);
lRaw->AddEntry(rawbgsubtracted,"subtracted raw spectrum");
rawbgsubtracted->Draw("samep");
lRaw->Draw("SAME");
+ gPad->SetGrid();
+ //rawspectra->SaveAs("rawspectra.eps");
+
+ if(fBeamType==1){
+ THnSparseF* sparseSubtracted = (THnSparseF *) spectrumSubtracted->GetGrid();
+ for(Int_t icent = 1; icent < kCentrality-1; icent++){
+ sparseSubtracted->GetAxis(0)->SetRange(icent,icent);
+ subtractedCent[icent-1] = (TH1D *) sparseSubtracted->Projection(ptpr);
+ CorrectFromTheWidth(subtractedCent[icent-1]);
+ }
+
+ TLegend *lCentRaw = new TLegend(0.55,0.55,0.85,0.85);
+ TCanvas *centRaw = new TCanvas("centRaw","centRaw",1000,800);
+ centRaw->Divide(3,3);
+ for(Int_t icent = 1; icent < kCentrality-1; icent++){
+ centRaw->cd(icent);
+ gPad->SetLogx();
+ gPad->SetLogy();
+ incElecCent[icent-1]->GetXaxis()->SetRangeUser(0.4,8.);
+ incElecCent[icent-1]->Draw("p");
+ incElecCent[icent-1]->SetMarkerColor(1);
+ incElecCent[icent-1]->SetMarkerStyle(20);
+ charmCent[icent-1]->Draw("samep");
+ charmCent[icent-1]->SetMarkerColor(3);
+ charmCent[icent-1]->SetMarkerStyle(20);
+ convCent[icent-1]->Draw("samep");
+ convCent[icent-1]->SetMarkerColor(4);
+ convCent[icent-1]->SetMarkerStyle(20);
+ nonHFECent[icent-1]->Draw("samep");
+ nonHFECent[icent-1]->SetMarkerColor(6);
+ nonHFECent[icent-1]->SetMarkerStyle(20);
+ subtractedCent[icent-1]->Draw("samep");
+ subtractedCent[icent-1]->SetMarkerStyle(24);
+ if(icent == 1){
+ lCentRaw->AddEntry(incElecCent[0],"inclusive electron spectrum");
+ lCentRaw->AddEntry(charmCent[0],"charm elecs");
+ lCentRaw->AddEntry(convCent[0],"conversion elecs");
+ lCentRaw->AddEntry(nonHFECent[0],"non-HF elecs");
+ lCentRaw->AddEntry(subtractedCent[0],"subtracted electron spectrum");
+ lCentRaw->Draw("SAME");
+ }
+ }
+ }
delete[] bins;
if(fBeamType==1)
{
- //charmbgaftertofpid = (TH1D *) charmBackgroundGrid->Project(0);
bins[0]=12;
charmbgaftertofpid = (TH1D *) charmBackgroundGrid->Project(1);
bins[1]=charmbgaftertofpid->GetNbinsX();
ipWeightedCharmContainer->SetGrid(GetPIDxIPEff(0)); // get charm efficiency
TH1D* parametrizedcharmpidipeff = (TH1D*)ipWeightedCharmContainer->Project(ptpr);
- if(fBeamType==0)charmBackgroundGrid->Multiply(ipWeightedCharmContainer,evtnorm);
- Double_t contents[2];
- AliCFDataGrid *eventTemp = new AliCFDataGrid("eventTemp","eventTemp",nDim,bins);
- if(fBeamType==1)
- {
- for(Int_t kCentr=0;kCentr<bins[0];kCentr++)
+ charmBackgroundGrid->Multiply(ipWeightedCharmContainer,1.);
+
+ Int_t *nBinpp=new Int_t[1];
+ Int_t* binspp=new Int_t[1];
+ binspp[0]=charmbgaftertofpid->GetNbinsX(); // number of pt bins
+
+ Int_t *nBinPbPb=new Int_t[2];
+ Int_t* binsPbPb=new Int_t[2];
+ binsPbPb[1]=charmbgaftertofpid->GetNbinsX(); // number of pt bins
+ binsPbPb[0]=12;
+
+ Int_t looppt=binspp[0];
+ if(fBeamType==1) looppt=binsPbPb[1];
+
+ for(Long_t iBin=1; iBin<= looppt;iBin++){
+ if(fBeamType==0)
{
- for(Int_t kpt=0;kpt<bins[1];kpt++)
- {
- Double_t evtnormPbPb=0;
- if(fNMCbgEvents[kCentr]) evtnormPbPb= double(fNEvents[kCentr])/double(fNMCbgEvents[kCentr]);
- contents[0]=kCentr;
- contents[1]=kpt;
- eventTemp->Fill(contents,evtnormPbPb);
- }
+ nBinpp[0]=iBin;
+ charmBackgroundGrid->SetElementError(nBinpp, charmBackgroundGrid->GetElementError(nBinpp)*evtnorm);
+ charmBackgroundGrid->SetElement(nBinpp,charmBackgroundGrid->GetElement(nBinpp)*evtnorm);
+ }
+ if(fBeamType==1)
+ {
+ // loop over centrality
+ for(Long_t iiBin=1; iiBin<= binsPbPb[0];iiBin++){
+ nBinPbPb[0]=iiBin;
+ nBinPbPb[1]=iBin;
+ Double_t evtnormPbPb=0;
+ if(fNMCbgEvents[iiBin-1]) evtnormPbPb= double(fNEvents[iiBin-1])/double(fNMCbgEvents[iiBin-1]);
+ charmBackgroundGrid->SetElementError(nBinPbPb, charmBackgroundGrid->GetElementError(nBinPbPb)*evtnormPbPb);
+ charmBackgroundGrid->SetElement(nBinPbPb,charmBackgroundGrid->GetElement(nBinPbPb)*evtnormPbPb);
+ }
}
- charmBackgroundGrid->Multiply(eventTemp,1);
}
+
TH1D* charmbgafteripcut = (TH1D*)charmBackgroundGrid->Project(ptpr);
AliCFDataGrid *weightedCharmContainer = new AliCFDataGrid("weightedCharmContainer","weightedCharmContainer",nDim,bins);
}
Int_t stepbackground = 3;
- if(TMath::Abs(fEtaRange[0]) < 0.5) stepbackground = 4;
AliCFDataGrid *backgroundGrid = new AliCFDataGrid("ConversionBgGrid","ConversionBgGrid",*backgroundContainer,stepbackground);
Int_t *nBinpp=new Int_t[1];
Int_t stepbackground = 3;
- if(TMath::Abs(fEtaRange[0]) < 0.5) stepbackground = 4;
AliCFDataGrid *backgroundGrid = new AliCFDataGrid("NonHFEBgGrid","NonHFEBgGrid",*backgroundContainer,stepbackground);
Int_t *nBinpp=new Int_t[1];
bins[0]=35;
AliCFEffGrid *beffContainer = new AliCFEffGrid("beffContainer","beffContainer",1,bins);
- beffContainer->SetGrid(GetBeautyIPEff());
+ beffContainer->SetGrid(GetBeautyIPEff(kTRUE));
efficiencyD->Multiply(beffContainer,1);
}
}
// Unfold
AliCFUnfolding unfolding("unfolding","",fNbDimensions,fCorrelation,efficiencyD->GetGrid(),dataGrid->GetGrid(),guessedTHnSparse);
- //unfolding.SetUseCorrelatedErrors();
if(fUnSetCorrelatedErrors) unfolding.UnsetCorrelatedErrors();
unfolding.SetMaxNumberOfIterations(fNumberOfIterations);
if(fSetSmoothing) unfolding.UseSmoothing();
bins[0]=35;
AliCFEffGrid *beffContainer = new AliCFEffGrid("beffContainer","beffContainer",1,bins);
- beffContainer->SetGrid(GetBeautyIPEff());
+ beffContainer->SetGrid(GetBeautyIPEff(kFALSE));
efficiencyD->Multiply(beffContainer,1);
}
}
return dynamic_cast<AliCFContainer *>(fCFContainers->At(type));
}
//____________________________________________________________
-AliCFContainer *AliHFEspectrum::GetSlicedContainer(AliCFContainer *container, Int_t nDim, Int_t *dimensions,Int_t source,Chargetype_t charge, Int_t centrality) {
+AliCFContainer *AliHFEspectrum::GetSlicedContainer(AliCFContainer *container, Int_t nDim, Int_t *dimensions,Int_t source,Chargetype_t charge, Int_t centralitylow, Int_t centralityhigh) {
//
// Slice bin for a given source of electron
// nDim is the number of dimension the corrections are done
AliInfo(Form("Normalization done in eta range [%f,%f]\n", fEtaRangeNorm[0], fEtaRangeNorm[0]));
}
if(ivar == 5){
- if((centrality>= 0) && (centrality<container->GetNBins(ivar))) {
- varMin[ivar] = binLimits[centrality];
- varMax[ivar] = binLimits[centrality];
+ if((centralitylow>= 0) && (centralitylow<container->GetNBins(ivar)) && (centralityhigh>= 0) && (centralityhigh<container->GetNBins(ivar))) {
+ varMin[ivar] = binLimits[centralitylow];
+ varMax[ivar] = binLimits[centralityhigh];
+
+ TAxis *axistest = container->GetAxis(5,0);
+ printf("GetSlicedContainer: Number of bin in centrality direction %d\n",axistest->GetNbins());
+ printf("Project from %f to %f\n",binLimits[centralitylow],binLimits[centralityhigh]);
+ Double_t lowcentrality = axistest->GetBinLowEdge(axistest->FindBin(binLimits[centralitylow]));
+ Double_t highcentrality = axistest->GetBinUpEdge(axistest->FindBin(binLimits[centralityhigh]));
+ printf("GetSlicedContainer: Low centrality %f and high centrality %f\n",lowcentrality,highcentrality);
+
}
+
}
}
//_________________________________________________________________________
-THnSparseF *AliHFEspectrum::GetSlicedCorrelation(THnSparseF *correlationmatrix, Int_t nDim, Int_t *dimensions) const {
+THnSparseF *AliHFEspectrum::GetSlicedCorrelation(THnSparseF *correlationmatrix, Int_t nDim, Int_t *dimensions, Int_t centralitylow, Int_t centralityhigh) const {
//
// Slice correlation
//
AliError("Problem in the dimensions");
return NULL;
}
+
+ // Cut in centrality is centrality > -1
+ if((centralitylow >=0) && (centralityhigh >=0)) {
+
+ TAxis *axiscentrality0 = correlationmatrix->GetAxis(5);
+ TAxis *axiscentrality1 = correlationmatrix->GetAxis(5+((Int_t)(ndimensions/2.)));
+
+ Int_t bins0 = axiscentrality0->GetNbins();
+ Int_t bins1 = axiscentrality1->GetNbins();
+
+ printf("Number of centrality bins: %d and %d\n",bins0,bins1);
+ if(bins0 != bins1) {
+ AliError("Problem in the dimensions");
+ return NULL;
+ }
+
+ if((centralitylow>= 0) && (centralitylow<bins0) && (centralityhigh>= 0) && (centralityhigh<bins0)) {
+ axiscentrality0->SetRangeUser(centralitylow,centralityhigh);
+ axiscentrality1->SetRangeUser(centralitylow,centralityhigh);
+
+ Double_t lowcentrality0 = axiscentrality0->GetBinLowEdge(axiscentrality0->FindBin(centralitylow));
+ Double_t highcentrality0 = axiscentrality0->GetBinUpEdge(axiscentrality0->FindBin(centralityhigh));
+ Double_t lowcentrality1 = axiscentrality1->GetBinLowEdge(axiscentrality1->FindBin(centralitylow));
+ Double_t highcentrality1 = axiscentrality1->GetBinUpEdge(axiscentrality1->FindBin(centralityhigh));
+ printf("GetCorrelation0: Low centrality %f and high centrality %f\n",lowcentrality0,highcentrality0);
+ printf("GetCorrelation1: Low centrality %f and high centrality %f\n",lowcentrality1,highcentrality1);
+
+ TCanvas * ctestcorrelation = new TCanvas("testcorrelationprojection","testcorrelationprojection",1000,700);
+ ctestcorrelation->cd(1);
+ TH2D* projection = (TH2D *) correlationmatrix->Projection(5,5+((Int_t)(ndimensions/2.)));
+ projection->Draw("colz");
+
+ }
+
+ }
+
+
Int_t ndimensionsContainer = (Int_t) ndimensions/2;
//printf("Number of dimension %d container\n",ndimensionsContainer);
loopcentr=nBinPbPb[0];
}
+ // Weighting factor for pp
Double_t weight[35]={0.859260, 0.872552, 0.847475, 0.823631, 0.839386, 0.874024, 0.916755, 0.942801, 0.965856, 0.933905, 0.933414, 0.931936, 0.847826, 0.810902, 0.796608, 0.727002, 0.659227, 0.583610, 0.549956, 0.512633, 0.472254, 0.412364, 0.353191, 0.319145, 0.305412, 0.290334, 0.269863, 0.254646, 0.230245, 0.200859, 0.275953, 0.276271, 0.227332, 0.197004, 0.474385};
- // Weighting factor for PbPb
- //Double_t weight[35]={0.645016, 0.643397, 0.617149, 0.641176, 0.752285, 0.838097, 0.996226, 1.152757, 1.243257, 1.441421, 1.526796, 1.755632, 1.731234, 1.900269, 1.859628, 1.945138, 1.943707, 1.739451, 1.640120, 1.318674, 1.041654, 0.826493, 0.704605, 0.662040, 0.572361, 0.644030, 0.479850, 0.559513, 0.504044, 0.449495, 0.227605, 0.186836, 0.038681, 0.000000, 0.000000};
+ // Weighting factor for PbPb (0-20%)
+ //Double_t weight[35]={0.641897, 0.640472, 0.615228, 0.650469, 0.737762, 0.847867, 1.009317, 1.158594, 1.307482, 1.476973, 1.551131, 1.677131, 1.785478, 1.888933, 2.017957, 2.074757, 1.926700, 1.869495, 1.546558, 1.222873, 1.160313, 0.903375, 0.799642, 0.706244, 0.705449, 0.599947, 0.719570, 0.499422, 0.703978, 0.477452, 0.325057, 0.093391, 0.096675, 0.000000, 0.000000};
+
+ // Weighting factor for PbPb (40-80%)
+ //Double_t weight[35]={0.181953, 0.173248, 0.166799, 0.182558, 0.206581, 0.236955, 0.279390, 0.329129, 0.365260, 0.423059, 0.452057, 0.482726, 0.462627, 0.537770, 0.584663, 0.579452, 0.587194, 0.499498, 0.443299, 0.398596, 0.376695, 0.322331, 0.260890, 0.374834, 0.249114, 0.310330, 0.287326, 0.243174, 0.758945, 0.138867, 0.170576, 0.107797, 0.011390, 0.000000, 0.000000};
//points
Double_t pt[1];
TF1 *fittofpid = new TF1("fittofpid","[0]*([1]+[2]*log(x)+[3]*log(x)*log(x))*tanh([4]*x-[5])",0.5,8.);
TF1 *fipfit = new TF1("fipfit","[0]*([1]+[2]*log(x)+[3]*log(x)*log(x))*tanh([4]*x-[5])",0.5,8.);
- TF1 *fipfit2 = new TF1("fipfit2","[0]*([1]+[2]*log(x)+[3]*log(x)*log(x))*tanh([4]*x-[5])",0.5,10.0);
+ TF1 *fipfitnonhfe = new TF1("fipfitnonhfe","[0]*([1]+[2]*log(x)+[3]*log(x)*log(x))*tanh([4]*x-[5])",0.3,10.0);
- TCanvas * cefficiencyParam = new TCanvas("efficiencyParam","efficiencyParam",1200,600);
- cefficiencyParam->Divide(2,1);
- cefficiencyParam->cd(1);
+ TCanvas * cefficiencyParamtof = new TCanvas("efficiencyParamtof","efficiencyParamtof",600,600);
+ cefficiencyParamtof->cd();
AliCFContainer *mccontainermcD = 0x0;
+ AliCFContainer *mccontaineresdD = 0x0;
TH1D* efficiencysigTOFPIDD[nCentralitybinning];
TH1D* efficiencyTOFPIDD[nCentralitybinning];
TH1D* efficiencysigesdTOFPIDD[nCentralitybinning];
legtofeff->Draw("same");
- cefficiencyParam->cd(2);
+ TCanvas * cefficiencyParamIP = new TCanvas("efficiencyParamIP","efficiencyParamIP",500,500);
+ cefficiencyParamIP->cd();
+ gStyle->SetOptStat(0);
+
// IP cut efficiencies
-
for(int icentr=0; icentr<loopcentr; icentr++) {
AliCFContainer *charmCombined = 0x0;
AliCFContainer *beautyCombined = 0x0;
+ AliCFContainer *beautyCombinedesd = 0x0;
printf("centrality printing %i \n",icentr);
beautyCombined = (AliCFContainer*)mccontainermcD->Clone("beautyCombined");
+ if(fBeamType==0) mccontaineresdD = GetSlicedContainer(containeresd, fNbDimensions, dimensions, source, fChargeChoosen);
+ else mccontaineresdD = GetSlicedContainer(containeresd, fNbDimensions, dimensions, source, fChargeChoosen, icentr+1);
+ AliCFEffGrid* efficiencyBeautySigesd = new AliCFEffGrid("efficiencyBeautySigesd","",*mccontaineresdD);
+ efficiencyBeautySigesd->CalculateEfficiency(AliHFEcuts::kNcutStepsMCTrack + fStepData,AliHFEcuts::kNcutStepsMCTrack + fStepData-1); // ip cut efficiency.
+
+ beautyCombinedesd = (AliCFContainer*)mccontaineresdD->Clone("beautyCombinedesd");
+
source = 0; //charm mb
if(fBeamType==0) mccontainermcD = GetSlicedContainer(containermb, fNbDimensions, dimensions, source, fChargeChoosen);
else mccontainermcD = GetSlicedContainer(containermb, fNbDimensions, dimensions, source, fChargeChoosen, icentr+1);
AliCFEffGrid* efficiencyBeautyCombined = new AliCFEffGrid("efficiencyBeautyCombined","",*beautyCombined);
efficiencyBeautyCombined->CalculateEfficiency(AliHFEcuts::kNcutStepsMCTrack + fStepData,AliHFEcuts::kNcutStepsMCTrack + fStepData-1);
+ if(fBeamType==0) mccontaineresdD = GetSlicedContainer(containeresdmb, fNbDimensions, dimensions, source, fChargeChoosen);
+ else mccontaineresdD = GetSlicedContainer(containeresdmb, fNbDimensions, dimensions, source, fChargeChoosen, icentr+1);
+ AliCFEffGrid* efficiencyBeautyesd = new AliCFEffGrid("efficiencyBeautyesd","",*mccontaineresdD);
+ efficiencyBeautyesd->CalculateEfficiency(AliHFEcuts::kNcutStepsMCTrack + fStepData,AliHFEcuts::kNcutStepsMCTrack + fStepData-1); // ip cut efficiency.
+
+ beautyCombinedesd->Add(mccontaineresdD);
+ AliCFEffGrid* efficiencyBeautyCombinedesd = new AliCFEffGrid("efficiencyBeautyCombinedesd","",*beautyCombinedesd);
+ efficiencyBeautyCombinedesd->CalculateEfficiency(AliHFEcuts::kNcutStepsMCTrack + fStepData,AliHFEcuts::kNcutStepsMCTrack + fStepData-1);
+
source = 2; //conversion mb
if(fBeamType==0) mccontainermcD = GetSlicedContainer(containermb, fNbDimensions, dimensions, source, fChargeChoosen);
else mccontainermcD = GetSlicedContainer(containermb, fNbDimensions, dimensions, source, fChargeChoosen, icentr+1);
if(fIPEffCombinedSamples){
fEfficiencyCharmSigD[icentr] = (TH1D*)efficiencyCharmCombined->Project(ptpr); //signal enhenced + mb
fEfficiencyBeautySigD[icentr] = (TH1D*)efficiencyBeautyCombined->Project(ptpr); //signal enhenced + mb
+ fEfficiencyBeautySigesdD[icentr] = (TH1D*)efficiencyBeautyCombinedesd->Project(ptpr); //signal enhenced + mb
}
else{
fEfficiencyCharmSigD[icentr] = (TH1D*)efficiencyCharmSig->Project(ptpr); //signal enhenced only
fEfficiencyBeautySigD[icentr] = (TH1D*)efficiencyBeautySig->Project(ptpr); //signal enhenced only
+ fEfficiencyBeautySigesdD[icentr] = (TH1D*)efficiencyBeautySigesd->Project(ptpr); //signal enhenced only
}
fCharmEff[icentr] = (TH1D*)efficiencyCharm->Project(ptpr); //mb only
fBeautyEff[icentr] = (TH1D*)efficiencyBeauty->Project(ptpr); //mb only
}
+ if(fBeamType==0){
+ AliCFEffGrid *nonHFEEffGrid = (AliCFEffGrid*) GetEfficiency(GetContainer(kMCWeightedContainerNonHFEESD),1,0);
+ fNonHFEEffbgc = (TH1D *) nonHFEEffGrid->Project(0);
+
+ AliCFEffGrid *conversionEffGrid = (AliCFEffGrid*) GetEfficiency(GetContainer(kMCWeightedContainerConversionESD),1,0);
+ fConversionEffbgc = (TH1D *) conversionEffGrid->Project(0);
+ }
+
for(int icentr=0; icentr<loopcentr; icentr++) {
fipfit->SetParameters(0.5,0.319,0.0157,0.00664,6.77,2.08);
fipfit->SetLineColor(2);
if(fBeauty2ndMethod)fEfficiencyIPBeautyD[icentr] = fEfficiencyBeautySigD[0]->GetFunction("fipfit"); //why do we need this line?
else fEfficiencyIPBeautyD[icentr] = fEfficiencyBeautySigD[icentr]->GetFunction("fipfit");
+ fipfit->SetParameters(0.5,0.319,0.0157,0.00664,6.77,2.08);
+ fipfit->SetLineColor(6);
+ fEfficiencyBeautySigesdD[icentr]->Fit(fipfit,"R");
+ fEfficiencyBeautySigesdD[icentr]->GetYaxis()->SetTitle("Efficiency");
+ if(fBeauty2ndMethod)fEfficiencyIPBeautyesdD[icentr] = fEfficiencyBeautySigesdD[0]->GetFunction("fipfit"); //why do we need this line?
+ else fEfficiencyIPBeautyesdD[icentr] = fEfficiencyBeautySigesdD[icentr]->GetFunction("fipfit");
+
fipfit->SetParameters(0.5,0.319,0.0157,0.00664,6.77,2.08);
fipfit->SetLineColor(1);
fEfficiencyCharmSigD[icentr]->Fit(fipfit,"R");
fEfficiencyIPCharmD[icentr] = fEfficiencyCharmSigD[icentr]->GetFunction("fipfit");
if(fIPParameterizedEff){
- fipfit2->SetParameters(0.5,0.319,0.0157,0.00664,6.77,2.08);
- fipfit2->SetLineColor(3);
- fConversionEff[icentr]->Fit(fipfit2,"R");
+ fipfitnonhfe->SetParameters(0.5,0.319,0.0157,0.00664,6.77,2.08);
+ fipfitnonhfe->SetLineColor(3);
+ fConversionEff[icentr]->Fit(fipfitnonhfe,"R");
fConversionEff[icentr]->GetYaxis()->SetTitle("Efficiency");
- fEfficiencyIPConversionD[icentr] = fConversionEff[icentr]->GetFunction("fipfit2");
+ fEfficiencyIPConversionD[icentr] = fConversionEff[icentr]->GetFunction("fipfitnonhfe");
- fipfit2->SetParameters(0.5,0.319,0.0157,0.00664,6.77,2.08);
- fipfit2->SetLineColor(4);
- fNonHFEEff[icentr]->Fit(fipfit2,"R");
+ fipfitnonhfe->SetParameters(0.5,0.319,0.0157,0.00664,6.77,2.08);
+ fipfitnonhfe->SetLineColor(4);
+ fNonHFEEff[icentr]->Fit(fipfitnonhfe,"R");
fNonHFEEff[icentr]->GetYaxis()->SetTitle("Efficiency");
- fEfficiencyIPNonhfeD[icentr] = fNonHFEEff[icentr]->GetFunction("fipfit2");
+ fEfficiencyIPNonhfeD[icentr] = fNonHFEEff[icentr]->GetFunction("fipfitnonhfe");
}
}
fEfficiencyBeautySigD[0]->SetMarkerStyle(21);
fEfficiencyBeautySigD[0]->SetMarkerColor(2);
fEfficiencyBeautySigD[0]->SetLineColor(2);
+ fEfficiencyBeautySigesdD[0]->SetStats(0);
+ fEfficiencyBeautySigesdD[0]->SetMarkerStyle(21);
+ fEfficiencyBeautySigesdD[0]->SetMarkerColor(6);
+ fEfficiencyBeautySigesdD[0]->SetLineColor(6);
fCharmEff[0]->SetMarkerStyle(24);
fCharmEff[0]->SetMarkerColor(1);
fCharmEff[0]->SetLineColor(1);
fNonHFEEff[0]->SetLineColor(4);
fEfficiencyCharmSigD[0]->Draw();
+ fEfficiencyCharmSigD[0]->GetXaxis()->SetRangeUser(0.0,7.9);
+ fEfficiencyCharmSigD[0]->GetYaxis()->SetRangeUser(0.0,0.5);
+
fEfficiencyBeautySigD[0]->Draw("same");
+ fEfficiencyBeautySigesdD[0]->Draw("same");
//fCharmEff[0]->Draw("same");
//fBeautyEff[0]->Draw("same");
- fConversionEff[0]->Draw("same");
- fNonHFEEff[0]->Draw("same");
+
+ if(fBeamType==0){
+ fConversionEffbgc->SetMarkerStyle(25);
+ fConversionEffbgc->SetMarkerColor(3);
+ fConversionEffbgc->SetLineColor(3);
+ fNonHFEEffbgc->SetMarkerStyle(25);
+ fNonHFEEffbgc->SetMarkerColor(4);
+ fNonHFEEffbgc->SetLineColor(4);
+ fConversionEffbgc->Draw("same");
+ fNonHFEEffbgc->Draw("same");
+ }
+ else{
+ fConversionEff[0]->Draw("same");
+ fNonHFEEff[0]->Draw("same");
+ }
fEfficiencyIPBeautyD[0]->Draw("same");
+ fEfficiencyIPBeautyesdD[0]->Draw("same");
fEfficiencyIPCharmD[0]->Draw("same");
if(fIPParameterizedEff){
fEfficiencyIPConversionD[0]->Draw("same");
fEfficiencyIPNonhfeD[0]->Draw("same");
}
- TLegend *legipeff = new TLegend(0.55,0.2,0.85,0.39);
+ TLegend *legipeff = new TLegend(0.58,0.2,0.88,0.39);
legipeff->AddEntry(fEfficiencyBeautySigD[0],"IP Step Efficiency","");
legipeff->AddEntry(fEfficiencyBeautySigD[0],"beauty e","p");
+ legipeff->AddEntry(fEfficiencyBeautySigesdD[0],"beauty e(esd pt)","p");
legipeff->AddEntry(fEfficiencyCharmSigD[0],"charm e","p");
- legipeff->AddEntry(fConversionEff[0],"conversion e","p");
- legipeff->AddEntry(fNonHFEEff[0],"Dalitz e","p");
+ legipeff->AddEntry(fConversionEffbgc,"conversion e(esd pt)","p");
+ legipeff->AddEntry(fNonHFEEffbgc,"Dalitz e(esd pt)","p");
+ //legipeff->AddEntry(fConversionEff[0],"conversion e","p");
+ //legipeff->AddEntry(fNonHFEEff[0],"Dalitz e","p");
legipeff->Draw("same");
+ gPad->SetGrid();
+ //cefficiencyParamIP->SaveAs("efficiencyParamIP.eps");
}
//____________________________________________________________________________
-THnSparse* AliHFEspectrum::GetBeautyIPEff(){
+THnSparse* AliHFEspectrum::GetBeautyIPEff(Bool_t isMCpt){
//
// Return beauty electron IP cut efficiency
//
}
Double_t pt[1];
Double_t weight[nCentralitybinning];
+ Double_t weightErr[nCentralitybinning];
Double_t contents[2];
for(Int_t a=0;a<11;a++)
{
weight[a] = 1.0;
+ weightErr[a] = 1.0;
}
Int_t looppt=nBin[0];
Int_t loopcentr=1;
+ Int_t ibin[2];
if(fBeamType==1)
{
loopcentr=nBinPbPb[0];
for(int i=0; i<looppt; i++)
{
pt[0]=(kPtRange[i]+kPtRange[i+1])/2.;
- if(fEfficiencyIPBeautyD[icentr])
- weight[icentr]=fEfficiencyIPBeautyD[icentr]->Eval(pt[0]);
+ if(isMCpt){
+ if(fEfficiencyIPBeautyD[icentr]){
+ weight[icentr]=fEfficiencyIPBeautyD[icentr]->Eval(pt[0]);
+ weightErr[icentr] = 0;
+ }
+ else{
+ printf("Fit failed on beauty IP cut efficiency for centrality %d. Contents in histo used!\n",icentr);
+ weight[icentr] = fEfficiencyBeautySigD[icentr]->GetBinContent(i+1);
+ weightErr[icentr] = fEfficiencyBeautySigD[icentr]->GetBinError(i+1);
+ }
+ }
else{
- printf("Fit failed on beauty IP cut efficiency for centrality %d. Contents in histo used!\n",icentr);
- weight[icentr] = fEfficiencyBeautySigD[icentr]->GetBinContent(i+1);
+ if(fEfficiencyIPBeautyesdD[icentr]){
+ weight[icentr]=fEfficiencyIPBeautyesdD[icentr]->Eval(pt[0]);
+ weightErr[icentr] = 0;
+ }
+ else{
+ printf("Fit failed on beauty IP cut efficiency for centrality %d. Contents in histo used!\n",icentr);
+ weight[icentr] = fEfficiencyBeautySigesdD[icentr]->GetBinContent(i+1);
+ weightErr[icentr] = fEfficiencyBeautySigD[icentr]->GetBinError(i+1);
+ }
}
if(fBeamType==1){
contents[0]=icentr;
contents[1]=pt[0];
+ ibin[0]=icentr;
+ ibin[1]=i+1;
}
if(fBeamType==0){
contents[0]=pt[0];
+ ibin[0]=i+1;
}
ipcut->Fill(contents,weight[icentr]);
+ ipcut->SetBinError(ibin,weightErr[icentr]);
}
}
-
Int_t nDimSparse = ipcut->GetNdimensions();
Int_t* binsvar = new Int_t[nDimSparse]; // number of bins for each variable
Long_t nBins = 1; // used to calculate the total number of bins in the THnSparse
binfill[iVar] = 1 + bintmp % binsvar[iVar] ;
bintmp /= binsvar[iVar] ;
}
- ipcut->SetBinError(binfill,0.); // put 0 everywhere
+ //ipcut->SetBinError(binfill,0.); // put 0 everywhere
}
delete[] binsvar;
{
weight[a] = 1.0;
weightErr[a] = 1.0;
-
}
Int_t looppt=nBin[0];
nDim=2;
}
- // const Int_t nPtbinning1 = 35;//number of pt bins, according to new binning
const Int_t nPtbinning1 = 18;//number of pt bins, according to new binning
const Int_t nCentralitybinning=11;//number of centrality bins
- // Double_t kPtRange[nPtbinning1+1] = { 0., 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1., 1.1, 1.2, 1.3, 1.4, 1.5, 1.75, 2., 2.25, 2.5, 2.75, 3., 3.5, 4., 4.5, 5., 5.5, 6., 7., 8., 10., 12., 14., 16., 18., 20.};//pt bin limits
Double_t kPtRange[19] = {0, 0.1, 0.3, 0.5, 0.7, 0.9, 1.1, 1.3, 1.5, 2, 2.5, 3, 4, 5, 6, 8, 12, 16, 20};
Double_t kCentralityRange[nCentralitybinning+1] = {0.,1.,2., 3., 4., 5., 6., 7.,8.,9., 10., 11.};
Int_t nBin[1] = {nPtbinning1};
Int_t nBinPbPb[2] = {nCentralitybinning,nPtbinning1};
- //fBSpectrum2ndMethod->GetNbinsX()
AliCFDataGrid *rawBeautyContainer;
if(fBeamType==0) rawBeautyContainer = new AliCFDataGrid("rawBeautyContainer","rawBeautyContainer",nDim,nBin);
void SetEtaRange(Double_t etamin, Double_t etamax) { fEtaRange[0] = etamin; fEtaRange[1] = etamax; fEtaSelected = kTRUE; }
void SetUnSetCorrelatedErrors(Bool_t unsetcorrelatederrors) {fUnSetCorrelatedErrors = unsetcorrelatederrors;};
void SetSmoothing(Bool_t setSmoothing) {fSetSmoothing = setSmoothing;};
+ void SetTestOneBinCentrality(Double_t centralitymin, Double_t centralitymax) { fTestCentralityLow = centralitymin; fTestCentralityHigh = centralitymax;}
void SetNCentralityBinAtTheEnd(Int_t nCentralityBinAtTheEnd) {fNCentralityBinAtTheEnd = nCentralityBinAtTheEnd; };
void SetLowHighBoundaryCentralityBinAtTheEnd(Int_t low, Int_t high, Int_t i) { fLowBoundaryCentralityBinAtTheEnd[i] = low; fHighBoundaryCentralityBinAtTheEnd[i] = high;};
AliCFDataGrid* GetConversionBackground();
AliCFDataGrid* GetNonHFEBackground();
THnSparse* GetCharmWeights();
- THnSparse* GetBeautyIPEff();
+ THnSparse* GetBeautyIPEff(Bool_t isMCpt);
THnSparse* GetPIDxIPEff(Int_t source);
void CalculateNonHFEsyst(Int_t centrality = 0);
protected:
AliCFContainer *GetContainer(AliHFEspectrum::CFContainer_t contt);
- AliCFContainer *GetSlicedContainer(AliCFContainer *cont, Int_t ndim, Int_t *dimensions,Int_t source=-1,Chargetype_t charge=kAllCharge,Int_t centrality=-1);
- THnSparseF *GetSlicedCorrelation(THnSparseF *correlationmatrix,Int_t nDim, Int_t *dimensions) const;
+ AliCFContainer *GetSlicedContainer(AliCFContainer *cont, Int_t ndim, Int_t *dimensions,Int_t source=-1,Chargetype_t charge=kAllCharge,Int_t centralitylow=-1, Int_t centralityhigh=-1);
+ THnSparseF *GetSlicedCorrelation(THnSparseF *correlationmatrix,Int_t nDim, Int_t *dimensions,Int_t centralitylow=-1, Int_t centralityhigh=-1) const;
TObject* GetSpectrum(const AliCFContainer * const c, Int_t step);
TObject* GetEfficiency(const AliCFContainer * const c, Int_t step, Int_t step0);
TF1 *fEfficiencyesdTOFPIDD[kCentrality]; // TOF PID efficiency parameterized
TF1 *fEfficiencyIPCharmD[kCentrality]; // IP efficiency parameterized for charm
TF1 *fEfficiencyIPBeautyD[kCentrality]; // IP efficiency parameterized for beauty
+ TF1 *fEfficiencyIPBeautyesdD[kCentrality]; // IP efficiency parameterized for beauty for esd
TF1 *fEfficiencyIPConversionD[kCentrality]; // IP efficiency parameterized for conversion
TF1 *fEfficiencyIPNonhfeD[kCentrality]; // IP efficiency parameterized for nonhfe
Int_t fNCentralityBinAtTheEnd;// Number of centrality class at the end
Int_t fLowBoundaryCentralityBinAtTheEnd[20]; // Boundary of the bins
Int_t fHighBoundaryCentralityBinAtTheEnd[20]; // Boundary of the bins
+ Int_t fTestCentralityLow; // To test one bin in centrality only
+ Int_t fTestCentralityHigh; // To test one bin in centrality only
THnSparseF *fHadronEffbyIPcut;// container for hadron efficiency by IP cut
TH1D *fEfficiencyCharmSigD[kCentrality]; // charm IP cut eff from signal enhanced MC
TH1D *fEfficiencyBeautySigD[kCentrality]; // beauty IP cut eff from signal enhanced MC
+ TH1D *fEfficiencyBeautySigesdD[kCentrality]; // beauty IP cut eff from signal enhanced MC for esd
TH1D *fConversionEff[kCentrality]; // conversion IP cut eff
TH1D *fNonHFEEff[kCentrality]; // nonhfe IP cut eff
TH1D *fCharmEff[kCentrality]; // charm IP cut eff
// Authors
// All authors of the HFE group
//
+#include <TArrayD.h>
#include <TMath.h>
#include <TParticle.h>
#include "TH1.h"
return bins;
}
+//__________________________________________
+void AliHFEtools::FillLinearBinning(TArrayD &bins, Int_t nBins, Double_t ymin, Double_t ymax){
+ //
+ // Helper function for linearly binned array
+ //
+ Double_t stepsize = (ymax - ymin) / static_cast<Double_t>(nBins);
+ bins[0] = ymin;
+ for(Int_t ibin = 1; ibin <= nBins; ibin++)
+ bins[ibin] = bins[ibin-1] + stepsize;
+}
+
//__________________________________________
Double_t *AliHFEtools::MakeLogarithmicBinning(Int_t nBins, Double_t ymin, Double_t ymax){
//
return bins;
}
+//__________________________________________
+void AliHFEtools::FillLogarithmicBinning(TArrayD &bins, Int_t nBins, Double_t ymin, Double_t ymax){
+ //
+ // Helper function for logartimically binned array
+ //
+ bins[0] = ymin;
+ Double_t factor = TMath::Power(ymax/ymin, 1./nBins);
+ for(Int_t ibin = 1; ibin <= nBins; ibin++)
+ bins[ibin] = factor * bins[ibin-1];
+}
+
//_________________________________________
Bool_t AliHFEtools::BinLogAxis(TObject *o, Int_t dim){
AliError(Form("Log binning not possible for this axis [min = %f]\n", from));
}
Double_t to = axis->GetXmax();
- Double_t *newBins = new Double_t[bins+1];
+ TArrayD newBins(bins+1);
newBins[0] = from;
Double_t factor = TMath::Power(to/from, 1./bins);
for(Int_t i=1; i<=bins; ++i){
newBins[i] = factor * newBins[i-1];
}
- axis->Set(bins, newBins);
- delete[] newBins;
+ axis->Set(bins, newBins.GetArray());
return kTRUE;
}
#include <TObject.h>
+class TArrayD;
class TParticle;
class AliAODMCParticle;
class AliPIDResponse;
static Double_t *MakeLinearBinning(Int_t nBins, Double_t ymin, Double_t ymax);
static Double_t *MakeLogarithmicBinning(Int_t nBins, Double_t ymin, Double_t ymax);
+ static void FillLinearBinning(TArrayD &bins, Int_t nBins, Double_t ymin, Double_t ymax);
+ static void FillLogarithmicBinning(TArrayD &bins, Int_t nBins, Double_t ymin, Double_t ymax);
Bool_t BinLogAxis(TObject *o, Int_t dim);
static Float_t GetRapidity(const TParticle *part);
static Float_t GetRapidity(const AliAODMCParticle *part); // return rapidity
-AliAnalysisTask *AddTaskHFEtpctofv2(Int_t tpcCls=110, Double_t tpcClsr=50, Int_t tpcClspid=60, Double_t tpcsharedfraction=95, Int_t itsCls=4, Double_t chi2peritscl=36, Int_t pixellayer=2, Double_t dcaxy=100,Double_t dcaz=200, Double_t tofsig=30., Double_t tpcdedx0=0.0, Double_t tpcdedx1=0.0, Double_t tpcdedx2=0.0, Double_t tpcdedx3=0.0, Double_t tpcdedx4=0.0, Int_t vzero=3, Int_t debuglevel=3){
+AliAnalysisTask *AddTaskHFEtpctofv2(Int_t tpcCls=110, Double_t tpcClsr=50, Int_t tpcClspid=60, Double_t tpcsharedfraction=10, Int_t itsCls=4, Double_t chi2peritscl=36, Int_t pixellayer=2, Double_t dcaxy=100,Double_t dcaz=200, Double_t tofsig=30., Double_t tpcdedx0=0.0, Double_t tpcdedx1=0.0, Double_t tpcdedx2=0.0, Double_t tpcdedx3=0.0, Double_t tpcdedx4=0.0, Int_t vzero=3, Int_t debuglevel=4){
// libraries in case
gSystem->Load("libANALYSIS.so");
printf("dcaxy %f\n",dcaxy*0.01);
printf("dcaz %f\n",dcaz*0.01);
printf("TOF sigma %f\n",tofsig*0.1);
- printf("TPC min sigma cut 0: %f\n",tpcdedx0*0.01);
- printf("TPC min sigma cut 1: %f\n",tpcdedx1*0.01);
- printf("TPC min sigma cut 2: %f\n",tpcdedx2*0.01);
- printf("TPC min sigma cut 3: %f\n",tpcdedx3*0.01);
- printf("TPC min sigma cut 4: %f\n",tpcdedx4*0.01);
+ printf("TPC min sigma cut 0: %f\n",tpcdedx0*0.001);
+ printf("TPC min sigma cut 1: %f\n",tpcdedx1*0.001);
+ printf("TPC min sigma cut 2: %f\n",tpcdedx2*0.001);
+ printf("TPC min sigma cut 3: %f\n",tpcdedx3*0.001);
+ printf("TPC min sigma cut 4: %f\n",tpcdedx4*0.001);
printf("VZERO event plane %d\n",vzero);
printf("Debug level %d\n",debuglevel);
hfecuts->SetTOFPIDStep(kTRUE);
+ // Cut HFE background
+ AliESDtrackCuts *hfeBackgroundCuts = new AliESDtrackCuts();
+ hfeBackgroundCuts->SetName("backgroundcuts");
+ hfeBackgroundCuts->SetAcceptKinkDaughters(kFALSE);
+ hfeBackgroundCuts->SetRequireTPCRefit(kTRUE);
+ hfeBackgroundCuts->SetEtaRange(-0.9,0.9);
+ hfeBackgroundCuts->SetRequireSigmaToVertex(kTRUE);
+ hfeBackgroundCuts->SetMaxChi2PerClusterTPC(4.0);
+ hfeBackgroundCuts->SetMinNClustersTPC(50);
+ hfeBackgroundCuts->SetPtRange(0.3,1e10);
+
// Name
TString appendix(TString::Format("TPC%dTPCr%dTPCpid%dTPCShared%dITScl%dChi2perITS%dPixelLayer%dDCAr%dz%dTOFsig%dTPCmindedx0%dTPCmindedx1%dTPCmindedx2%dTPCmindedx3%dTPCmindedx4%dVZERO%dDebugLevel%d",tpcCls,(Int_t)tpcClsr,tpcClspid,(Int_t) tpcsharedfraction,itsCls,(Int_t) chi2peritscl,(Int_t) pixellayer,(Int_t)dcaxy,(Int_t)dcaz,(Int_t)tofsig,(Int_t)tpcdedx0,(Int_t)tpcdedx1,(Int_t)tpcdedx2,(Int_t)tpcdedx3,(Int_t)tpcdedx4,vzero,debuglevel));
printf("appendix %s\n", appendix.Data());
task->SetDebugLevel(1);
task->GetPIDQAManager()->SetHighResolutionHistos();
task->SetHFECuts(hfecuts);
+ task->SetHFEBackgroundCuts(hfeBackgroundCuts);
if(useMC) {
task->SetMCPID(kTRUE);
- //task->SetUseMCReactionPlane(kTRUE);
+ task->SetUseMCReactionPlane(kTRUE);
task->SetAfterBurnerOn(kTRUE);
task->SetV1V2V3V4V5(0.0,0.2,0.0,0.0,0.0);
}
Double_t params_centr_10_20[1];
Double_t params_centr_20_30[1];
Double_t params_centr_per[1];
- params_centr_0_5[0]=tpcdedx0*0.01; // cut tuned for 0-10%
- params_centr_5_10[0]=tpcdedx1*0.01; // cut tuned for 0-10%
- params_centr_10_20[0]=tpcdedx2*0.01;
- params_centr_20_30[0]=tpcdedx3*0.01;
- params_centr_per[0]=tpcdedx4*0.01;
+ params_centr_0_5[0]=tpcdedx0*0.001; // cut tuned for 0-10%
+ params_centr_5_10[0]=tpcdedx1*0.001; // cut tuned for 0-10%
+ params_centr_10_20[0]=tpcdedx2*0.001;
+ params_centr_20_30[0]=tpcdedx3*0.001;
+ params_centr_per[0]=tpcdedx4*0.001;
char *cutmodel;
cutmodel="pol0";
}
pid->ConfigureTOF(tofsig*0.1);
+
+ // Define PID background
+ AliHFEpid *pidbackground = task->GetPIDBackground();
+ if(useMC) pidbackground->SetHasMCData(kTRUE);
+ pidbackground->AddDetector("TOF", 0);
+ pidbackground->AddDetector("TPC", 1);
+ pidbackground->ConfigureTOF(3.0);
+ pidbackground->ConfigureTPCasymmetric(0.1,20.,-2.0,5.0);
+ task->SetMaxopeningtheta(9999.0);
+ task->SetMaxopeningphi(9999.0);
+ task->SetMaxopening3D(0.5);
printf("*************************************\n");
printf("Configuring standard Task:\n");
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff