Adding the track-by-track correction for both the trigger and the associated particle...

diff --git a/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskBFPsi.cxx b/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskBFPsi.cxx
index 5a2b626edeeccba27ba35c3971d9e58d7ca79d5e..66ac8a1f78c4644b79891acc9df9133266d4ce06 100755 (executable)
--- a/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskBFPsi.cxx
+++ b/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskBFPsi.cxx
@@ -99,6 +99,8 @@ AliAnalysisTaskBFPsi::AliAnalysisTaskBFPsi(const char *name)
   fHistProbTPCTOFvsPtafterPID(NULL),\r
   fHistNSigmaTPCvsPtafterPID(NULL), \r
   fHistNSigmaTOFvsPtafterPID(NULL),  \r
+  //fHistMatrixCorrectionPlus(0), //======================================================correction\r
+  //fHistMatrixCorrectionMinus(0), //=====================================================correction\r
   fPIDResponse(0x0),\r
   fPIDCombined(0x0),\r
   fParticleOfInterest(kPion),\r
@@ -126,8 +128,19 @@ AliAnalysisTaskBFPsi::AliAnalysisTaskBFPsi(const char *name)
   nAODtrackCutBit(128),\r
   fPtMin(0.3),\r
   fPtMax(1.5),\r
+  fPtMinForCorrections(0.3),\r
+  fPtMaxForCorrections(1.5),\r
+  fPtBinForCorrections(36), //=================================correction\r
   fEtaMin(-0.8),\r
   fEtaMax(-0.8),\r
+  fEtaMinForCorrections(-0.8),\r
+  fEtaMaxForCorrections(-0.8),\r
+  fEtaBinForCorrections(36), //=================================correction\r
+  fPhiMin(0.),//=================================correction\r
+  fPhiMax(360.),//=================================correction\r
+  fPhiMinForCorrections(0.),//=================================correction\r
+  fPhiMaxForCorrections(360.),//=================================correction\r
+  fPhiBinForCorrections(100), //=================================correction\r
   fDCAxyCut(-1),\r
   fDCAzCut(-1),\r
   fTPCchi2Cut(-1),\r
@@ -138,10 +151,15 @@ AliAnalysisTaskBFPsi::AliAnalysisTaskBFPsi(const char *name)
   fExcludeResonancesInMC(kFALSE),\r
   fUseMCPdgCode(kFALSE),\r
   fPDGCodeToBeAnalyzed(-1),\r
-  fEventClass("EventPlane") {\r
+  fEventClass("EventPlane") \r
+{\r
   // Constructor\r
   // Define input and output slots here\r
   // Input slot #0 works with a TChain\r
+  for (Int_t i=0; i<=kCENTRALITY; i++){\r
+    fHistMatrixCorrectionPlus[i] = NULL; //======================================================correction\r
+    fHistMatrixCorrectionMinus[i] = NULL; //=====================================================correction\r
+  }\r
   DefineInput(0, TChain::Class());\r
   // Output slot #0 writes into a TH1 container\r
   DefineOutput(1, TList::Class());\r
@@ -450,6 +468,76 @@ void AliAnalysisTaskBFPsi::UserCreateOutputObjects() {
 \r
 }\r
 \r
+//====================================correction=============================================//\r
+//________________________________________________________________________\r
+void AliAnalysisTaskBFPsi::SetInputCorrection(const char* filename, const char* gCollSystem) {\r
+  //Open files that will be used for correction\r
+  TString gCollidingSystem = gCollSystem;\r
+\r
+  //Open the input file\r
+  TFile *f = TFile::Open(filename);\r
+  if(!f->IsOpen()) {\r
+    Printf("File not found!!!");\r
+     return;\r
+    //break;\r
+  }\r
+  \r
+  //Get the TDirectoryFile and TList\r
+  TDirectoryFile *dir = dynamic_cast<TDirectoryFile *>(f->Get("PWGCFEbyE.outputBalanceFunctionEfficiencyAnalysis"));\r
+  if(!dir) {\r
+    Printf("TDirectoryFile not found!!!");\r
+     return;\r
+    //break;\r
+  }\r
+  \r
+  TString listEffName = "";\r
+  for (Int_t iCent = 0; iCent < kCENTRALITY; iCent++) {\r
+    listEffName = "listEfficiencyBF_";  \r
+    listEffName += (TString)((Int_t)(centralityArrayForPbPb[iCent]));\r
+    listEffName += "_";\r
+    listEffName += (TString)((Int_t)(centralityArrayForPbPb[iCent+1]));\r
+    TList *list = (TList*)dir->Get(listEffName.Data());\r
+    if(!list) {\r
+      Printf("TList Efficiency not found!!!");\r
+      return;\r
+      // break;\r
+    }\r
+    \r
+    TString histoName = "fHistMatrixCorrectionPlus";\r
+    fHistMatrixCorrectionPlus[iCent]= dynamic_cast<TH3D *>(list->FindObject(histoName.Data()));\r
+    //fHistMatrixCorrectionPlus[iCent]->SetName(Form(histoName.Data()));\r
+    // histoName.Data() = fHistMatrixCorrectionPlus->Clone;\r
+    if(!fHistMatrixCorrectionPlus[iCent]) {\r
+      Printf("fHist not found!!!");\r
+      return;\r
+      //break;\r
+    }\r
+    \r
+    histoName = "fHistMatrixCorrectionMinus";\r
+    fHistMatrixCorrectionMinus[iCent] = dynamic_cast<TH3D *>(list->FindObject(histoName.Data()));  \r
+    //fHistMatrixCorrectionMinus[iCent]->SetName(Form(histoName.Data()));\r
+    //histoName.Data() = fHistMatrixCorrectionPlus->Clone;\r
+    if(!fHistMatrixCorrectionMinus[iCent]) {\r
+      Printf("fHist not found!!!");\r
+      return; \r
+      //break;\r
+    }\r
+  }//loop over centralities: ONLY the PbPb case is covered\r
+\r
+  fEtaMinForCorrections = fHistMatrixCorrectionPlus[0]->GetXaxis()->GetXmin();\r
+  fEtaMaxForCorrections = fHistMatrixCorrectionPlus[0]->GetXaxis()->GetXmax();\r
+  fEtaBinForCorrections = fHistMatrixCorrectionPlus[0]->GetNbinsX();\r
+\r
+  fPtMinForCorrections = fHistMatrixCorrectionPlus[0]->GetYaxis()->GetXmin();\r
+  fPtMaxForCorrections = fHistMatrixCorrectionPlus[0]->GetYaxis()->GetXmax();\r
+  fPtBinForCorrections = fHistMatrixCorrectionPlus[0]->GetNbinsY();\r
+\r
+  fPhiMinForCorrections = fHistMatrixCorrectionPlus[0]->GetZaxis()->GetXmin();\r
+  fPhiMaxForCorrections = fHistMatrixCorrectionPlus[0]->GetZaxis()->GetXmax();\r
+  fPhiBinForCorrections = fHistMatrixCorrectionPlus[0]->GetNbinsZ();\r
+}\r
+//====================================correction=============================================//\r
+\r
 //________________________________________________________________________\r
 void AliAnalysisTaskBFPsi::UserExec(Option_t *) {\r
   // Main loop\r
@@ -457,10 +545,10 @@ void AliAnalysisTaskBFPsi::UserExec(Option_t *) {
 \r
   TString gAnalysisLevel = fBalance->GetAnalysisLevel();\r
   Int_t gNumberOfAcceptedTracks = 0;\r
-  Double_t fCentrality          = -1.;\r
+  Double_t gCentrality          = -1.;\r
   Double_t gReactionPlane       = -1.; \r
   Float_t bSign = 0.;\r
-\r
+  \r
   // get the event (for generator level: MCEvent())\r
   AliVEvent* eventMain = NULL;\r
   if(gAnalysisLevel == "MC") {\r
@@ -484,7 +572,7 @@ void AliAnalysisTaskBFPsi::UserExec(Option_t *) {
   }\r
   \r
   // check event cuts and fill event histograms\r
-  if((fCentrality = IsEventAccepted(eventMain)) < 0){\r
+  if((gCentrality = IsEventAccepted(eventMain)) < 0){\r
     return;\r
   }\r
   \r
@@ -493,17 +581,17 @@ void AliAnalysisTaskBFPsi::UserExec(Option_t *) {
 \r
   // get the reaction plane\r
   gReactionPlane = GetEventPlane(eventMain);\r
-  fHistEventPlane->Fill(gReactionPlane,fCentrality);\r
+  fHistEventPlane->Fill(gReactionPlane,gCentrality);\r
   if(gReactionPlane < 0){\r
     return;\r
   }\r
   \r
   // get the accepted tracks in main event\r
-  TObjArray *tracksMain = GetAcceptedTracks(eventMain,fCentrality,gReactionPlane);\r
+  TObjArray *tracksMain = GetAcceptedTracks(eventMain,gCentrality,gReactionPlane);\r
   gNumberOfAcceptedTracks = tracksMain->GetEntriesFast();\r
 \r
   //multiplicity cut (used in pp)\r
-  fHistNumberOfAcceptedTracks->Fill(gNumberOfAcceptedTracks,fCentrality);\r
+  fHistNumberOfAcceptedTracks->Fill(gNumberOfAcceptedTracks,gCentrality);\r
   if(fUseMultiplicity) {\r
     if((gNumberOfAcceptedTracks < fNumberOfAcceptedTracksMin)||(gNumberOfAcceptedTracks > fNumberOfAcceptedTracksMax))\r
       return;\r
@@ -512,7 +600,7 @@ void AliAnalysisTaskBFPsi::UserExec(Option_t *) {
   // store charges of all accepted tracks, shuffle and reassign (two extra loops!)\r
   TObjArray* tracksShuffled = NULL;\r
   if(fRunShuffling){\r
-    tracksShuffled = GetShuffledTracks(tracksMain);\r
+    tracksShuffled = GetShuffledTracks(tracksMain,gCentrality);\r
   }\r
   \r
   // Event mixing \r
@@ -534,10 +622,10 @@ void AliAnalysisTaskBFPsi::UserExec(Option_t *) {
       //    FillCorrelations(). Also nMix should be passed in, so a weight\r
       //    of 1./nMix can be applied.\r
       \r
-      AliEventPool* pool = fPoolMgr->GetEventPool(fCentrality, eventMain->GetPrimaryVertex()->GetZ(),gReactionPlane);\r
+      AliEventPool* pool = fPoolMgr->GetEventPool(gCentrality, eventMain->GetPrimaryVertex()->GetZ(),gReactionPlane);\r
       \r
       if (!pool){\r
-       AliFatal(Form("No pool found for centrality = %f, zVtx = %f, psi = %f", fCentrality, eventMain->GetPrimaryVertex()->GetZ(),gReactionPlane));\r
+       AliFatal(Form("No pool found for centrality = %f, zVtx = %f, psi = %f", gCentrality, eventMain->GetPrimaryVertex()->GetZ(),gReactionPlane));\r
       }\r
       else{\r
        \r
@@ -584,10 +672,10 @@ Double_t AliAnalysisTaskBFPsi::IsEventAccepted(AliVEvent *event){
   // Fills Event statistics histograms\r
   \r
   Bool_t isSelectedMain = kTRUE;\r
-  Float_t fCentrality = -1.;\r
+  Float_t gCentrality = -1.;\r
   TString gAnalysisLevel = fBalance->GetAnalysisLevel();\r
 \r
-  fHistEventStats->Fill(1,fCentrality); //all events\r
+  fHistEventStats->Fill(1,gCentrality); //all events\r
 \r
   // Event trigger bits\r
   fHistTriggerStats->Fill(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected());\r
@@ -595,14 +683,14 @@ Double_t AliAnalysisTaskBFPsi::IsEventAccepted(AliVEvent *event){
     isSelectedMain = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();\r
   \r
   if(isSelectedMain) {\r
-    fHistEventStats->Fill(2,fCentrality); //triggered events\r
+    fHistEventStats->Fill(2,gCentrality); //triggered events\r
     \r
     //Centrality stuff \r
     if(fUseCentrality) {\r
       if(gAnalysisLevel == "AOD") { //centrality in AOD header\r
        AliAODHeader *header = (AliAODHeader*) event->GetHeader();\r
        if(header){\r
-         fCentrality = header->GetCentralityP()->GetCentralityPercentile(fCentralityEstimator.Data());\r
+         gCentrality = header->GetCentralityP()->GetCentralityPercentile(fCentralityEstimator.Data());\r
          \r
          // QA for centrality estimators\r
          fHistCentStats->Fill(0.,header->GetCentralityP()->GetCentralityPercentile("V0M"));\r
@@ -633,7 +721,7 @@ Double_t AliAnalysisTaskBFPsi::IsEventAccepted(AliVEvent *event){
 \r
       else if(gAnalysisLevel == "ESD" || gAnalysisLevel == "MCESD"){ // centrality class for ESDs or MC-ESDs\r
        AliCentrality *centrality = event->GetCentrality();\r
-       fCentrality = centrality->GetCentralityPercentile(fCentralityEstimator.Data());\r
+       gCentrality = centrality->GetCentralityPercentile(fCentralityEstimator.Data());\r
 \r
        // QA for centrality estimators\r
        fHistCentStats->Fill(0.,centrality->GetCentralityPercentile("V0M"));\r
@@ -654,11 +742,11 @@ Double_t AliAnalysisTaskBFPsi::IsEventAccepted(AliVEvent *event){
        AliGenHijingEventHeader* headerH = dynamic_cast<AliGenHijingEventHeader*>(dynamic_cast<AliMCEvent*>(event)->GenEventHeader());\r
        if(headerH){\r
          gImpactParameter = headerH->ImpactParameter();\r
-         fCentrality      = gImpactParameter;\r
+         gCentrality      = gImpactParameter;\r
        }//MC header\r
       }//MC\r
       else{\r
-       fCentrality = -1.;\r
+       gCentrality = -1.;\r
       }\r
     }\r
     \r
@@ -672,18 +760,18 @@ Double_t AliAnalysisTaskBFPsi::IsEventAccepted(AliVEvent *event){
        //gVertexArray.At(0),\r
        //gVertexArray.At(1),\r
        //gVertexArray.At(2));\r
-       fHistEventStats->Fill(3,fCentrality); //events with a proper vertex\r
+       fHistEventStats->Fill(3,gCentrality); //events with a proper vertex\r
        if(TMath::Abs(gVertexArray.At(0)) < fVxMax) {\r
          if(TMath::Abs(gVertexArray.At(1)) < fVyMax) {\r
            if(TMath::Abs(gVertexArray.At(2)) < fVzMax) {\r
-             fHistEventStats->Fill(4,fCentrality); //analayzed events\r
+             fHistEventStats->Fill(4,gCentrality); //analayzed events\r
              fHistVx->Fill(gVertexArray.At(0));\r
              fHistVy->Fill(gVertexArray.At(1));\r
-             fHistVz->Fill(gVertexArray.At(2),fCentrality);\r
+             fHistVz->Fill(gVertexArray.At(2),gCentrality);\r
              \r
              // take only events inside centrality class\r
-             if((fImpactParameterMin < fCentrality) && (fImpactParameterMax > fCentrality)){\r
-               return fCentrality;         \r
+             if((fImpactParameterMin < gCentrality) && (fImpactParameterMax > gCentrality)){\r
+               return gCentrality;         \r
              }//centrality class\r
            }//Vz cut\r
          }//Vy cut\r
@@ -700,19 +788,19 @@ Double_t AliAnalysisTaskBFPsi::IsEventAccepted(AliVEvent *event){
        vertex->GetCovarianceMatrix(fCov);\r
        if(vertex->GetNContributors() > 0) {\r
          if(fCov[5] != 0) {\r
-           fHistEventStats->Fill(3,fCentrality); //events with a proper vertex\r
+           fHistEventStats->Fill(3,gCentrality); //events with a proper vertex\r
            if(TMath::Abs(vertex->GetX()) < fVxMax) {\r
              if(TMath::Abs(vertex->GetY()) < fVyMax) {\r
                if(TMath::Abs(vertex->GetZ()) < fVzMax) {\r
-                 fHistEventStats->Fill(4,fCentrality); //analyzed events\r
+                 fHistEventStats->Fill(4,gCentrality); //analyzed events\r
                  fHistVx->Fill(vertex->GetX());\r
                  fHistVy->Fill(vertex->GetY());\r
-                 fHistVz->Fill(vertex->GetZ(),fCentrality);\r
+                 fHistVz->Fill(vertex->GetZ(),gCentrality);\r
                  \r
                  // take only events inside centrality class\r
-                 if((fCentrality > fCentralityPercentileMin) && (fCentrality < fCentralityPercentileMax)){\r
-                   fHistEventStats->Fill(5,fCentrality); //events with correct centrality\r
-                   return fCentrality;         \r
+                 if((gCentrality > fCentralityPercentileMin) && (gCentrality < fCentralityPercentileMax)){\r
+                   fHistEventStats->Fill(5,gCentrality); //events with correct centrality\r
+                   return gCentrality;         \r
                  }//centrality class\r
                }//Vz cut\r
              }//Vy cut\r
@@ -733,7 +821,7 @@ Double_t AliAnalysisTaskBFPsi::GetRefMultiOrCentrality(AliVEvent *event){
     // Checks the Event cuts\r
     // Fills Event statistics histograms\r
   \r
-  Float_t fCentrality = -1.;\r
+  Float_t gCentrality = -1.;\r
   Double_t fMultiplicity = -100.;\r
   TString gAnalysisLevel = fBalance->GetAnalysisLevel();\r
   if(fEventClass == "Centrality"){\r
@@ -742,24 +830,24 @@ Double_t AliAnalysisTaskBFPsi::GetRefMultiOrCentrality(AliVEvent *event){
     if(gAnalysisLevel == "AOD") { //centrality in AOD header\r
       AliAODHeader *header = (AliAODHeader*) event->GetHeader();\r
       if(header){\r
-        fCentrality = header->GetCentralityP()->GetCentralityPercentile(fCentralityEstimator.Data());\r
+        gCentrality = header->GetCentralityP()->GetCentralityPercentile(fCentralityEstimator.Data());\r
       }//AOD header\r
     }//AOD\r
     \r
     else if(gAnalysisLevel == "ESD" || gAnalysisLevel == "MCESD"){ // centrality class for ESDs or MC-ESDs\r
       AliCentrality *centrality = event->GetCentrality();\r
-      fCentrality = centrality->GetCentralityPercentile(fCentralityEstimator.Data());\r
+      gCentrality = centrality->GetCentralityPercentile(fCentralityEstimator.Data());\r
     }//ESD\r
     else if(gAnalysisLevel == "MC"){\r
       Double_t gImpactParameter = 0.;\r
       AliGenHijingEventHeader* headerH = dynamic_cast<AliGenHijingEventHeader*>(dynamic_cast<AliMCEvent*>(event)->GenEventHeader());\r
       if(headerH){\r
         gImpactParameter = headerH->ImpactParameter();\r
-        fCentrality      = gImpactParameter;\r
+        gCentrality      = gImpactParameter;\r
       }//MC header\r
     }//MC\r
     else{\r
-      fCentrality = -1.;\r
+      gCentrality = -1.;\r
     }\r
   }//End if "Centrality"\r
   if(fEventClass=="Multiplicity"&&gAnalysisLevel == "ESD"){\r
@@ -775,7 +863,7 @@ Double_t AliAnalysisTaskBFPsi::GetRefMultiOrCentrality(AliVEvent *event){
   if(fEventClass=="Multiplicity"){\r
     lReturnVal = fMultiplicity;\r
   }else if(fEventClass=="Centrality"){\r
-    lReturnVal = fCentrality;\r
+    lReturnVal = gCentrality;\r
   }\r
   return lReturnVal;\r
 }\r
@@ -815,8 +903,46 @@ Double_t AliAnalysisTaskBFPsi::GetEventPlane(AliVEvent *event){
   return gReactionPlane;\r
 }\r
 \r
+//========================correction=============================//\r
+Double_t AliAnalysisTaskBFPsi::GetTrackbyTrackCorrectionMatrix( Double_t vEta, Double_t vPhi, \r
+                                                               Double_t vPt, Short_t vCharge, Double_t gCentrality) {\r
+  // -- Get efficiency correction of particle dependent on (eta, phi, pt, charge, centrality)  \r
+  \r
+  Double_t correction = 1.;\r
+  //Double_t dimBin[3] = {vEta, vPhi, vPt, gCentrality}; // eta, phi, pt, centrality\r
+ \r
+  Double_t widthEta = (fEtaMaxForCorrections - fEtaMinForCorrections)/fEtaBinForCorrections;\r
+  Int_t binEta = (Int_t)(vEta/widthEta)+1;\r
+  Double_t widthPt = (fPtMaxForCorrections - fPtMinForCorrections)/fPtBinForCorrections;\r
+  Int_t binPt = (Int_t)(vPt/widthPt) + 1;\r
+  Double_t widthPhi = (fPhiMaxForCorrections - fPhiMinForCorrections)/fPhiBinForCorrections;\r
+  Int_t binPhi = (Int_t)(vPhi/widthPhi)+ 1;\r
+\r
+  Int_t gCentralityInt = 1;\r
+  for (Int_t i=1; i<=kCENTRALITY; i++){\r
+    if ((Int_t)(centralityArrayForPbPb[i]) <= gCentrality <= (Int_t)(centralityArrayForPbPb[i])){\r
+      gCentralityInt = i;\r
+    }\r
+  }\r
+    \r
+  if (vCharge > 0) {\r
+    correction = fHistMatrixCorrectionPlus[gCentralityInt-1]->GetBinContent(fHistMatrixCorrectionPlus[gCentralityInt-1]->GetBin(binEta, binPt, binPhi));\r
+  }\r
+  if (vCharge < 0) {\r
+    //correction = fHistMatrixCorrectionMinus->GetBinContent(fHistMatrixCorrectionMinus->GetBin(dimBin));\r
+    correction = fHistMatrixCorrectionMinus[gCentralityInt-1]->GetBinContent(fHistMatrixCorrectionMinus[gCentralityInt-1]->GetBin(binEta, binPt, binPhi));\r
+  }\r
+  if (correction == 0.) { \r
+    AliError(Form("Should not happen : bin content = 0. >> eta: %.2f | phi : %.2f | pt : %.2f | cent %d",vEta, vPhi, vPt, gCentralityInt)); \r
+    correction = 1.; \r
+  } \r
+\r
+  return correction;\r
+}\r
+//========================correction=============================//\r
+\r
 //________________________________________________________________________\r
-TObjArray* AliAnalysisTaskBFPsi::GetAcceptedTracks(AliVEvent *event, Double_t fCentrality, Double_t gReactionPlane){\r
+TObjArray* AliAnalysisTaskBFPsi::GetAcceptedTracks(AliVEvent *event, Double_t gCentrality, Double_t gReactionPlane){\r
   // Returns TObjArray with tracks after all track cuts (only for AOD!)\r
   // Fills QA histograms\r
 \r
@@ -826,7 +952,7 @@ TObjArray* AliAnalysisTaskBFPsi::GetAcceptedTracks(AliVEvent *event, Double_t fC
   TObjArray* tracksAccepted = new TObjArray;\r
   tracksAccepted->SetOwner(kTRUE);\r
 \r
-  Double_t vCharge;\r
+  Short_t vCharge;\r
   Double_t vEta;\r
   Double_t vY;\r
   Double_t vPhi;\r
@@ -884,16 +1010,19 @@ TObjArray* AliAnalysisTaskBFPsi::GetAcceptedTracks(AliVEvent *event, Double_t fC
       // fill QA histograms\r
       fHistClus->Fill(aodTrack->GetITSNcls(),aodTrack->GetTPCNcls());\r
       fHistDCA->Fill(dcaZ,dcaXY);\r
-      fHistChi2->Fill(aodTrack->Chi2perNDF(),fCentrality);\r
-      fHistPt->Fill(vPt,fCentrality);\r
-      fHistEta->Fill(vEta,fCentrality);\r
-      fHistRapidity->Fill(vY,fCentrality);\r
-      if(vCharge > 0) fHistPhiPos->Fill(vPhi,fCentrality);\r
-      else if(vCharge < 0) fHistPhiNeg->Fill(vPhi,fCentrality);\r
-      fHistPhi->Fill(vPhi,fCentrality);\r
+      fHistChi2->Fill(aodTrack->Chi2perNDF(),gCentrality);\r
+      fHistPt->Fill(vPt,gCentrality);\r
+      fHistEta->Fill(vEta,gCentrality);\r
+      fHistRapidity->Fill(vY,gCentrality);\r
+      if(vCharge > 0) fHistPhiPos->Fill(vPhi,gCentrality);\r
+      else if(vCharge < 0) fHistPhiNeg->Fill(vPhi,gCentrality);\r
+      fHistPhi->Fill(vPhi,gCentrality);\r
+      \r
+      //=======================================correction\r
+      Double_t correction = GetTrackbyTrackCorrectionMatrix(vEta, vPhi, vPt, vCharge, gCentrality);  \r
       \r
       // add the track to the TObjArray\r
-      tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, 1.*vCharge));\r
+      tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge, correction));  \r
     }//track loop\r
   }// AOD analysis\r
 \r
@@ -1064,16 +1193,19 @@ TObjArray* AliAnalysisTaskBFPsi::GetAcceptedTracks(AliVEvent *event, Double_t fC
       vPt     = trackTPC->Pt();\r
       fHistClus->Fill(trackTPC->GetITSclusters(0),nClustersTPC);\r
       fHistDCA->Fill(b[1],b[0]);\r
-      fHistChi2->Fill(chi2PerClusterTPC,fCentrality);\r
-      fHistPt->Fill(vPt,fCentrality);\r
-      fHistEta->Fill(vEta,fCentrality);\r
-      fHistPhi->Fill(vPhi,fCentrality);\r
-      fHistRapidity->Fill(vY,fCentrality);\r
-      if(vCharge > 0) fHistPhiPos->Fill(vPhi,fCentrality);\r
-      else if(vCharge < 0) fHistPhiNeg->Fill(vPhi,fCentrality);\r
+      fHistChi2->Fill(chi2PerClusterTPC,gCentrality);\r
+      fHistPt->Fill(vPt,gCentrality);\r
+      fHistEta->Fill(vEta,gCentrality);\r
+      fHistPhi->Fill(vPhi,gCentrality);\r
+      fHistRapidity->Fill(vY,gCentrality);\r
+      if(vCharge > 0) fHistPhiPos->Fill(vPhi,gCentrality);\r
+      else if(vCharge < 0) fHistPhiNeg->Fill(vPhi,gCentrality);\r
       \r
+      //=======================================correction\r
+      Double_t correction = GetTrackbyTrackCorrectionMatrix(vEta, vPhi, vPt, vCharge, gCentrality);  \r
+\r
       // add the track to the TObjArray\r
-      tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge));      \r
+      tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge, correction));   \r
 \r
       delete trackTPC;\r
     }//track loop\r
@@ -1151,15 +1283,15 @@ TObjArray* AliAnalysisTaskBFPsi::GetAcceptedTracks(AliVEvent *event, Double_t fC
       vPhi    = track->Phi();\r
       //Printf("phi (before): %lf",vPhi);\r
       \r
-      fHistPt->Fill(vPt,fCentrality);\r
-      fHistEta->Fill(vEta,fCentrality);\r
-      fHistPhi->Fill(vPhi,fCentrality);\r
-      //fHistPhi->Fill(vPhi*TMath::RadToDeg(),fCentrality);\r
-      fHistRapidity->Fill(vY,fCentrality);\r
-      //if(vCharge > 0) fHistPhiPos->Fill(vPhi*TMath::RadToDeg(),fCentrality);\r
-      //else if(vCharge < 0) fHistPhiNeg->Fill(vPhi*TMath::RadToDeg(),fCentrality);\r
-      if(vCharge > 0) fHistPhiPos->Fill(vPhi,fCentrality);\r
-      else if(vCharge < 0) fHistPhiNeg->Fill(vPhi,fCentrality);\r
+      fHistPt->Fill(vPt,gCentrality);\r
+      fHistEta->Fill(vEta,gCentrality);\r
+      fHistPhi->Fill(vPhi,gCentrality);\r
+      //fHistPhi->Fill(vPhi*TMath::RadToDeg(),gCentrality);\r
+      fHistRapidity->Fill(vY,gCentrality);\r
+      //if(vCharge > 0) fHistPhiPos->Fill(vPhi*TMath::RadToDeg(),gCentrality);\r
+      //else if(vCharge < 0) fHistPhiNeg->Fill(vPhi*TMath::RadToDeg(),gCentrality);\r
+      if(vCharge > 0) fHistPhiPos->Fill(vPhi,gCentrality);\r
+      else if(vCharge < 0) fHistPhiNeg->Fill(vPhi,gCentrality);\r
       \r
       //Flow after burner\r
       if(fUseFlowAfterBurner) {\r
@@ -1179,16 +1311,20 @@ TObjArray* AliAnalysisTaskBFPsi::GetAcceptedTracks(AliVEvent *event, Double_t fC
        //Printf("phi (after): %lf\n",vPhi);\r
        Double_t vDeltaphiBefore = phi0 - gReactionPlane*TMath::DegToRad();\r
        if(vDeltaphiBefore < 0) vDeltaphiBefore += 2*TMath::Pi();\r
-       fHistPhiBefore->Fill(vDeltaphiBefore,fCentrality);\r
+       fHistPhiBefore->Fill(vDeltaphiBefore,gCentrality);\r
        \r
        Double_t vDeltaphiAfter = vPhi - gReactionPlane*TMath::DegToRad();\r
        if(vDeltaphiAfter < 0) vDeltaphiAfter += 2*TMath::Pi();\r
-       fHistPhiAfter->Fill(vDeltaphiAfter,fCentrality);\r
+       fHistPhiAfter->Fill(vDeltaphiAfter,gCentrality);\r
+      \r
       }\r
       \r
       //vPhi *= TMath::RadToDeg();\r
-                \r
-      tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge));\r
+\r
+      //=======================================correction\r
+      Double_t correction = GetTrackbyTrackCorrectionMatrix(vEta, vPhi, vPt, vCharge, gCentrality);  \r
+      \r
+      tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge, correction)); \r
       \r
     } //track loop\r
   }//MC\r
@@ -1196,7 +1332,7 @@ TObjArray* AliAnalysisTaskBFPsi::GetAcceptedTracks(AliVEvent *event, Double_t fC
   return tracksAccepted;  \r
 }\r
 //________________________________________________________________________\r
-TObjArray* AliAnalysisTaskBFPsi::GetShuffledTracks(TObjArray *tracks){\r
+TObjArray* AliAnalysisTaskBFPsi::GetShuffledTracks(TObjArray *tracks, Double_t gCentrality){\r
   // Clones TObjArray and returns it with tracks after shuffling the charges\r
 \r
   TObjArray* tracksShuffled = new TObjArray;\r
@@ -1214,7 +1350,9 @@ TObjArray* AliAnalysisTaskBFPsi::GetShuffledTracks(TObjArray *tracks){
   \r
   for(Int_t i = 0; i < tracks->GetEntriesFast(); i++){\r
     AliVParticle* track = (AliVParticle*) tracks->At(i);\r
-    tracksShuffled->Add(new AliBFBasicParticle(track->Eta(), track->Phi(), track->Pt(),chargeVector->at(i)));\r
+    //==============================correction\r
+    Double_t correction = GetTrackbyTrackCorrectionMatrix(track->Eta(), track->Phi(),track->Pt(), chargeVector->at(i), gCentrality);\r
+    tracksShuffled->Add(new AliBFBasicParticle(track->Eta(), track->Phi(), track->Pt(),chargeVector->at(i), correction));\r
   }\r
 \r
   delete chargeVector;\r

diff --git a/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskBFPsi.h b/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskBFPsi.h
index 346ba27d305fbeb662b7e48c4384739fed759185..24f72e99ba736ffa7b91d743f8ee28f4a074d59b 100755 (executable)
--- a/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskBFPsi.h
+++ b/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskBFPsi.h
@@ -8,6 +8,7 @@ class TList;
 class TH1F;\r
 class TH2F;\r
 class TF1;\r
+class TH3D;\r
 \r
 class AliBalancePsi;\r
 class AliESDtrackCuts;\r
@@ -21,18 +22,25 @@ class AliEventPoolManager;
 #include "AliPIDResponse.h"\r
 #include "AliPIDCombined.h"\r
  \r
+//================================correction\r
+#define kCENTRALITY 9  \r
+const Double_t centralityArrayForPbPb[kCENTRALITY+1] = {0.,5.,10.,20.,30.,40.,50.,60.,70.,80.};\r
+//================================correction\r
 \r
 class AliAnalysisTaskBFPsi : public AliAnalysisTaskSE {\r
  public:\r
   AliAnalysisTaskBFPsi(const char *name = "AliAnalysisTaskBFPsi");\r
   virtual ~AliAnalysisTaskBFPsi(); \r
-  \r
-  \r
+   \r
   virtual void   UserCreateOutputObjects();\r
   virtual void   UserExec(Option_t *option);\r
   virtual void   FinishTaskOutput();\r
   virtual void   Terminate(Option_t *);\r
 \r
+  //========================correction\r
+  virtual void   SetInputCorrection(const char* filename, const char* gSystem); \r
+  //========================correction\r
+\r
   void SetAnalysisObject(AliBalancePsi *const analysis) {\r
     fBalance         = analysis;\r
     }\r
@@ -137,13 +145,19 @@ class AliAnalysisTaskBFPsi : public AliAnalysisTaskSE {
         fEventClass = receivedEventClass;\r
     }\r
 \r
+    //========================correction=============================//\r
+    // void SetInputCorrection(const char* filename = "efficiencyCorrection.root"){}\r
+\r
  private:\r
   Double_t    IsEventAccepted(AliVEvent* event);\r
   Double_t    GetRefMultiOrCentrality(AliVEvent* event);\r
   Double_t    GetEventPlane(AliVEvent* event);\r
-  TObjArray* GetAcceptedTracks(AliVEvent* event, Double_t fCentrality, Double_t gReactionPlane);\r
-  TObjArray* GetShuffledTracks(TObjArray* tracks);\r
-\r
+  //===============================correction\r
+  Double_t    GetTrackbyTrackCorrectionMatrix(Double_t vEta, Double_t vPhi, Double_t vPt, Short_t vCharge, Double_t gCentrality);\r
+  //===============================correction\r
+  TObjArray* GetAcceptedTracks(AliVEvent* event, Double_t gCentrality, Double_t gReactionPlane);\r
+  TObjArray* GetShuffledTracks(TObjArray* tracks, Double_t gCentrality);\r
+  \r
   AliBalancePsi *fBalance; //BF object\r
   Bool_t fRunShuffling;//run shuffling or not\r
   AliBalancePsi *fShuffledBalance; //BF object (shuffled)\r
@@ -198,6 +212,9 @@ class AliAnalysisTaskBFPsi : public AliAnalysisTaskSE {
   TH2D *fHistProbTPCTOFvsPtafterPID;//TOF/TPC probability vs pT after PID cuts (QA histogram)\r
   TH2D *fHistNSigmaTPCvsPtafterPID;//TPC nsigma vs pT after PID cuts (QA histogram)\r
   TH2D *fHistNSigmaTOFvsPtafterPID;//TOF nsigma vs pT after PID cuts (QA histogram)\r
+  \r
+  TH3D *fHistMatrixCorrectionPlus[kCENTRALITY]; //======================================================correction\r
+  TH3D *fHistMatrixCorrectionMinus[kCENTRALITY]; //=====================================================correction\r
 \r
   AliPIDResponse *fPIDResponse;     //! PID response object\r
   AliPIDCombined       *fPIDCombined;     //! combined PID object\r
@@ -236,8 +253,19 @@ class AliAnalysisTaskBFPsi : public AliAnalysisTaskSE {
 \r
   Double_t fPtMin;//only used for AODs\r
   Double_t fPtMax;//only used for AODs\r
+  Double_t fPtMinForCorrections;//only used for AODs\r
+  Double_t fPtMaxForCorrections;//only used for AODs\r
+  Double_t fPtBinForCorrections; //=================================correction\r
   Double_t fEtaMin;//only used for AODs\r
   Double_t fEtaMax;//only used for AODs\r
+  Double_t fEtaMinForCorrections;//only used for AODs\r
+  Double_t fEtaMaxForCorrections;//only used for AODs\r
+  Double_t fEtaBinForCorrections; //=================================correction\r
+  Double_t fPhiMin; //=================================correction \r
+  Double_t fPhiMax; //=================================correction\r
+  Double_t fPhiMinForCorrections;//only used for AODs\r
+  Double_t fPhiMaxForCorrections;//only used for AODs\r
+  Double_t fPhiBinForCorrections; //=================================correction\r
 \r
   Double_t fDCAxyCut;//only used for AODs\r
   Double_t fDCAzCut;//only used for AODs\r

diff --git a/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskTriggeredBF.cxx b/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskTriggeredBF.cxx
index 7b13b4afd5120351529d32622b59033d5d676156..ae49008dd30e759d92b999137935ebbbd45ea196 100755 (executable)
--- a/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskTriggeredBF.cxx
+++ b/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskTriggeredBF.cxx
@@ -707,7 +707,7 @@ TObjArray* AliAnalysisTaskTriggeredBF::GetAcceptedTracks(AliVEvent *event){
   TObjArray* tracksAccepted = new TObjArray;\r
   tracksAccepted->SetOwner(kTRUE);\r
 \r
-  Double_t vCharge = 0.;\r
+  Short_t vCharge = 0;\r
   Double_t vEta    = 0.;\r
   Double_t vPhi    = 0.;\r
   Double_t vPt     = 0.;\r
@@ -764,7 +764,7 @@ TObjArray* AliAnalysisTaskTriggeredBF::GetAcceptedTracks(AliVEvent *event){
     fHistPhi->Fill(vPhi);\r
     \r
     // add the track to the TObjArray\r
-    tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge));\r
+    tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge,1.));\r
   }\r
 \r
   return tracksAccepted;\r
@@ -779,7 +779,7 @@ TObjArray* AliAnalysisTaskTriggeredBF::GetAcceptedV0s(AliVEvent *event){
   TObjArray* tracksAccepted = new TObjArray;\r
   tracksAccepted->SetOwner(kTRUE);\r
 \r
-  Double_t vCharge = 0.;\r
+  Short_t vCharge = 0;\r
   Double_t vEta    = 0.;\r
   Double_t vPhi    = 0.;\r
   Double_t vPt     = 0.;\r
@@ -1040,7 +1040,7 @@ TObjArray* AliAnalysisTaskTriggeredBF::GetAcceptedV0s(AliVEvent *event){
                    fHistPhi->Fill(vPhi);\r
                    \r
                    // add the track to the TObjArray\r
-                   tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge));\r
+                   tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge,1.));\r
                  }\r
                }\r
              }\r
@@ -1072,7 +1072,7 @@ TObjArray* AliAnalysisTaskTriggeredBF::GetShuffledTracks(TObjArray *tracks){
   \r
   for(Int_t i = 0; i < tracks->GetEntriesFast(); i++){\r
     AliVParticle* track = (AliVParticle*) tracks->At(i);\r
-    tracksShuffled->Add(new AliBFBasicParticle(track->Eta(), track->Phi(), track->Pt(),chargeVector->at(i)));\r
+    tracksShuffled->Add(new AliBFBasicParticle(track->Eta(), track->Phi(), track->Pt(),chargeVector->at(i),1.));\r
   }\r
 \r
   delete chargeVector;\r

diff --git a/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskTriggeredBF.h b/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskTriggeredBF.h
index 42743225c15389cc33f17160816a5be5f030b3ba..1cb244a3ff70421e803190e8d57dd29802898928 100755 (executable)
--- a/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskTriggeredBF.h
+++ b/PWGCF/EBYE/BalanceFunctions/AliAnalysisTaskTriggeredBF.h
@@ -206,8 +206,8 @@ class AliAnalysisTaskTriggeredBF : public AliAnalysisTaskSE {
 class AliBFBasicParticle : public AliVParticle\r
 { \r
   public:\r
-    AliBFBasicParticle(Float_t eta, Float_t phi, Float_t pt, Short_t charge)\r
-      : fEta(eta), fPhi(phi), fpT(pt), fCharge(charge)\r
+  AliBFBasicParticle(Float_t eta, Float_t phi, Float_t pt, Short_t charge, Double_t correction)\r
+    : fEta(eta), fPhi(phi), fpT(pt), fCharge(charge), fCorrection(correction)\r
     {\r
     }\r
     ~AliBFBasicParticle() {}\r
@@ -238,6 +238,9 @@ class AliBFBasicParticle : public AliVParticle
     \r
     virtual Short_t Charge()      const { return fCharge; }\r
     virtual Int_t   GetLabel()    const { AliFatal("Not implemented"); return 0; }\r
+\r
+    virtual Double_t Correction()        const { return fCorrection; } //=============================correction\r
+    \r
     // PID\r
     virtual Int_t   PdgCode()     const { AliFatal("Not implemented"); return 0; }      \r
     virtual const Double_t *PID() const { AliFatal("Not implemented"); return 0; }\r
@@ -247,7 +250,8 @@ class AliBFBasicParticle : public AliVParticle
     Float_t fPhi;      // phi\r
     Float_t fpT;       // pT\r
     Short_t fCharge;   // charge\r
-    \r
+    Double_t fCorrection; //============================correction\r
+\r
     ClassDef( AliBFBasicParticle, 1); // class which contains only quantities requires for this analysis to reduce memory consumption for event mixing\r
 };\r
 \r

diff --git a/PWGCF/EBYE/BalanceFunctions/AliBalancePsi.cxx b/PWGCF/EBYE/BalanceFunctions/AliBalancePsi.cxx
index a307c6edab600422b94f39c87d14a4736ce6902b..7015d4f600441834534e616e813bb2062e3df6f1 100644 (file)
--- a/PWGCF/EBYE/BalanceFunctions/AliBalancePsi.cxx
+++ b/PWGCF/EBYE/BalanceFunctions/AliBalancePsi.cxx
@@ -37,6 +37,7 @@
 #include "AliESDtrack.h"
 #include "AliAODTrack.h"
 #include "AliTHn.h"
+#include "AliAnalysisTaskTriggeredBF.h"
 
 #include "AliBalancePsi.h"
 
@@ -340,23 +341,27 @@ void AliBalancePsi::CalculateBalance(Double_t gReactionPlane,
   TArrayF secondPhi(jMax);
   TArrayF secondPt(jMax);
   TArrayS secondCharge(jMax);
+  TArrayD secondCorrection(jMax);
 
   for (Int_t i=0; i<jMax; i++){
     secondEta[i] = ((AliVParticle*) particlesSecond->At(i))->Eta();
     secondPhi[i] = ((AliVParticle*) particlesSecond->At(i))->Phi();
     secondPt[i]  = ((AliVParticle*) particlesSecond->At(i))->Pt();
     secondCharge[i]  = (Short_t)((AliVParticle*) particlesSecond->At(i))->Charge();
+    secondCorrection[i]  = (Double_t)((AliBFBasicParticle*) particlesSecond->At(i))->Correction();   //==========================correction
   }
   
   // 1st particle loop
   for (Int_t i = 0; i < iMax; i++) {
-    AliVParticle* firstParticle = (AliVParticle*) particles->At(i);
+    //AliVParticle* firstParticle = (AliVParticle*) particles->At(i);
+    AliBFBasicParticle* firstParticle = (AliBFBasicParticle*) particles->At(i); //==========================correction
     
     // some optimization
     Float_t firstEta = firstParticle->Eta();
     Float_t firstPhi = firstParticle->Phi();
     Float_t firstPt  = firstParticle->Pt();
-    
+    Float_t firstCorrection  = firstParticle->Correction();//==========================correction
+
     // Event plane (determine psi bin)
     Double_t gPsiMinusPhi    =   0.;
     Double_t gPsiMinusPhiBin = -10.;
@@ -388,8 +393,8 @@ void AliBalancePsi::CalculateBalance(Double_t gReactionPlane,
       if(fEventClass=="Multiplicity" || fEventClass == "Centrality" ) trackVariablesSingle[0] = kMultorCent;
     
     //fill single particle histograms
-    if(charge1 > 0)      fHistP->Fill(trackVariablesSingle,0,1.); 
-    else if(charge1 < 0) fHistN->Fill(trackVariablesSingle,0,1.);  
+    if(charge1 > 0)      fHistP->Fill(trackVariablesSingle,0,firstCorrection); //==========================correction
+    else if(charge1 < 0) fHistN->Fill(trackVariablesSingle,0,firstCorrection);  //==========================correction
     
     // 2nd particle loop
     for(Int_t j = 0; j < jMax; j++) {   
@@ -502,10 +507,10 @@ void AliBalancePsi::CalculateBalance(Double_t gReactionPlane,
        }
       }//conversion cut
       
-      if( charge1 > 0 && charge2 < 0)  fHistPN->Fill(trackVariablesPair,0,1.); 
-      else if( charge1 < 0 && charge2 > 0)  fHistNP->Fill(trackVariablesPair,0,1.); 
-      else if( charge1 > 0 && charge2 > 0)  fHistPP->Fill(trackVariablesPair,0,1.); 
-      else if( charge1 < 0 && charge2 < 0)  fHistNN->Fill(trackVariablesPair,0,1.); 
+      if( charge1 > 0 && charge2 < 0)  fHistPN->Fill(trackVariablesPair,0,firstCorrection*secondCorrection[j]); //==========================correction
+      else if( charge1 < 0 && charge2 > 0)  fHistNP->Fill(trackVariablesPair,0,firstCorrection*secondCorrection[j]);//==========================correction 
+      else if( charge1 > 0 && charge2 > 0)  fHistPP->Fill(trackVariablesPair,0,firstCorrection*secondCorrection[j]);//==========================correction 
+      else if( charge1 < 0 && charge2 < 0)  fHistNN->Fill(trackVariablesPair,0,firstCorrection*secondCorrection[j]);//==========================correction 
       else {
        //AliWarning(Form("Wrong charge combination: charge1 = %d and charge2 = %d",charge,charge2));
        continue;