From 587d006ad0d7a8f2f8082e8c82076baac57ad43d Mon Sep 17 00:00:00 2001 From: iseliouj Date: Fri, 2 Mar 2012 09:19:20 +0000 Subject: [PATCH] From Francesco: Added AliFlowVZEROResults, PID-VZERO class restructuring --- PWG/CMakelibPWGflowTasks.pkg | 1 + PWG/FLOW/Tasks/AliAnalysisTaskVnV0.cxx | 659 ++++++++++++++++--------- PWG/FLOW/Tasks/AliAnalysisTaskVnV0.h | 22 +- PWG/FLOW/Tasks/AliFlowBayesianPID.cxx | 5 +- PWG/FLOW/Tasks/AliFlowBayesianPID.h | 5 + PWG/PWGflowTasksLinkDef.h | 1 + PWGCF/FLOW/macros/AddTaskVZERO.C | 17 +- PWGCF/FLOW/macros/extractFlowVZERO.C | 302 ++--------- 8 files changed, 507 insertions(+), 505 deletions(-) diff --git a/PWG/CMakelibPWGflowTasks.pkg b/PWG/CMakelibPWGflowTasks.pkg index 176f139aa7f..135c0230603 100644 --- a/PWG/CMakelibPWGflowTasks.pkg +++ b/PWG/CMakelibPWGflowTasks.pkg @@ -52,6 +52,7 @@ set ( SRCS FLOW/Tasks/AliAnalysisTaskPhiFlow.cxx FLOW/Tasks/AliAnalysisTaskFilterFE.cxx FLOW/Tasks/AliAnalysisTaskVnV0.cxx + FLOW/Tasks/AliFlowVZEROResults.cxx ) string ( REPLACE ".cxx" ".h" HDRS "${SRCS}" ) diff --git a/PWG/FLOW/Tasks/AliAnalysisTaskVnV0.cxx b/PWG/FLOW/Tasks/AliAnalysisTaskVnV0.cxx index 10c018f8f2e..885fe1c3581 100644 --- a/PWG/FLOW/Tasks/AliAnalysisTaskVnV0.cxx +++ b/PWG/FLOW/Tasks/AliAnalysisTaskVnV0.cxx @@ -15,6 +15,10 @@ #include "AliOADBContainer.h" #include "TH2F.h" #include "TF1.h" +#include "AliGenHijingEventHeader.h" +#include "AliMCEvent.h" +#include "AliAODMCHeader.h" +#include "AliAODMCParticle.h" // STL includes //#include @@ -22,48 +26,39 @@ ClassImp(AliAnalysisTaskVnV0) -/* -fMeanQ[nCentrBin][2][2]; // and recentering - Float_t fWidthQ[nCentrBin][2][2]; // ... - Float_t fMeanQv3[nCentrBin][2][2]; // also for v3 - Float_t fWidthQv3[nCentrBin][2][2]; // ... -*/ - //_____________________________________________________________________________ AliAnalysisTaskVnV0::AliAnalysisTaskVnV0(): AliAnalysisTaskSE(), fAOD(0), - fVtxCut(10), // cut on |vertex| < fVtxCut + fVtxCut(10.0), // cut on |vertex| < fVtxCut fEtaCut(0.8), // cut on |eta| < fEtaCut fMinPt(0.15), // cut on pt > fMinPt fRun(-1), fList(new TList()), fList2(new TList()), - fMultV0(0), + fList3(new TList()), + fList4(new TList()), + fMultV0(NULL), fV0Cpol(100), fV0Apol(100), - fContAllChargesV0A(0), - fContAllChargesV0C(0), - fContAllChargesV0Av3(0), - fContAllChargesV0Cv3(0), - fHResTPCv0A2(0), - fHResTPCv0C2(0), - fHResv0Cv0A2(0), - fHResTPCv0A3(0), - fHResTPCv0C3(0), - fHResv0Cv0A3(0), - fPhiRPv0A(0), - fPhiRPv0C(0), - fPhiRPv0Av3(0), - fPhiRPv0Cv3(0), - fPhiTracks(0), - fQA(0), - fQA2(0), - fQAv3(0), - fQA2v3(0), + fHResTPCv0A2(NULL), + fHResTPCv0C2(NULL), + fHResv0Cv0A2(NULL), + fHResTPCv0A3(NULL), + fHResTPCv0C3(NULL), + fHResv0Cv0A3(NULL), + fPhiRPv0A(NULL), + fPhiRPv0C(NULL), + fPhiRPv0Av3(NULL), + fPhiRPv0Cv3(NULL), + fPhiTracks(NULL), + fQA(NULL), + fQA2(NULL), + fQAv3(NULL), + fQA2v3(NULL), fPID(new AliFlowBayesianPID()), - fTree(0), - fCentrality(0), + fTree(NULL), + fCentrality(-1), evPlAngV0ACor2(0), evPlAngV0CCor2(0), evPlAng2(0), @@ -71,54 +66,63 @@ AliAnalysisTaskVnV0::AliAnalysisTaskVnV0(): evPlAngV0CCor3(0), evPlAng3(0), fV2(kTRUE), - fV3(kTRUE) + fV3(kTRUE), + fContAllChargesV0A(NULL), + fContAllChargesV0C(NULL), + fContAllChargesV0Av3(NULL), + fContAllChargesV0Cv3(NULL), + fContAllChargesMC(NULL), + fContAllChargesMCv3(NULL), + fIsMC(kFALSE), + fQAsw(kFALSE) { DefineOutput(1, TList::Class()); DefineOutput(2, TList::Class()); + DefineOutput(3, TList::Class()); + DefineOutput(4, TList::Class()); // Default constructor (should not be used) fList->SetName("resultsV2"); fList2->SetName("resultsV3"); + fList3->SetName("resultsMC"); + fList4->SetName("QA"); fPID->SetNewTrackParam(); // Better tuning for TOF PID tracking effect in LHC10h - } //______________________________________________________________________________ AliAnalysisTaskVnV0::AliAnalysisTaskVnV0(const char *name): AliAnalysisTaskSE(name), fAOD(0), - fVtxCut(10), // cut on |vertex| < fVtxCut + fVtxCut(10.0), // cut on |vertex| < fVtxCut fEtaCut(0.8), // cut on |eta| < fEtaCut fMinPt(0.15), // cut on pt > fMinPt fRun(-1), fList(new TList()), fList2(new TList()), - fMultV0(0), + fList3(new TList()), + fList4(new TList()), + fMultV0(NULL), fV0Cpol(100), fV0Apol(100), - fContAllChargesV0A(0), - fContAllChargesV0C(0), - fContAllChargesV0Av3(0), - fContAllChargesV0Cv3(0), - fHResTPCv0A2(0), - fHResTPCv0C2(0), - fHResv0Cv0A2(0), - fHResTPCv0A3(0), - fHResTPCv0C3(0), - fHResv0Cv0A3(0), - fPhiRPv0A(0), - fPhiRPv0C(0), - fPhiRPv0Av3(0), - fPhiRPv0Cv3(0), - fPhiTracks(0), - fQA(0), - fQA2(0), - fQAv3(0), - fQA2v3(0), + fHResTPCv0A2(NULL), + fHResTPCv0C2(NULL), + fHResv0Cv0A2(NULL), + fHResTPCv0A3(NULL), + fHResTPCv0C3(NULL), + fHResv0Cv0A3(NULL), + fPhiRPv0A(NULL), + fPhiRPv0C(NULL), + fPhiRPv0Av3(NULL), + fPhiRPv0Cv3(NULL), + fPhiTracks(NULL), + fQA(NULL), + fQA2(NULL), + fQAv3(NULL), + fQA2v3(NULL), fPID(new AliFlowBayesianPID()), - fTree(0), - fCentrality(0), + fTree(NULL), + fCentrality(-1), evPlAngV0ACor2(0), evPlAngV0CCor2(0), evPlAng2(0), @@ -126,15 +130,27 @@ AliAnalysisTaskVnV0::AliAnalysisTaskVnV0(const char *name): evPlAngV0CCor3(0), evPlAng3(0), fV2(kTRUE), - fV3(kTRUE) + fV3(kTRUE), + fContAllChargesV0A(NULL), + fContAllChargesV0C(NULL), + fContAllChargesV0Av3(NULL), + fContAllChargesV0Cv3(NULL), + fContAllChargesMC(NULL), + fContAllChargesMCv3(NULL), + fIsMC(kFALSE), + fQAsw(kFALSE) { DefineOutput(1, TList::Class()); DefineOutput(2, TList::Class()); + DefineOutput(3, TList::Class()); + DefineOutput(4, TList::Class()); // Output slot #1 writes into a TTree fList->SetName("resultsV2"); fList2->SetName("resultsV3"); + fList3->SetName("resultsMC"); + fList4->SetName("QA"); fPID->SetNewTrackParam(); // Better tuning for TOF PID tracking effect in LHC10h } @@ -149,6 +165,9 @@ AliAnalysisTaskVnV0::~AliAnalysisTaskVnV0() void AliAnalysisTaskVnV0::UserCreateOutputObjects() { + if(fIsMC) fPID->SetMC(kTRUE); + + // Tree for EP debug (comment the adding to v2 list id not needed) fTree = new TTree("tree","tree"); fTree->Branch("evPlAngV0ACor2",&evPlAngV0ACor2,"evPlAngV0ACor2/F"); @@ -165,99 +184,147 @@ void AliAnalysisTaskVnV0::UserCreateOutputObjects() Double_t binsPtTOF[nPtBinsTOF+1] = {0., 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0, 2.25, 2.5, 2.75,3.0,3.25,3.5,3.75,4.0,4.5,5,5.5,6,6.5,7,8,9,10,12,15,20}; const Int_t nChargeBinsTOF = 3; Double_t binChargeTOF[nChargeBinsTOF+1] = {-1.5,-0.5,0.5,1.5}; - const Int_t nv_2BinsTOF = 50; - Double_t binV_2TOF[nv_2BinsTOF+1]; - for(Int_t i=0;iSetBinLimits(0,binCentrTOF); - fContAllChargesV0A->SetBinLimits(1,binV_2TOF); - fContAllChargesV0A->SetBinLimits(2,binChargeTOF); - fContAllChargesV0A->SetBinLimits(3,binsPtTOF); - fContAllChargesV0A->SetBinLimits(4,binProbTOF); - fContAllChargesV0A->SetBinLimits(5,binPsiTOF); - fContAllChargesV0A->SetBinLimits(6,binMaskPID); - fContAllChargesV0C->SetBinLimits(0,binCentrTOF); - fContAllChargesV0C->SetBinLimits(1,binV_2TOF); - fContAllChargesV0C->SetBinLimits(2,binChargeTOF); - fContAllChargesV0C->SetBinLimits(3,binsPtTOF); - fContAllChargesV0C->SetBinLimits(4,binProbTOF); - fContAllChargesV0C->SetBinLimits(5,binPsiTOF); - fContAllChargesV0C->SetBinLimits(6,binMaskPID); - - fContAllChargesV0A->SetVarTitle(0,"centrality"); - fContAllChargesV0A->SetVarTitle(1,"v_{2}"); - fContAllChargesV0A->SetVarTitle(2,"charge"); - fContAllChargesV0A->SetVarTitle(3,"p_{t} (GeV/c)"); - fContAllChargesV0A->SetVarTitle(4,"Bayesian Probability"); - fContAllChargesV0A->SetVarTitle(5,"#Psi V0A"); - fContAllChargesV0A->SetVarTitle(6,"TOF PID"); - fContAllChargesV0C->SetVarTitle(0,"centrality"); - fContAllChargesV0C->SetVarTitle(1,"v_{2}"); - fContAllChargesV0C->SetVarTitle(2,"charge"); - fContAllChargesV0C->SetVarTitle(3,"p_{t} (GeV/c)"); - fContAllChargesV0C->SetVarTitle(4,"Bayesian Probability"); - fContAllChargesV0C->SetVarTitle(5,"#Psi V0C"); - fContAllChargesV0C->SetVarTitle(6,"TOF PID"); + fContAllChargesV0A = new AliFlowVZEROResults("v2A",5,binsTOF); + fContAllChargesV0A->SetVarRange(0,-0.5,8.5); // centrality + fContAllChargesV0A->SetVarRange(1,-1.5,1.5); // charge + fContAllChargesV0A->SetVarRange(2,0.6,1.0001);// prob + fContAllChargesV0A->SetVarRange(3,-TMath::Pi()/2,TMath::Pi()/2); // Psi + fContAllChargesV0A->SetVarRange(4,-0.5,2.5); // pid mask + fContAllChargesV0A->SetVarName(0,"centrality"); + fContAllChargesV0A->SetVarName(1,"charge"); + fContAllChargesV0A->SetVarName(2,"prob"); + fContAllChargesV0A->SetVarName(3,"#Psi"); + fContAllChargesV0A->SetVarName(4,"PIDmask"); + fContAllChargesV0A->AddSpecies("all",nPtBinsTOF,binsPtTOF); + fContAllChargesV0A->AddSpecies("pi",nPtBinsTOF,binsPtTOF); + fContAllChargesV0A->AddSpecies("k",nPtBinsTOF,binsPtTOF); + fContAllChargesV0A->AddSpecies("pr",nPtBinsTOF,binsPtTOF); + fContAllChargesV0A->AddSpecies("e",nPtBinsTOF,binsPtTOF); + fContAllChargesV0A->AddSpecies("d",nPtBinsTOF,binsPtTOF); + fContAllChargesV0A->AddSpecies("t",nPtBinsTOF,binsPtTOF); + fContAllChargesV0A->AddSpecies("he3",nPtBinsTOF,binsPtTOF); + fContAllChargesV0A->AddSpecies("mu",nPtBinsTOF,binsPtTOF); + + fContAllChargesV0C = new AliFlowVZEROResults("v2C",5,binsTOF); + fContAllChargesV0C->SetVarRange(0,-0.5,8.5); // centrality + fContAllChargesV0C->SetVarRange(1,-1.5,1.5); // charge + fContAllChargesV0C->SetVarRange(2,0.6,1.0001);// prob + fContAllChargesV0C->SetVarRange(3,-TMath::Pi()/2,TMath::Pi()/2); // Psi + fContAllChargesV0C->SetVarRange(4,-0.5,2.5); // pid mask + fContAllChargesV0C->SetVarName(0,"centrality"); + fContAllChargesV0C->SetVarName(1,"charge"); + fContAllChargesV0C->SetVarName(2,"prob"); + fContAllChargesV0C->SetVarName(3,"#Psi"); + fContAllChargesV0C->SetVarName(4,"PIDmask"); + fContAllChargesV0C->AddSpecies("all",nPtBinsTOF,binsPtTOF); + fContAllChargesV0C->AddSpecies("pi",nPtBinsTOF,binsPtTOF); + fContAllChargesV0C->AddSpecies("k",nPtBinsTOF,binsPtTOF); + fContAllChargesV0C->AddSpecies("pr",nPtBinsTOF,binsPtTOF); + fContAllChargesV0C->AddSpecies("e",nPtBinsTOF,binsPtTOF); + fContAllChargesV0C->AddSpecies("d",nPtBinsTOF,binsPtTOF); + fContAllChargesV0C->AddSpecies("t",nPtBinsTOF,binsPtTOF); + fContAllChargesV0C->AddSpecies("he3",nPtBinsTOF,binsPtTOF); + fContAllChargesV0C->AddSpecies("mu",nPtBinsTOF,binsPtTOF); fList->Add(fContAllChargesV0A); fList->Add(fContAllChargesV0C); + if(fIsMC){ + fContAllChargesMC = new AliFlowVZEROResults("v2mc",5,binsTOFmc); + fContAllChargesMC->SetVarRange(0,-0.5,8.5); // centrality + fContAllChargesMC->SetVarRange(1,-1.5,1.5); // charge + fContAllChargesMC->SetVarRange(2,0.6,1.0001);// prob + fContAllChargesMC->SetVarRange(3,-TMath::Pi()/2,TMath::Pi()/2); // Psi + fContAllChargesMC->SetVarRange(4,-0.5,1.5); // pid mask + fContAllChargesMC->SetVarName(0,"centrality"); + fContAllChargesMC->SetVarName(1,"charge"); + fContAllChargesMC->SetVarName(2,"prob"); + fContAllChargesMC->SetVarName(3,"#Psi"); + fContAllChargesMC->SetVarName(4,"PIDmask"); + fContAllChargesMC->AddSpecies("all",nPtBinsTOF,binsPtTOF); + fContAllChargesMC->AddSpecies("pi",nPtBinsTOF,binsPtTOF); + fContAllChargesMC->AddSpecies("k",nPtBinsTOF,binsPtTOF); + fContAllChargesMC->AddSpecies("pr",nPtBinsTOF,binsPtTOF); + fContAllChargesMC->AddSpecies("e",nPtBinsTOF,binsPtTOF); + fContAllChargesMC->AddSpecies("mu",nPtBinsTOF,binsPtTOF); + fList3->Add(fContAllChargesMC); + } + // v3 container - fContAllChargesV0Av3 = new AliCFContainer("fContAllChargesV0Av3","centr:v3:charge:pt:prob:Psi:PIDmask",1+7,7,binsTOF); - fContAllChargesV0Cv3 = new AliCFContainer("fContAllChargesV0Cv3","centr:v3:charge:pt:prob:Psi:PIDmask",1+7,7,binsTOF); - fContAllChargesV0Av3->SetBinLimits(0,binCentrTOF); - fContAllChargesV0Av3->SetBinLimits(1,binV_2TOF); - fContAllChargesV0Av3->SetBinLimits(2,binChargeTOF); - fContAllChargesV0Av3->SetBinLimits(3,binsPtTOF); - fContAllChargesV0Av3->SetBinLimits(4,binProbTOF); - fContAllChargesV0Av3->SetBinLimits(5,binPsiTOFv3); - fContAllChargesV0Av3->SetBinLimits(6,binMaskPID); - fContAllChargesV0Cv3->SetBinLimits(0,binCentrTOF); - fContAllChargesV0Cv3->SetBinLimits(1,binV_2TOF); - fContAllChargesV0Cv3->SetBinLimits(2,binChargeTOF); - fContAllChargesV0Cv3->SetBinLimits(3,binsPtTOF); - fContAllChargesV0Cv3->SetBinLimits(4,binProbTOF); - fContAllChargesV0Cv3->SetBinLimits(5,binPsiTOFv3); - fContAllChargesV0Cv3->SetBinLimits(6,binMaskPID); - - fContAllChargesV0Av3->SetVarTitle(0,"centrality"); - fContAllChargesV0Av3->SetVarTitle(1,"v_{2}"); - fContAllChargesV0Av3->SetVarTitle(2,"charge"); - fContAllChargesV0Av3->SetVarTitle(3,"p_{t} (GeV/c)"); - fContAllChargesV0Av3->SetVarTitle(4,"Bayesian Probability"); - fContAllChargesV0Av3->SetVarTitle(5,"#Psi V0A"); - fContAllChargesV0Av3->SetVarTitle(6,"TOF PID"); - fContAllChargesV0Cv3->SetVarTitle(0,"centrality"); - fContAllChargesV0Cv3->SetVarTitle(1,"v_{3}"); - fContAllChargesV0Cv3->SetVarTitle(2,"charge"); - fContAllChargesV0Cv3->SetVarTitle(3,"p_{t} (GeV/c)"); - fContAllChargesV0Cv3->SetVarTitle(4,"Bayesian Probability"); - fContAllChargesV0Cv3->SetVarTitle(5,"#Psi V0C"); - fContAllChargesV0Cv3->SetVarTitle(6,"TOF PID"); + fContAllChargesV0Av3 = new AliFlowVZEROResults("v3A",5,binsTOF); + fContAllChargesV0Av3->SetVarRange(0,-0.5,8.5); // centrality + fContAllChargesV0Av3->SetVarRange(1,-1.5,1.5); // charge + fContAllChargesV0Av3->SetVarRange(2,0.6,1.0001);// prob + fContAllChargesV0Av3->SetVarRange(3,-TMath::Pi()/3,TMath::Pi()/3); // Psi + fContAllChargesV0Av3->SetVarRange(4,-0.5,2.5); // pid mask + fContAllChargesV0Av3->SetVarName(0,"centrality"); + fContAllChargesV0Av3->SetVarName(1,"charge"); + fContAllChargesV0Av3->SetVarName(2,"prob"); + fContAllChargesV0Av3->SetVarName(3,"#Psi"); + fContAllChargesV0Av3->SetVarName(4,"PIDmask"); + fContAllChargesV0Av3->AddSpecies("all",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Av3->AddSpecies("pi",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Av3->AddSpecies("k",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Av3->AddSpecies("pr",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Av3->AddSpecies("e",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Av3->AddSpecies("d",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Av3->AddSpecies("t",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Av3->AddSpecies("he3",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Av3->AddSpecies("mu",nPtBinsTOF,binsPtTOF); + + fContAllChargesV0Cv3 = new AliFlowVZEROResults("v3C",5,binsTOF); + fContAllChargesV0Cv3->SetVarRange(0,-0.5,8.5); // centrality + fContAllChargesV0Cv3->SetVarRange(1,-1.5,1.5); // charge + fContAllChargesV0Cv3->SetVarRange(2,0.6,1.0001);// prob + fContAllChargesV0Cv3->SetVarRange(3,-TMath::Pi()/3,TMath::Pi()/3); // Psi + fContAllChargesV0Cv3->SetVarRange(4,-0.5,2.5); // pid mask + fContAllChargesV0Cv3->SetVarName(0,"centrality"); + fContAllChargesV0Cv3->SetVarName(1,"charge"); + fContAllChargesV0Cv3->SetVarName(2,"prob"); + fContAllChargesV0Cv3->SetVarName(3,"#Psi"); + fContAllChargesV0Cv3->SetVarName(4,"PIDmask"); + fContAllChargesV0Cv3->AddSpecies("all",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Cv3->AddSpecies("pi",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Cv3->AddSpecies("k",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Cv3->AddSpecies("pr",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Cv3->AddSpecies("e",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Cv3->AddSpecies("d",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Cv3->AddSpecies("t",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Cv3->AddSpecies("he3",nPtBinsTOF,binsPtTOF); + fContAllChargesV0Cv3->AddSpecies("mu",nPtBinsTOF,binsPtTOF); fList2->Add(fContAllChargesV0Av3); fList2->Add(fContAllChargesV0Cv3); @@ -316,6 +383,11 @@ void AliAnalysisTaskVnV0::UserCreateOutputObjects() for(Int_t i=0;iAdd(fPhiRPv0A); fList->Add(fPhiRPv0C); - // fList->Add(fPhiTracks); // comment if not needed - fList->Add(fQA); - fList->Add(fQA2); + if(fQAsw) + fList4->Add(fPhiTracks); // comment if not needed + + if(fQAsw && fV2){ + fList4->Add(fQA); + fList4->Add(fQA2); + } fList2->Add(fPhiRPv0Av3); fList2->Add(fPhiRPv0Cv3); - fList2->Add(fQAv3); - fList2->Add(fQA2v3); - // fList->Add(fTree); // comment if not needed + if(fQAsw && fV3){ + fList4->Add(fQAv3); + fList4->Add(fQA2v3); + } + + fList->Add(fTree); // comment if not needed printf("Output creation ok!!\n\n\n\n"); // Post output data. if(fV2) PostData(1, fList); if(fV3) PostData(2, fList2); + if(fIsMC) PostData(3, fList3); + if(fQAsw) PostData(4, fList4); } //______________________________________________________________________________ @@ -442,6 +523,36 @@ void AliAnalysisTaskVnV0::UserExec(Option_t *) Float_t zvtx = GetVertex(fAOD); + + + //Get the MC object + + AliAODMCHeader *mcHeader = dynamic_cast(fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName())); + if (!mcHeader) { + AliError("Could not find MC Header in AOD"); + return; + } + + /* + AliMCEvent* mcEvent = MCEvent(); + if (!mcEvent) { + Printf("ERROR: Could not retrieve MC event"); + return; + } + + Double_t gReactionPlane = -999., gImpactParameter = -999.; + //Get the MC header + AliGenHijingEventHeader* headerH = dynamic_cast(mcEvent->GenEventHeader()); + if (headerH) { + //Printf("====================================================="); + //Printf("Reaction plane angle: %lf",headerH->ReactionPlaneAngle()); + //Printf("====================================================="); + gReactionPlane = headerH->ReactionPlaneAngle(); + gImpactParameter = headerH->ImpactParameter(); + } + +*/ + if (TMath::Abs(zvtx) < fVtxCut) { //Centrality Float_t v0Centr = -10.; @@ -490,7 +601,23 @@ void AliAnalysisTaskVnV0::Analyze(AliAODEvent* aodEvent, Float_t v0Centr) Double_t Qxc2 = 0, Qyc2 = 0; Double_t Qxa3 = 0, Qya3 = 0; Double_t Qxc3 = 0, Qyc3 = 0; - + + Int_t nAODTracks = aodEvent->GetNumberOfTracks(); + + AliAODMCHeader *mcHeader = NULL; + TClonesArray *mcArray = NULL; + Float_t evplaneMC = 0; + if(fIsMC){ + mcHeader = dynamic_cast(fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName())); + evplaneMC = mcHeader->GetReactionPlaneAngle(); + if(evplaneMC > TMath::Pi()/2 && evplaneMC <= TMath::Pi()*3/2) evplaneMC-=TMath::Pi(); + else if(evplaneMC > TMath::Pi()*3/2) evplaneMC-=2*TMath::Pi(); + + if (mcHeader) { + mcArray = (TClonesArray*)fAOD->GetList()->FindObject(AliAODMCParticle::StdBranchName()); + } + } + //V0 info AliAODVZERO* aodV0 = aodEvent->GetVZEROData(); @@ -543,7 +670,6 @@ void AliAnalysisTaskVnV0::Analyze(AliAODEvent* aodEvent, Float_t v0Centr) evPlAngV0ACor3 = TMath::ATan2(QyaCor3, QxaCor3)/3.; evPlAngV0CCor3 = TMath::ATan2(QycCor3, QxcCor3)/3.; - Int_t nAODTracks = aodEvent->GetNumberOfTracks(); //loop track and get pid for(Int_t iT = 0; iT < nAODTracks; iT++) { // loop on the tracks AliAODTrack* aodTrack = aodEvent->GetTrack(iT); @@ -569,36 +695,66 @@ void AliAnalysisTaskVnV0::Analyze(AliAODEvent* aodEvent, Float_t v0Centr) // re-map the container in an array to do the analysis for V0A and V0C within a loop Float_t evPlAngV0[2] = {evPlAngV0ACor2,evPlAngV0CCor2}; - AliCFContainer *contV0[2] = {fContAllChargesV0A,fContAllChargesV0C}; + AliFlowVZEROResults *contV0[2] = {fContAllChargesV0A,fContAllChargesV0C}; AliCFContainer *QA[2] = {fQA,fQA2}; Float_t evPlAngV0v3[2] = {evPlAngV0ACor3,evPlAngV0CCor3}; - AliCFContainer *contV0v3[2] = {fContAllChargesV0Av3,fContAllChargesV0Cv3}; + AliFlowVZEROResults *contV0v3[2] = {fContAllChargesV0Av3,fContAllChargesV0Cv3}; AliCFContainer *QAv3[2] = {fQAv3,fQA2v3}; + // Fill MC results + if(fIsMC && mcArray){ + fPID->ComputeProb(aodTrack,fAOD); // compute Bayesian probabilities + Float_t tofMismProbMC = fPID->GetTOFMismProb(); // TOF mismatch probability requested to be lower than 50% for TOF analysis + + Float_t xMC[5] = {iC,aodTrack->Charge(),1,evplaneMC,fPID->GetCurrentMask(1)&&tofMismProbMC < 0.5}; // to fill analysis v2 container + + + Float_t v2mc = TMath::Cos(2*(aodTrack->Phi() - evplaneMC)); + + fContAllChargesMC->Fill(0,aodTrack->Pt(),v2mc,xMC); + + Int_t iS = TMath::Abs(((AliAODMCParticle*)mcArray->At(TMath::Abs(aodTrack->GetLabel())))->GetPdgCode()); + if(iS==11){ + fContAllChargesMC->Fill(4,aodTrack->Pt(),v2mc,xMC); + } + else if(iS==13){ + fContAllChargesMC->Fill(5,aodTrack->Pt(),v2mc,xMC); + } + else if(iS==211){ + fContAllChargesMC->Fill(1,aodTrack->Pt(),v2mc,xMC); + } + else if(iS==321){ + fContAllChargesMC->Fill(2,aodTrack->Pt(),v2mc,xMC); + } + else if(iS==2212){ + fContAllChargesMC->Fill(3,aodTrack->Pt(),v2mc,xMC); + } + } + for(Int_t iV0=0;iV0<2;iV0++){ // loop on A and C side fPID->SetPsiCorrectionDeDx(evPlAngV0[iV0],evPlRes[iV0*8+iC]); // set the PID dE/dx correction as a function of the v2-EP (resolution is needed) - Double_t v2V0 = TMath::Cos(2*(aodTrack->Phi() - evPlAngV0[iV0])); - Double_t v3V0 = TMath::Cos(3*(aodTrack->Phi() - evPlAngV0v3[iV0])); + Float_t v2V0 = TMath::Cos(2*(aodTrack->Phi() - evPlAngV0[iV0])); + Float_t v3V0 = TMath::Cos(3*(aodTrack->Phi() - evPlAngV0v3[iV0])); - fPID->ComputeProb(aodTrack); // compute Bayesian probabilities + fPID->ComputeProb(aodTrack,fAOD); // compute Bayesian probabilities + Float_t dedx = fPID->GetDeDx();//aodTrack->GetTPCsignal(); Float_t *probRead = fPID->GetProb(); Float_t prob[8] = {probRead[0],probRead[1],probRead[2],probRead[3],probRead[4],probRead[5],probRead[6],probRead[7]}; Float_t tofMismProb = fPID->GetTOFMismProb(); // TOF mismatch probability requested to be lower than 50% for TOF analysis - Double_t x[7] = {iC,v2V0,aodTrack->Charge(),aodTrack->Pt(),1,evPlAngV0[iV0],fPID->GetCurrentMask(1)&&tofMismProb < 0.5}; // to fill analysis v2 container - Double_t x3[7] = {iC,v3V0,aodTrack->Charge(),aodTrack->Pt(),1,evPlAngV0v3[iV0],fPID->GetCurrentMask(1)&&tofMismProb < 0.5}; // to fill analysis v3 container + Float_t x[5] = {iC,aodTrack->Charge(),1,evPlAngV0[iV0],fPID->GetCurrentMask(1)&&tofMismProb < 0.5}; // to fill analysis v2 container + Float_t x3[5] = {iC,aodTrack->Charge(),1,evPlAngV0v3[iV0],fPID->GetCurrentMask(1)&&tofMismProb < 0.5}; // to fill analysis v3 container Double_t phi[5] = {iC,aodTrack->Charge(),aodTrack->Pt(),1,aodTrack->Phi()}; // to fill track container // Fill no PID - if(iV0 && fV2) fPhiTracks->Fill(phi,0); - if(fV2) contV0[iV0]->Fill(x,0); - if(fV3) contV0v3[iV0]->Fill(x3,0); - + if(iV0 && fQAsw) fPhiTracks->Fill(phi,0); + if(fV2) contV0[iV0]->Fill(0,aodTrack->Pt(),v2V0,x); + if(fV3) contV0v3[iV0]->Fill(0,aodTrack->Pt(),v3V0,x3); + - Float_t dedx = 0; Double_t dedxExp[8]; Float_t tof = -1; Double_t inttimes[8] = {0.,0.,0.,0.,0.,0.,0.,0.}; @@ -606,9 +762,7 @@ void AliAnalysisTaskVnV0::Analyze(AliAODEvent* aodEvent, Float_t v0Centr) if(aodTrack->GetDetPid()){ // check the PID object is available for(Int_t iS=0;iS < 8;iS++) dedxExp[iS] = fPID->GetExpDeDx(aodTrack,iS); - - dedx = aodTrack->GetTPCsignal(); - + if(fPID->GetCurrentMask(1)){ // if TOF is present Float_t ptrack = aodTrack->P(); tof = aodTrack->GetTOFsignal() - fPID->GetESDpid()->GetTOFResponse().GetStartTime(ptrack); @@ -641,169 +795,186 @@ void AliAnalysisTaskVnV0::Analyze(AliAODEvent* aodEvent, Float_t v0Centr) if(!(fPID->GetCurrentMask(0)) || !aodTrack->GetDetPid()){} // TPC PID and PID object strictly required (very important!!!!) else if(prob[2] > 0.6){ // pi phi[3] = prob[2]; // set probability in the container variables - x[4] = prob[2]; + x[2] = prob[2]; xQA[2] = prob[2]; - x3[4] = x[4]; + x3[2] = x[2]; xQA3[2] = xQA[2]; if(dedx > 10.){ // set TPC signal in the QA container variables xQA[3] = (dedx - dedxExp[2])/(dedxExp[2]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } - if(x[6] > 0.5){ // set TOF signal in the QA container variables if present + if(x[4] > 0.5){ // set TOF signal in the QA container variables if present xQA[4] = (tof - inttimes[2])/expTOFsigma[2]; // TOF xQA3[4] = xQA[4]; // TOF } if(TMath::Abs(xQA[3]) < 5){ // TPC 5 sigma extra cut to accept the track - if(iV0 && fV2) fPhiTracks->Fill(phi,1); - if(fV2) contV0[iV0]->Fill(x,1); + if(iV0 && fQAsw) fPhiTracks->Fill(phi,1); + if(fV2) contV0[iV0]->Fill(1,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(1,aodTrack->Pt(),v3V0,x3); if(fV2) QA[iV0]->Fill(xQA,0); - if(fV3) contV0v3[iV0]->Fill(x3,1); if(fV3) QAv3[iV0]->Fill(xQA3,0); } } else if(prob[3] > 0.6){ // K phi[3] = prob[3]; - x[4] = prob[3]; + x[2] = prob[3]; xQA[2] = prob[3]; - x3[4] = x[4]; + x3[2] = x[2]; xQA3[2] = xQA[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[3])/(dedxExp[3]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } - if(x[6] > 0.5){ + if(x[4] > 0.5){ xQA[4] = (tof - inttimes[3])/expTOFsigma[3]; // TOF xQA3[4] = xQA[4]; // TOF } if(TMath::Abs(xQA[3]) < 5){ - if(iV0 && fV2) fPhiTracks->Fill(phi,2); - if(fV2) contV0[iV0]->Fill(x,2); + if(iV0 && fQAsw) fPhiTracks->Fill(phi,2); + if(fV2) contV0[iV0]->Fill(2,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(2,aodTrack->Pt(),v3V0,x3); if(fV2) QA[iV0]->Fill(xQA,1); - if(fV3) contV0v3[iV0]->Fill(x3,2); if(fV3) QAv3[iV0]->Fill(xQA3,1); } } else if(prob[4] > 0.6){ // p phi[3] = prob[4]; - x[4] = prob[4]; + x[2] = prob[4]; xQA[2] = prob[4]; - x3[4] = x[4]; + x3[2] = x[2]; xQA3[2] = xQA[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[4])/(dedxExp[4]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } - if(x[6] > 0.5){ + if(x[4] > 0.5){ xQA[4] = (tof - inttimes[4])/expTOFsigma[4]; // TOF xQA3[4] = xQA[4]; // TOF } if(TMath::Abs(xQA[3]) < 5){ - if(iV0 && fV2) fPhiTracks->Fill(phi,3); - if(fV2) contV0[iV0]->Fill(x,3); + if(iV0 && fQAsw) fPhiTracks->Fill(phi,3); + if(fV2) contV0[iV0]->Fill(3,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(3,aodTrack->Pt(),v3V0,x3); if(fV2) QA[iV0]->Fill(xQA,2); - if(fV3) contV0v3[iV0]->Fill(x3,3); if(fV3) QAv3[iV0]->Fill(xQA3,2); } } else if(prob[0] > 0.6){ // e phi[3] = prob[0]; - x[4] = prob[0]; + x[2] = prob[0]; xQA[2] = prob[0]; - x3[4] = x[4]; + x3[2] = x[2]; xQA3[2] = xQA[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[0])/(dedxExp[0]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } - if(x[6] > 0.5){ + if(x[4] > 0.5){ xQA[4] = (tof - inttimes[0])/expTOFsigma[0]; // TOF xQA3[4] = xQA[4]; // TOF } if(TMath::Abs(xQA[3]) < 5){ - if(iV0 && fV2) fPhiTracks->Fill(phi,4); - if(fV2) contV0[iV0]->Fill(x,4); + if(iV0 && fQAsw) fPhiTracks->Fill(phi,4); + if(fV2) contV0[iV0]->Fill(4,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(4,aodTrack->Pt(),v3V0,x3); if(fV2) QA[iV0]->Fill(xQA,3); - if(fV3) contV0v3[iV0]->Fill(x3,4); if(fV3) QAv3[iV0]->Fill(xQA3,3); } } + else if(prob[1] > 0.6){ // mu + phi[3] = prob[1]; + x[2] = prob[1]; + xQA[2] = prob[1]; + x3[2] = x[2]; + xQA3[2] = xQA[2]; + if(dedx > 10.){ + xQA[3] = (dedx - dedxExp[1])/(dedxExp[1]*0.07); // TPC + xQA3[3] = xQA[3]; // TPC + } + if(x[4] > 0.5){ + xQA[4] = (tof - inttimes[1])/expTOFsigma[1]; // TOF + xQA3[4] = xQA[4]; // TOF + } + if(TMath::Abs(xQA[3]) < 5){ + if(fV2) contV0[iV0]->Fill(8,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(8,aodTrack->Pt(),v3V0,x3); + } + } else if(prob[5] > 0.6){ // d phi[3] = prob[5]; - x[4] = prob[5]; + x[2] = prob[5]; xQA[2] = prob[5]; - x3[4] = x[4]; + x3[2] = x[2]; xQA3[2] = xQA[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[5])/(dedxExp[5]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } - if(x[6] > 0.5){ + if(x[4] > 0.5){ xQA[4] = (tof - inttimes[5])/expTOFsigma[5]; // TOF xQA3[4] = xQA[4]; // TOF } if(TMath::Abs(xQA[3]) < 5){ - if(iV0 && fV2) fPhiTracks->Fill(phi,5); - if(fV2) contV0[iV0]->Fill(x,5); + if(iV0 && fQAsw) fPhiTracks->Fill(phi,5); + if(fV2) contV0[iV0]->Fill(5,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(5,aodTrack->Pt(),v3V0,x3); if(fV2) QA[iV0]->Fill(xQA,4); - if(fV3) contV0v3[iV0]->Fill(x3,5); if(fV3) QAv3[iV0]->Fill(xQA3,4); } } else if(prob[6] > 0.6){ // t phi[3] = prob[6]; - x[4] = prob[6]; + x[2] = prob[6]; xQA[2] = prob[6]; - x3[4] = x[4]; + x3[2] = x[2]; xQA3[2] = xQA[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[6])/(dedxExp[6]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } - if(x[6] > 0.5){ + if(x[4] > 0.5){ xQA[4] = (tof - inttimes[6])/expTOFsigma[6]; // TOF xQA3[4] = xQA[4]; // TOF } if(TMath::Abs(xQA[3]) < 5){ - if(iV0 && fV2) fPhiTracks->Fill(phi,6); - if(fV2) contV0[iV0]->Fill(x,6); + if(iV0 && fQAsw) fPhiTracks->Fill(phi,6); + if(fV2) contV0[iV0]->Fill(6,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(6,aodTrack->Pt(),v3V0,x3); if(fV2) QA[iV0]->Fill(xQA,5); - if(fV3) contV0v3[iV0]->Fill(x3,6); if(fV3) QAv3[iV0]->Fill(xQA3,5); } } else if(prob[7] > 0.6){ // He3 phi[3] = prob[7]; phi[1] *= 2; - x[4] = prob[7]; - x[3] *= 2; + x[2] = prob[7]; xQA[2] = prob[7]; - x3[3] = x[3]; - x3[4] = x[4]; + x3[2] = x[2]; xQA3[2] = xQA[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[7])/(dedxExp[7]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } - if(x[6] > 0.5){ + if(x[4] > 0.5){ xQA[4] = (tof - inttimes[7])/expTOFsigma[7]; // TOF xQA3[4] = xQA[4]; // TOF } if(TMath::Abs(xQA[3]) < 5){ - if(iV0 && fV2) fPhiTracks->Fill(phi,7); - if(fV2) contV0[iV0]->Fill(x,7); + if(iV0 && fQAsw) fPhiTracks->Fill(phi,7); + if(fV2) contV0[iV0]->Fill(7,aodTrack->Pt()*2,v2V0,x); + if(fV3) contV0[iV0]->Fill(7,aodTrack->Pt()*2,v3V0,x3); if(fV2) QA[iV0]->Fill(xQA,6); - if(fV3) contV0v3[iV0]->Fill(x3,7); if(fV3) QAv3[iV0]->Fill(xQA3,6); } - x[3] *= 0.5; - x3[3] = x[3]; + phi[1] *= 0.5; } - if(x[6] > 0.5){ // if TOF was present redo TPC stand alone PID to check the PID in the same acceptance (PID mask = 2) + if(x[4] > 0.5){ // if TOF was present redo TPC stand alone PID to check the PID in the same acceptance (PID mask = 2) fPID->ResetDetOR(1); // exclude TOF from PID tofMismProb = 0; - fPID->ComputeProb(aodTrack); + fPID->ComputeProb(aodTrack,fAOD); + dedx = fPID->GetDeDx();//aodTrack->GetTPCsignal(); probRead = fPID->GetProb(); fPID->SetDetOR(1); // include TOF for PID @@ -811,95 +982,103 @@ void AliAnalysisTaskVnV0::Analyze(AliAODEvent* aodEvent, Float_t v0Centr) Float_t probTPC[8] = {probRead[0],probRead[1],probRead[2],probRead[3],probRead[4],probRead[5],probRead[6],probRead[7]}; // TPC stand alone prbabilities //pid selection TPC S.A. with TOF matching - x[6]*=2; // set the mask to 2 id TOF is present - if(x[6]<1 || !(fPID->GetCurrentMask(0)) || !aodTrack->GetDetPid()){} // TPC PID S.A. PID in TOF acceptance + x[4]*=2; // set the mask to 2 id TOF is present + if(x[4]<1 || !(fPID->GetCurrentMask(0)) || !aodTrack->GetDetPid()){} // TPC PID S.A. PID in TOF acceptance else if(probTPC[2] > 0.6){ // pi - x[4] = probTPC[2]; - x3[4] = x[4]; + x[2] = probTPC[2]; + x3[2] = x[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[2])/(dedxExp[2]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } if(TMath::Abs(xQA[3]) < 5){ - if(fV2) contV0[iV0]->Fill(x,1); - if(fV3) contV0v3[iV0]->Fill(x3,1); + if(fV2) contV0[iV0]->Fill(1,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(1,aodTrack->Pt(),v3V0,x3); } } else if(probTPC[3] > 0.6){ // K - x[4] = probTPC[3]; - x3[4] = x[4]; + x[2] = probTPC[3]; + x3[2] = x[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[3])/(dedxExp[3]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } if(TMath::Abs(xQA[3]) < 5){ - if(fV2) contV0[iV0]->Fill(x,2); - if(fV3) contV0v3[iV0]->Fill(x3,2); + if(fV2) contV0[iV0]->Fill(2,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(2,aodTrack->Pt(),v3V0,x3); } } else if(probTPC[4] > 0.6){ // p - x[4] = probTPC[4]; - x3[4] = x[4]; + x[2] = probTPC[4]; + x3[2] = x[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[4])/(dedxExp[4]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } if(TMath::Abs(xQA[3]) < 5){ - if(fV2) contV0[iV0]->Fill(x,3); - if(fV3) contV0v3[iV0]->Fill(x3,3); + if(fV2) contV0[iV0]->Fill(3,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(3,aodTrack->Pt(),v3V0,x3); } } else if(probTPC[0] > 0.6){ // e - x[4] = probTPC[0]; - x3[4] = x[4]; + x[2] = probTPC[0]; + x3[2] = x[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[0])/(dedxExp[0]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } if(TMath::Abs(xQA[3]) < 5){ - if(fV2) contV0[iV0]->Fill(x,4); - if(fV3) contV0v3[iV0]->Fill(x3,4); + if(fV2) contV0[iV0]->Fill(4,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(4,aodTrack->Pt(),v3V0,x3); + } + } + else if(probTPC[1] > 0.6){ // mu + x[2] = probTPC[1]; + x3[2] = x[2]; + if(dedx > 10.){ + xQA[3] = (dedx - dedxExp[1])/(dedxExp[1]*0.07); // TPC + xQA3[3] = xQA[3]; // TPC + } + if(TMath::Abs(xQA[3]) < 5){ + if(fV2) contV0[iV0]->Fill(8,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(8,aodTrack->Pt(),v3V0,x3); } } else if(probTPC[5] > 0.6){ // d - x[4] = probTPC[5]; - x3[4] = x[4]; + x[2] = probTPC[5]; + x3[2] = x[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[5])/(dedxExp[5]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } if(TMath::Abs(xQA[3]) < 5){ - if(fV2) contV0[iV0]->Fill(x,5); - if(fV3) contV0v3[iV0]->Fill(x3,5); + if(fV2) contV0[iV0]->Fill(5,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(5,aodTrack->Pt(),v3V0,x3); } } else if(probTPC[6] > 0.6){ // t - x[4] = probTPC[6]; - x3[4] = x[4]; + x[2] = probTPC[6]; + x3[2] = x[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[6])/(dedxExp[6]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } if(TMath::Abs(xQA[3]) < 5){ - if(fV2) contV0[iV0]->Fill(x,6); - if(fV3) contV0v3[iV0]->Fill(x3,6); + if(fV2) contV0[iV0]->Fill(6,aodTrack->Pt(),v2V0,x); + if(fV3) contV0[iV0]->Fill(6,aodTrack->Pt(),v3V0,x3); } } else if(probTPC[7] > 0.6){ // He3 - x[4] = probTPC[7]; - x3[4] = x[4]; - x[3] *= 2; - x3[3] = x[3]; + x[2] = probTPC[7]; + x3[2] = x[2]; if(dedx > 10.){ xQA[3] = (dedx - dedxExp[7])/(dedxExp[7]*0.07); // TPC xQA3[3] = xQA[3]; // TPC } if(TMath::Abs(xQA[3]) < 5){ - if(fV2) contV0[iV0]->Fill(x,7); - if(fV3) contV0v3[iV0]->Fill(x3,7); + if(fV2) contV0[iV0]->Fill(7,aodTrack->Pt()*2,v2V0,x); + if(fV3) contV0[iV0]->Fill(7,aodTrack->Pt()*2,v3V0,x3); } - x[3] *= 0.5; - x3[3] = x[3]; } } // end side loop } // end track loop diff --git a/PWG/FLOW/Tasks/AliAnalysisTaskVnV0.h b/PWG/FLOW/Tasks/AliAnalysisTaskVnV0.h index 80dab1d9253..3e7243c2079 100644 --- a/PWG/FLOW/Tasks/AliAnalysisTaskVnV0.h +++ b/PWG/FLOW/Tasks/AliAnalysisTaskVnV0.h @@ -13,6 +13,7 @@ #include #include #include "AliFlowBayesianPID.h" +#include "AliFlowVZEROResults.h" class TH2F; @@ -37,11 +38,14 @@ class AliAnalysisTaskVnV0 : public AliAnalysisTaskSE { virtual void SetV2(Bool_t val){fV2 = val;}; virtual void SetV3(Bool_t val){fV3 = val;}; + virtual void SetMC(Bool_t flag = kTRUE){fIsMC = flag;}; + virtual void SetQA(Bool_t flag = kTRUE){fQAsw = flag;}; + void OpenInfoCalbration(Int_t run); private: - AliAnalysisTaskVnV0(const AliAnalysisTaskVnV0 &); - AliAnalysisTaskVnV0 & operator=(const AliAnalysisTaskVnV0 &); + AliAnalysisTaskVnV0(const AliAnalysisTaskVnV0 &old); + AliAnalysisTaskVnV0& operator=(const AliAnalysisTaskVnV0 &source); virtual Float_t GetVertex(AliAODEvent* aod) const; virtual void Analyze(AliAODEvent* aodEvent, Float_t v0Centr); @@ -59,7 +63,7 @@ class AliAnalysisTaskVnV0 : public AliAnalysisTaskSE { Int_t fRun; // current run checked to load VZERO calibrations - TList *fList,*fList2; // List for output objects + TList *fList,*fList2,*fList3,*fList4; // List for output objects // // Output objects TProfile *fMultV0; // object containing VZERO calibration information @@ -69,11 +73,6 @@ class AliAnalysisTaskVnV0 : public AliAnalysisTaskSE { Float_t fMeanQv3[nCentrBin][2][2]; // also for v3 Float_t fWidthQv3[nCentrBin][2][2]; // ... - AliCFContainer *fContAllChargesV0A; // containers for v2 (A nd C) and v3 (A and C) - AliCFContainer *fContAllChargesV0C; // ... - AliCFContainer *fContAllChargesV0Av3; // ... - AliCFContainer *fContAllChargesV0Cv3; // ... - TProfile *fHResTPCv0A2,*fHResTPCv0C2,*fHResv0Cv0A2; // TProfile for subevent resolution (output) TProfile *fHResTPCv0A3,*fHResTPCv0C3,*fHResv0Cv0A3; // also for v3 @@ -96,7 +95,12 @@ class AliAnalysisTaskVnV0 : public AliAnalysisTaskSE { Bool_t fV2,fV3; // swith to set the armonics - ClassDef(AliAnalysisTaskVnV0, 1); //Analysis task v2 and v3 analysis on AOD + AliFlowVZEROResults *fContAllChargesV0A,*fContAllChargesV0C,*fContAllChargesV0Av3,*fContAllChargesV0Cv3,*fContAllChargesMC,*fContAllChargesMCv3; + + Bool_t fIsMC; // if MC + Bool_t fQAsw; // if QA + + ClassDef(AliAnalysisTaskVnV0, 2); //Analysis task v2 and v3 analysis on AOD }; #endif diff --git a/PWG/FLOW/Tasks/AliFlowBayesianPID.cxx b/PWG/FLOW/Tasks/AliFlowBayesianPID.cxx index 1d023c23599..b7c24d2dcbf 100644 --- a/PWG/FLOW/Tasks/AliFlowBayesianPID.cxx +++ b/PWG/FLOW/Tasks/AliFlowBayesianPID.cxx @@ -34,7 +34,7 @@ AliTOFGeometry* AliFlowBayesianPID::fgTofGeo = NULL; // TOF geometry needed to r //________________________________________________________________________ AliFlowBayesianPID::AliFlowBayesianPID(AliESDpid *esdpid) - : AliPIDResponse(), fPIDesd(NULL), fDB(TDatabasePDG::Instance()), fNewTrackParam(0), fTOFresolution(84.0), fTOFResponseF(NULL), fTPCResponseF(NULL), fTOFmaker(NULL),fWTofMism(0.0), fProbTofMism(0.0), fZ(0) ,fMassTOF(0), fBBdata(NULL),fCurrCentrality(100),fPsi(999),fPsiRes(999),fIsMC(kFALSE) + : AliPIDResponse(), fPIDesd(NULL), fDB(TDatabasePDG::Instance()), fNewTrackParam(0), fTOFresolution(84.0), fTOFResponseF(NULL), fTPCResponseF(NULL), fTOFmaker(NULL),fWTofMism(0.0), fProbTofMism(0.0), fZ(0) ,fMassTOF(0), fBBdata(NULL),fCurrCentrality(100),fPsi(999),fPsiRes(999),fIsMC(kFALSE),fDedx(0.0) { // Constructor Bool_t redopriors = kFALSE; @@ -534,6 +534,8 @@ void AliFlowBayesianPID::ComputeWeights(const AliESDtrack *t){ } } + fDedx = dedx; + if(t->GetStatus() & AliESDtrack::kTPCout && dedx > 40 && fMaskOR[0]){ // if TPC PID available for(Int_t iS=0;iSSetV2(kV2); task->SetV3(kV3); + if(ismc) task->SetMC(); + if(qa) task->SetQA(); mgr->AddTask(task); @@ -39,7 +45,16 @@ AliAnalysisTask *AddTaskVZERO(AliAnalysisManager *mgr,Bool_t kV2=kTRUE,Bool_t kV AliAnalysisDataContainer *cOutputL2= mgr->CreateContainer("contVZEROv3",TList::Class(), AliAnalysisManager::kOutputContainer, fileout2); mgr->ConnectOutput(task, 2, cOutputL2); } + if(ismc){ + AliAnalysisDataContainer *cOutputL3= mgr->CreateContainer("contVZEROmc",TList::Class(), AliAnalysisManager::kOutputContainer, fileout3); + mgr->ConnectOutput(task, 3, cOutputL3); + } + if(qa){ + AliAnalysisDataContainer *cOutputL4= mgr->CreateContainer("contVZEROqa",TList::Class(), AliAnalysisManager::kOutputContainer, fileout4); + mgr->ConnectOutput(task, 4, cOutputL4); + } printf("task really added\n"); return task; } + diff --git a/PWGCF/FLOW/macros/extractFlowVZERO.C b/PWGCF/FLOW/macros/extractFlowVZERO.C index 94ec7c179df..3257afd7d3d 100644 --- a/PWGCF/FLOW/macros/extractFlowVZERO.C +++ b/PWGCF/FLOW/macros/extractFlowVZERO.C @@ -1,242 +1,54 @@ -extractFlowVZERO(Int_t itech=0,Int_t ibin = 3,Int_t step = 0,Int_t arm=2,Float_t pTh = 0.9,Int_t charge=0,Int_t addbins=0){ - TCanvas *can1 = new TCanvas("cV0Acheck","cV0Acheck"); - -/* -read outVZEROv2.root or outVZEROv3.root files - -itech=0 TPC stand alone PID and selection -itech=1 TPC stand alone PID but TOF track selection -itech=2 TPC&TOF PID (TOF strictly required) -itech=3 TPC|TOF PID (the same of before but using TPC stand alone when TOF not available) - -ibin = centrality bin - -step = species (0=all charges, 1=pi, 2=K, 3=p, 4=e, 5=d, 6=t, 7=3He - -arm = armonic 2 (v2) or 3 (v3) - -pTh = probability threshold = 0.6 , 0.8 , 0.9 , 0.95 - -charge = 1(pos), -1(neg), 0(both) - -addbins = to merge more centrality bin (i.e. ibin = 2, addbins = 3 merge 2,3,4,5 centrality bins) - preliminary version!!!!!! - -*/ - - // load libraries - gSystem->Load("libANALYSIS"); - gSystem->Load("libANALYSISalice"); - gSystem->Load("libCORRFW.so"); - - // get objects - char filein[100]; - sprintf(filein,"outVZEROv%i.root",arm); - TFile *f=new TFile(filein); - sprintf(filein,"contVZEROv%i",arm); - TList *l = f->Get(filein); - l->ls(); - AliCFContainer *c = l->At(0); - AliCFContainer *c3 = l->At(1); - TProfile *p1 = l->At(2); - TProfile *p2 = l->At(3); - TProfile *p3 = l->At(4); - Int_t iVar[2] = {1,3}; - - Int_t iVarPsi[1] = {5}; - /* 0=centrality , 1= cos(arm*deltaphi) , 2=charge , 3=pt , 4=probability , 5=EP , 6=PIDmask*/ - Double_t iMin[7] = {ibin,-1,-1.5,0,pTh+0.00001,-3.14,0}; - Double_t iMax[7] = {ibin+addbins,1,1.5,20,1.1,3.14,1}; - - char *iTechName[4] = {"tpc","tpcInTof","tof","tpcORtof"}; - - Bool_t kOnlyTPC = kFALSE; - - if(itech==0 && (step!=0)){ - kOnlyTPC = kTRUE; - } - else if(itech==1 && (step!=0)){ - iMin[6] = 2; - iMax[6] = 2; - } - else if(itech==2 || itech==1){ - iMin[6] = 1; - iMax[6] = 1; - } - else{ - iMin[6] = 0; - iMax[6] = 1; - } +extractFlowVZERO(Int_t icentr,Int_t spec,Int_t arm=2,Bool_t isMC=kFALSE){ + // NUA correction currently are missing + char name[100]; + sprintf(name,"outVZEROv%i.root",arm); + TFile *fo = new TFile(name); + sprintf(name,"contVZEROv%i",arm); + TList *cont = (TList *) fo->Get(name); - if(charge==-1) iMax[2]=0; - if(charge==1) iMin[2]=0; + Float_t xMin[5] = {icentr,-1,0,-10,0}; + Float_t xMax[5] = {icentr,1,1,10,1.5}; - // EP distribution needed for NUA corrections - TH2F *hPhiA = l->At(5); - TH2F *hPhiC = l->At(6); + cont->ls(); - TH1D *hPhiPA = hPhiA->ProjectionY("_py",iMin[0]+1,iMax[0]+1); - TH1D *hPhiPC = hPhiC->ProjectionY("_py",iMin[0]+1,iMax[0]+1); + TProfile *p1 = cont->At(2); + TProfile *p2 = cont->At(3); + TProfile *p3 = cont->At(4); - // EP resolution variables Float_t res1=0,res2=0,res3=0; Float_t eres1=0,eres2=0,eres3=0; - for(Int_t i=iMin[0];i<=iMax[0];i++){ // in case of more centrality bins weighting with the errors - if(p1->GetBinError(i+1)){ - eres1 += 1./p1->GetBinError(i+1)/p1->GetBinError(i+1); - res1 += p1->GetBinContent(i+1)/p1->GetBinError(i+1)/p1->GetBinError(i+1); - } - if(p2->GetBinError(i+1)){ - eres2 += 1./p2->GetBinError(i+1)/p2->GetBinError(i+1); - res2 += p2->GetBinContent(i+1)/p2->GetBinError(i+1)/p2->GetBinError(i+1); - } - if(p3->GetBinError(i+1)){ - eres3 += 1./p3->GetBinError(i+1)/p3->GetBinError(i+1); - res3 += p3->GetBinContent(i+1)/p3->GetBinError(i+1)/p3->GetBinError(i+1); - } - - if(eres1) res1 /= eres1; - if(eres2) res2 /= eres2; - if(eres3) res3 /= eres3; - - if(eres1) eres1 = TMath::Sqrt(1./eres1); - if(eres1) eres2 = TMath::Sqrt(1./eres2); - if(eres1) eres3 = TMath::Sqrt(1./eres3); + Int_t i = icentr; + if(p1->GetBinError(i+1)){ + eres1 += 1./p1->GetBinError(i+1)/p1->GetBinError(i+1); + res1 += p1->GetBinContent(i+1)/p1->GetBinError(i+1)/p1->GetBinError(i+1); } - - // NUA correction (fit to EP distribution) - TF1 *fpol = new TF1("fPol","pol0"); - hPhiPA->Fit("fPol","","",-TMath::Pi()/arm,TMath::Pi()/arm); - Float_t scalA = fPol->GetParameter(0); - hPhiPC->Fit("fPol","","",-TMath::Pi()/arm,TMath::Pi()/arm); - Float_t scalC = fPol->GetParameter(0); - - AliCFContainer *c2[20]; - AliCFContainer *c4[20]; - TH2D *h[20],*h2[20]; - - AliCFContainer *c2bis[20]; - AliCFContainer *c4bis[20]; - - AliCFContainer *cPsi2[20]; - AliCFContainer *cPsi4[20]; - TH1D *hPsi[20],*hPsi2[20]; - AliCFContainer *cPsi2bis[20]; - AliCFContainer *cPsi4bis[20]; - - Float_t intA = 0; - Float_t intC = 0; - - for(Int_t i=5;i<15;i++){ - printf("%i\n",i); - iMin[5] = hPhiPA->GetBinCenter(i+1); - iMax[5] = hPhiPA->GetBinCenter(i+1); - if(!kOnlyTPC){ - c2[i] = c->MakeSlice(2,iVar,iMin,iMax); - c4[i]= c3->MakeSlice(2,iVar,iMin,iMax); + if(p2->GetBinError(i+1)){ + eres2 += 1./p2->GetBinError(i+1)/p2->GetBinError(i+1); + res2 += p2->GetBinContent(i+1)/p2->GetBinError(i+1)/p2->GetBinError(i+1); } - else{ // merge to maskPID bins (TPC standalone without TOF + TPC stand alone with TOF) - iMin[6] = 0; - iMax[6] = 0; - c2[i] = c->MakeSlice(2,iVar,iMin,iMax); - c4[i]= c3->MakeSlice(2,iVar,iMin,iMax); - - iMin[6] = 2; - iMax[6] = 2; - c2bis[i] = c->MakeSlice(2,iVar,iMin,iMax); - c4bis[i]= c3->MakeSlice(2,iVar,iMin,iMax); - } - - - if(!kOnlyTPC){ - cPsi2[i] = c->MakeSlice(1,iVarPsi,iMin,iMax); - cPsi4[i]= c3->MakeSlice(1,iVarPsi,iMin,iMax); - } - else{ // merge to maskPID bins (TPC standalone without TOF + TPC stand alone with TOF) - iMin[6] = 0; - iMax[6] = 0; - cPsi2[i] = c->MakeSlice(1,iVarPsi,iMin,iMax); - cPsi4[i]= c3->MakeSlice(1,iVarPsi,iMin,iMax); - - iMin[6] = 2; - iMax[6] = 2; - cPsi2bis[i] = c->MakeSlice(1,iVarPsi,iMin,iMax); - cPsi4bis[i]= c3->MakeSlice(1,iVarPsi,iMin,iMax); - } - - h[i] = (TH2D *) c2[i]->Project(step,0,1); - h2[i] = (TH2D *) c4[i]->Project(step,0,1); - if(kOnlyTPC){ - TH2D *htemp = (TH2D *) c2bis[i]->Project(step,0,1); - h[i]->Add(htemp); - htemp = (TH2D *) c4bis[i]->Project(step,0,1); - h2[i]->Add(htemp); - } - - if(hPhiPA->GetBinContent(i)) h[i]->Scale(scalA/hPhiPA->GetBinContent(i+1)); // reweighting for NUA correction - if(hPhiPC->GetBinContent(i)) h2[i]->Scale(scalC/hPhiPC->GetBinContent(i+1)); // reweighting for NUA correction - - hPsi[i] = (TH1D *) cPsi2[i]->Project(step,0); - hPsi2[i] = (TH1D *) cPsi4[i]->Project(step,0); - - if(kOnlyTPC){ // merge to maskPID bins (TPC standalone without TOF + TPC stand alone with TOF) - TH1D *htemp2 = (TH1D *) cPsi2bis[i]->Project(step,0); - hPsi[i]->Add(htemp2); - htemp2 = (TH1D *) cPsi4bis[i]->Project(step,0); - hPsi2[i]->Add(htemp2); - } - - if(hPhiPA->GetBinContent(i+1)) hPsi[i]->Scale(scalA/hPhiPA->GetBinContent(i+1)); // check of reweighting for NUA correction - if(hPhiPC->GetBinContent(i+1)) hPsi2[i]->Scale(scalC/hPhiPC->GetBinContent(i+1)); // check of reweighting for NUA correction - - intA += hPsi[i]->Integral(); - intC += hPsi2[i]->Integral(); - + if(p3->GetBinError(i+1)){ + eres3 += 1./p3->GetBinError(i+1)/p3->GetBinError(i+1); + res3 += p3->GetBinContent(i+1)/p3->GetBinError(i+1)/p3->GetBinError(i+1); } - -// hPhiPA->Scale(intA / (scalA)/10); -// hPhiPC->Scale(intC / (scalC)/10); - - // NUA correction check V0A - hPhiPA->Draw(); - hPsi[5]->Draw("SAME"); - hPsi[5]->Scale(scalA / intA * 10); - for(Int_t i=6;i<15;i++){ - h[5]->Add(h[i]); - h2[5]->Add(h2[i]); - hPsi[i]->Draw("SAME"); - hPsi[i]->Scale(scalA / intA * 10); - } - - // NUA correction check V0C - TCanvas *can2 = new TCanvas("cV0Ccheck","cV0Ccheck"); - hPhiPC->Draw(); - hPsi2[5]->Draw("SAME"); - for(Int_t i=6;i<15;i++){ - h[5]->Add(h[i]); - h2[5]->Add(h2[i]); - hPsi2[i]->Draw("SAME"); - hPsi2[i]->Scale(scalC /intC * 10); - } - hPsi2[5]->Scale(scalC /intC * 10); - - // Flow for V0A and V0C separately - TCanvas *can3 = new TCanvas("cFlowComp","cFlowComp"); - TProfile *pp = h[5]->ProfileY(); - pp->Draw(); - TProfile *pp2 = h2[5]->ProfileY(); - pp2->Draw("SAME"); - - printf("nev (selected within 0-80%) = %i\n",p1->GetEntries()); - - - // correction for resoltion + + res1 /= eres1; + res2 /= eres2; + res3 /= eres3; + + AliFlowVZEROResults *a = (AliFlowVZEROResults *) cont->At(0); + AliFlowVZEROResults *b = (AliFlowVZEROResults *) cont->At(1); + TProfile *pp = a->GetV2(spec,xMin,xMax); + TProfile *pp2 = b->GetV2(spec,xMin,xMax); + + Float_t scaling = sqrt(res1*res3/res2); pp->Scale(1./scaling); + printf("resolution V0A = %f\n",scaling); Float_t err1_2 = eres1*eres1/res1/res1/4 + - eres2*eres2/res2/res2/4 + - eres3*eres3/res3/res3/4; + eres2*eres2/res2/res2/4 + + eres3*eres3/res3/res3/4; Float_t err2_2 = err1_2; err1_2 /= scaling*scaling; scaling = sqrt(res2*res3/res1); @@ -244,36 +56,18 @@ addbins = to merge more centrality bin (i.e. ibin = 2, addbins = 3 merge 2,3,4,5 pp2->Scale(1./scaling); printf("resolution V0C =%f\n",scaling); - char title[100]; - sprintf(title,"VZERO EP;p_{t} (GeV/c);v_{%i}",arm); - pp->SetTitle(title); + pp->SetName("V0A"); + pp2->SetName("V0C"); - // Average V0A-V0C - TH1D *pAll = pp->ProjectionX(); + pp->Draw(); + pp2->Draw("SAME"); - for(Int_t i=1;i <= pAll->GetNbinsX();i++){ - Float_t e1 = err1_2*pp->GetBinContent(i)*pp->GetBinContent(i) + pp->GetBinError(i)*pp->GetBinError(i); - Float_t e2 = err2_2*pp2->GetBinContent(i)*pp2->GetBinContent(i) + pp2->GetBinError(i)*pp2->GetBinError(i); - Float_t xval = 0,exval = 0; - if(e1 >0 && e2>0){ - xval = (pp->GetBinContent(i)/e1 + pp2->GetBinContent(i)/e2)/(1/e1 + 1/e2); - exval = 1./sqrt(1/e1 + 1/e2); - } - pAll->SetBinContent(i,xval); - pAll->SetBinError(i,exval); + if(arm == 2 && isMC){ + sprintf(name,"outVZEROmc.root"); + fo = new TFile(name); + sprintf(name,"contVZEROmc"); + cont = (TList *) fo->Get(name); + AliFlowVZEROResults *c = (AliFlowVZEROResults *) cont->At(0); + c->GetV2(spec,xMin,xMax)->Draw("SAME"); } - // combined measurement - TCanvas *can4 = new TCanvas("cFlowCombined","cFlowCombined"); - pAll->Draw(); - - char name[100]; - sprintf(name,"out%i-%i_%i_%4.2f_%i%sv%i.root",iMin[0],iMax[0],step,pTh,charge,iTechName[itech],arm); - TFile *fout = new TFile(name,"RECREATE"); - - pAll->SetName("histo"); - pAll->Write(); - can1->Write(); - can2->Write(); - can3->Write(); - fout->Close(); } -- 2.43.0