]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Add histo to study fakes versus Rabs. Add histo and modify counters to study fakes...
authormartinez <martinez@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 14 Jul 2011 17:25:27 +0000 (17:25 +0000)
committermartinez <martinez@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 14 Jul 2011 17:25:27 +0000 (17:25 +0000)
PWG3/muondep/AddTaskMuonFakes.C
PWG3/muondep/AliAnalysisTaskMuonFakes.cxx
PWG3/muondep/AliAnalysisTaskMuonFakes.h

index 12d07d023d86d6343395962951a4195cf221e4e0..bd075119f528719037d228c5f8b0a645fcb8310e 100644 (file)
@@ -1,4 +1,4 @@
-AliAnalysisTaskMuonFakes* AddTaskMuonFakes(Bool_t useMCLabels = kFALSE)
+AliAnalysisTaskMuonFakes* AddTaskMuonFakes(Bool_t useMCLabels = kFALSE, Bool_t matchTrig = kFALSE, Bool_t applyAccCut = kFALSE)
 {
   /// Add AliAnalysisTaskMuonFakes to the train (Philippe Pillot)
   
@@ -23,6 +23,8 @@ AliAnalysisTaskMuonFakes* AddTaskMuonFakes(Bool_t useMCLabels = kFALSE)
     return NULL;
   }
   task->UseMCLabels(useMCLabels);
+  task->MatchTrigger(matchTrig);
+  task->ApplyAccCut(applyAccCut);
   
   // Add task to analysis manager
   mgr->AddTask(task);
@@ -44,11 +46,15 @@ AliAnalysisTaskMuonFakes* AddTaskMuonFakes(Bool_t useMCLabels = kFALSE)
   AliAnalysisDataContainer *cout_fakeTrack = mgr->CreateContainer("fakeTrackCounters", AliCounterCollection::Class(), AliAnalysisManager::kOutputContainer, outputfile.Data());
   AliAnalysisDataContainer *cout_matchTrack = mgr->CreateContainer("matchedTrackCounters", AliCounterCollection::Class(), AliAnalysisManager::kOutputContainer, outputfile.Data());
   AliAnalysisDataContainer *cout_event = mgr->CreateContainer("eventCounters", AliCounterCollection::Class(), AliAnalysisManager::kOutputContainer, outputfile.Data());
+  AliAnalysisDataContainer *cout_histo2 = mgr->CreateContainer("histos2", TObjArray::Class(), AliAnalysisManager::kOutputContainer, outputfile.Data());
+  AliAnalysisDataContainer *cout_pair = mgr->CreateContainer("pairCounters", AliCounterCollection::Class(), AliAnalysisManager::kOutputContainer, outputfile.Data());
   mgr->ConnectOutput(task, 1, cout_histo);
   mgr->ConnectOutput(task, 2, cout_track);
   mgr->ConnectOutput(task, 3, cout_fakeTrack);
   mgr->ConnectOutput(task, 4, cout_matchTrack);
   mgr->ConnectOutput(task, 5, cout_event);
+  mgr->ConnectOutput(task, 6, cout_histo2);
+  mgr->ConnectOutput(task, 7, cout_pair);
   
   return task;
 }
index 4fa4439ab759f854ea352fce9081f61031194080..12c441d7474cc61debe84788f655aee67804e15d 100644 (file)
 #include "AliESDMuonTrack.h"
 #include "AliESDInputHandler.h"
 #include "AliMCEventHandler.h"
+#include "AliCDBManager.h"
+#include "AliCentrality.h"
 
 // ANALYSIS includes
 #include "AliAnalysisDataSlot.h"
 #include "AliAnalysisDataContainer.h"
 #include "AliAnalysisManager.h"
-#include "AliCDBManager.h"
 
 // MUON includes
 #include "AliMUONCDB.h"
@@ -55,16 +56,20 @@ ClassImp(AliAnalysisTaskMuonFakes)
 AliAnalysisTaskMuonFakes::AliAnalysisTaskMuonFakes() :
 AliAnalysisTaskSE(), 
 fList(0x0),
+fList2(0x0),
 fCanvases(0x0),
 fTrackCounters(0x0),
 fFakeTrackCounters(0x0),
 fMatchedTrackCounters(0x0),
 fEventCounters(0x0),
+fPairCounters(0x0),
 fCurrentFileName(""),
 fRequestedStationMask(0),
 fRequest2ChInSameSt45(kFALSE),
 fSigmaCut(-1.),
 fUseLabel(kFALSE),
+fMatchTrig(kFALSE),
+fApplyAccCut(kFALSE),
 fRecoParamLocation("alien://folder=/alice/simulation/2008/v4-15-Release/Full")
 {
   /// Default constructor.
@@ -74,16 +79,20 @@ fRecoParamLocation("alien://folder=/alice/simulation/2008/v4-15-Release/Full")
 AliAnalysisTaskMuonFakes::AliAnalysisTaskMuonFakes(const char *name) :
 AliAnalysisTaskSE(name), 
 fList(0x0),
+fList2(0x0),
 fCanvases(0x0),
 fTrackCounters(0x0),
 fFakeTrackCounters(0x0),
 fMatchedTrackCounters(0x0),
 fEventCounters(0x0),
+fPairCounters(0x0),
 fCurrentFileName(""),
 fRequestedStationMask(0),
 fRequest2ChInSameSt45(kFALSE),
 fSigmaCut(-1.),
 fUseLabel(kFALSE),
+fMatchTrig(kFALSE),
+fApplyAccCut(kFALSE),
 fRecoParamLocation("alien://folder=/alice/simulation/2008/v4-15-Release/Full")
 {
   /// Constructor.
@@ -97,6 +106,10 @@ fRecoParamLocation("alien://folder=/alice/simulation/2008/v4-15-Release/Full")
   DefineOutput(4,AliCounterCollection::Class());
   // Output slot #5 writes into an AliCounterCollection container
   DefineOutput(5,AliCounterCollection::Class());
+  // Output slot #6 writes into a TObjArray container
+  DefineOutput(6,TObjArray::Class());
+  // Output slot #7 writes into an AliCounterCollection container
+  DefineOutput(7,AliCounterCollection::Class());
 }
 
 //________________________________________________________________________
@@ -105,10 +118,12 @@ AliAnalysisTaskMuonFakes::~AliAnalysisTaskMuonFakes()
   /// Destructor.
   if (!AliAnalysisManager::GetAnalysisManager()->IsProofMode()) {
     delete fList;
+    delete fList2;
     delete fTrackCounters;
     delete fFakeTrackCounters;
     delete fMatchedTrackCounters;
     delete fEventCounters;
+    delete fPairCounters;
   }
   delete fCanvases;
 }
@@ -118,6 +133,7 @@ void AliAnalysisTaskMuonFakes::UserCreateOutputObjects()
 {
   /// Create histograms and counters.
   
+  // single track histograms
   fList = new TObjArray(100);
   fList->SetOwner();
   
@@ -150,11 +166,11 @@ void AliAnalysisTaskMuonFakes::UserCreateOutputObjects()
   fList->AddAtAndExpand(hNumberOfChamberHitF, kNumberOfChamberHitF);
   
   // chi2
-  TH1F *hChi2PerDof = new TH1F("hChi2PerDof", "track chi2/d.o.f.", 100, 0., 20.);
+  TH1F *hChi2PerDof = new TH1F("hChi2PerDof", "track chi2/d.o.f.", 200, 0., 20.);
   fList->AddAtAndExpand(hChi2PerDof, kChi2PerDof);
-  TH1F *hChi2PerDofM = new TH1F("hChi2PerDofM", "matched track chi2/d.o.f.", 100, 0., 20.);
+  TH1F *hChi2PerDofM = new TH1F("hChi2PerDofM", "matched track chi2/d.o.f.", 200, 0., 20.);
   fList->AddAtAndExpand(hChi2PerDofM, kChi2PerDofM);
-  TH1F *hChi2PerDofF = new TH1F("hChi2PerDofF", "fake track chi2/d.o.f.", 100, 0., 20.);
+  TH1F *hChi2PerDofF = new TH1F("hChi2PerDofF", "fake track chi2/d.o.f.", 200, 0., 20.);
   fList->AddAtAndExpand(hChi2PerDofF, kChi2PerDofF);
   
   // chi2 versus number of clusters
@@ -186,11 +202,11 @@ void AliAnalysisTaskMuonFakes::UserCreateOutputObjects()
   fList->AddAtAndExpand(hPtM, kPtM);
   TH1F *hPtF = new TH1F("hPtF", "fake track Pt distribution (GeV/c)", 100, 0., 20.);
   fList->AddAtAndExpand(hPtF, kPtF);
-  TH1F *hEta = new TH1F("hEta", "Muon pseudo-rapidity distribution", 100, -10., 0.);
+  TH1F *hEta = new TH1F("hEta", "Muon pseudo-rapidity distribution", 200, -10., 0.);
   fList->AddAtAndExpand(hEta , kEta );
-  TH1F *hEtaM = new TH1F("hEtaM", "matched track pseudo-rapidity distribution", 100, -10., 0.);
+  TH1F *hEtaM = new TH1F("hEtaM", "matched track pseudo-rapidity distribution", 200, -10., 0.);
   fList->AddAtAndExpand(hEtaM, kEtaM);
-  TH1F *hEtaF = new TH1F("hEtaF", "fake track pseudo-rapidity distribution", 100, -10., 0.);
+  TH1F *hEtaF = new TH1F("hEtaF", "fake track pseudo-rapidity distribution", 200, -10., 0.);
   fList->AddAtAndExpand(hEtaF, kEtaF);
   TH1F *hPhi = new TH1F("hPhi", "Muon phi distribution", 100, -1., 9.);
   fList->AddAtAndExpand(hPhi, kPhi);
@@ -198,12 +214,41 @@ void AliAnalysisTaskMuonFakes::UserCreateOutputObjects()
   fList->AddAtAndExpand(hPhiM, kPhiM);
   TH1F *hPhiF = new TH1F("hPhiF", "fake track phi distribution", 100, -1., 9.);
   fList->AddAtAndExpand(hPhiF, kPhiF);
-  TH1F *hDCA = new TH1F("hDCA", "Muon DCA distribution", 125, 0., 500.);
+  TH1F *hDCA = new TH1F("hDCA", "Muon DCA distribution", 250, 0., 500.);
   fList->AddAtAndExpand(hDCA, kDCA);
-  TH1F *hDCAM = new TH1F("hDCAM", "matched track DCA distribution", 125, 0., 500.);
+  TH1F *hDCAM = new TH1F("hDCAM", "matched track DCA distribution", 250, 0., 500.);
   fList->AddAtAndExpand(hDCAM, kDCAM);
-  TH1F *hDCAF = new TH1F("hDCAF", "fake track DCA distribution", 125, 0., 500.);
+  TH1F *hDCAF = new TH1F("hDCAF", "fake track DCA distribution", 250, 0., 500.);
   fList->AddAtAndExpand(hDCAF, kDCAF);
+  TH1F *hRAbs = new TH1F("hRAbs", "Muon R_{Abs} distribution", 200, 0., 100.);
+  fList->AddAtAndExpand(hRAbs, kRAbs);
+  TH1F *hRAbsM = new TH1F("hRAbsM", "matched track R_{Abs} distribution", 300, 0., 150.);
+  fList->AddAtAndExpand(hRAbsM, kRAbsM);
+  TH1F *hRAbsF = new TH1F("hRAbsF", "fake track R_{Abs} distribution", 200, 0., 100.);
+  fList->AddAtAndExpand(hRAbsF, kRAbsF);
+  
+  // track pair histograms
+  fList2 = new TObjArray(100);
+  fList2->SetOwner();
+  
+  // physics quantities of opposite-sign track pairs
+  TH1F* h = 0x0;
+  TString nameSuffix[4] = {"", "M", "F1", "F2"};
+  TString titlePrefix[4] = {"Dimuon", "matched-matched pair", "matched-fake pair", "fake-fake pair"};
+  for (Int_t i = 0; i < 4; i++) {
+    h = new TH1F(Form("h2Mass%s",nameSuffix[i].Data()), Form("%s mass distribution (GeV/c^{2})",titlePrefix[i].Data()), 300, 0., 15.);
+    fList2->AddAtAndExpand(h, k2Mass+i);
+    h = new TH1F(Form("h2P%s",nameSuffix[i].Data()), Form("%s P distribution (GeV/c)",titlePrefix[i].Data()), 100, 0., 200.);
+    fList2->AddAtAndExpand(h, k2P+i);
+    h = new TH1F(Form("h2Pt%s",nameSuffix[i].Data()), Form("%s Pt distribution (GeV/c)",titlePrefix[i].Data()), 100, 0., 20.);
+    fList2->AddAtAndExpand(h, k2Pt+i);
+    h = new TH1F(Form("h2Y%s",nameSuffix[i].Data()), Form("%s rapidity distribution",titlePrefix[i].Data()), 200, -10., 0.);
+    fList2->AddAtAndExpand(h, k2Y+i);
+    h = new TH1F(Form("h2Eta%s",nameSuffix[i].Data()), Form("%s pseudo-rapidity distribution",titlePrefix[i].Data()), 200, -10., 0.);
+    fList2->AddAtAndExpand(h, k2Eta+i);
+    h = new TH1F(Form("h2Phi%s",nameSuffix[i].Data()), Form("%s phi distribution",titlePrefix[i].Data()), 100, -1., 9.);
+    fList2->AddAtAndExpand(h, k2Phi+i);
+  }
   
   // global counters of tracks:
   // - reconstructible = number of reconstructible tracks
@@ -219,6 +264,8 @@ void AliAnalysisTaskMuonFakes::UserCreateOutputObjects()
   fTrackCounters->AddRubric("trig", "yes/no/unknown");
   fTrackCounters->AddRubric("selected", "yes/no");
   fTrackCounters->AddRubric("acc", "in/out/unknown");
+  TString centralityClasses = "5/10/15/20/25/30/35/40/45/50/55/60/65/70/75/80/85/90/95/100";
+  fTrackCounters->AddRubric("cent", centralityClasses.Data());
   fTrackCounters->Init();
   
   // detailled counters of fake tracks:
@@ -230,6 +277,7 @@ void AliAnalysisTaskMuonFakes::UserCreateOutputObjects()
   fFakeTrackCounters->AddRubric("trig", "yes/no/unknown");
   fFakeTrackCounters->AddRubric("selected", "yes/no");
   fFakeTrackCounters->AddRubric("acc", "in/out/unknown");
+  fFakeTrackCounters->AddRubric("cent", centralityClasses.Data());
   fFakeTrackCounters->Init();
   
   // counters of tracks matched by position or by using MC labels
@@ -240,6 +288,7 @@ void AliAnalysisTaskMuonFakes::UserCreateOutputObjects()
   fMatchedTrackCounters->AddRubric("trig", "yes/no");
   fMatchedTrackCounters->AddRubric("selected", "yes/no");
   fMatchedTrackCounters->AddRubric("acc", "in/out");
+  fMatchedTrackCounters->AddRubric("cent", centralityClasses.Data());
   fMatchedTrackCounters->Init();
   
   // global counters of events
@@ -254,8 +303,24 @@ void AliAnalysisTaskMuonFakes::UserCreateOutputObjects()
   fEventCounters->AddRubric("run", 1000000);
   fEventCounters->AddRubric("trig", "any/yes");
   fEventCounters->AddRubric("selected", "yes/no");
+  fEventCounters->AddRubric("cent", centralityClasses.Data());
   fEventCounters->Init();
   
+  // global counters of track pairs:
+  // - reconstructible = number of reconstructible track pairs
+  // - reconstructed   = number of reconstructed track pairs
+  // - matched         = number of reconstructed track pairs fully matched with a simulated one (reconstructible or not)
+  // - 1fake           = number of reconstructed track pairs made of one matched track and one fake
+  // - 2fakes          = number of reconstructed track pairs fully fake
+  fPairCounters = new AliCounterCollection(GetOutputSlot(7)->GetContainer()->GetName());
+  fPairCounters->AddRubric("pair", "reconstructible/reconstructed/matched/1fake/2fakes");
+  fPairCounters->AddRubric("run", 1000000);
+  fPairCounters->AddRubric("trig", "0/1/2");
+  fPairCounters->AddRubric("selected", "yes/no");
+  fPairCounters->AddRubric("acc", "in/out/unknown");
+  fPairCounters->AddRubric("cent", centralityClasses.Data());
+  fPairCounters->Init();
+  
   // Disable printout of AliMCEvent
   AliLog::SetClassDebugLevel("AliMCEvent",-1);
   
@@ -265,6 +330,8 @@ void AliAnalysisTaskMuonFakes::UserCreateOutputObjects()
   PostData(3, fFakeTrackCounters);
   PostData(4, fMatchedTrackCounters);
   PostData(5, fEventCounters);
+  PostData(6, fList2);
+  PostData(7, fPairCounters);
 }
 
 //________________________________________________________________________
@@ -291,6 +358,21 @@ void AliAnalysisTaskMuonFakes::UserExec(Option_t *)
     return;
   }      
   
+  // current centrality class
+  TString centrality = "cent:";
+  Double_t centralityValue = esd->GetCentrality()->GetCentralityPercentile("V0M");
+  TObjArray* centralylimits = fTrackCounters->GetKeyWords("cent").Tokenize(",");
+  TObjString* limit = 0x0;
+  TIter nextLimit(centralylimits);
+  while ((limit = static_cast<TObjString*>(nextLimit()))) {
+    if (centralityValue < limit->String().Atoi()) {
+      centrality += limit->GetName();
+      break;
+    }
+  }
+  if (!limit) centrality += static_cast<TObjString*>(centralylimits->Last())->GetName();
+  delete centralylimits;
+  
   // Load MC event 
   AliMCEventHandler *mcH = 0;
   if(MCEvent()) mcH = static_cast<AliMCEventHandler*>((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler()); 
@@ -302,16 +384,35 @@ void AliAnalysisTaskMuonFakes::UserExec(Option_t *)
   AliMUONVTriggerTrackStore* triggerTrackRefStore = rc.TriggerableTracks(-1);
   if (!muonTrackStore || !trackRefStore) return;
   
-  // count the number of reconstructible tracks
+  // loop over trackRefs
+  Int_t nMuPlus[2] = {0, 0};
+  Int_t nMuMinus[2] = {0, 0};
   TIter next(trackRefStore->CreateIterator());
   AliMUONTrack* trackRef;
   while ( ( trackRef = static_cast<AliMUONTrack*>(next()) ) ) {
-    if (trackRef->IsValid(fRequestedStationMask, fRequest2ChInSameSt45)) {
-      TString trig = (triggerTrackRefStore->FindObject(trackRef->GetUniqueID())) ? "trig:yes" : "trig:no";
-      fTrackCounters->Count(Form("track:reconstructible/run:%d/%s/%s/acc:unknown", fCurrentRunNumber, trig.Data(), selected.Data()));
-    }
+    
+    // skip trackRefs that are not reconstructible
+    if (!trackRef->IsValid(fRequestedStationMask, fRequest2ChInSameSt45)) continue;
+    
+    // trigger condition
+    Bool_t trigger = (triggerTrackRefStore->FindObject(trackRef->GetUniqueID()));
+    Int_t iTrig = trigger ? 1 : 0;
+    TString trig = trigger ? "trig:yes" : "trig:no";
+    
+    // count muons
+    if (trackRef->GetTrackParamAtVertex()->GetCharge() > 0) nMuPlus[iTrig]++;
+    else nMuMinus[iTrig]++;
+    
+    // fill global counters
+    fTrackCounters->Count(Form("track:reconstructible/run:%d/%s/%s/acc:unknown/%s", fCurrentRunNumber, trig.Data(), selected.Data(), centrality.Data()));
+    
   }
   
+  // fill global counters
+  fPairCounters->Count(Form("pair:reconstructible/run:%d/trig:0/%s/acc:unknown/%s", fCurrentRunNumber, selected.Data(), centrality.Data()), nMuPlus[0]*nMuMinus[0]);
+  fPairCounters->Count(Form("pair:reconstructible/run:%d/trig:1/%s/acc:unknown/%s", fCurrentRunNumber, selected.Data(), centrality.Data()), nMuPlus[1]*nMuMinus[0]+nMuPlus[0]*nMuMinus[1]);
+  fPairCounters->Count(Form("pair:reconstructible/run:%d/trig:2/%s/acc:unknown/%s", fCurrentRunNumber, selected.Data(), centrality.Data()), nMuPlus[1]*nMuMinus[1]);
+  
   // loop over ESD tracks
   Int_t nTrackerTracks = 0;
   Bool_t containTrack[2] = {kFALSE, kFALSE};
@@ -319,32 +420,29 @@ void AliAnalysisTaskMuonFakes::UserExec(Option_t *)
   Bool_t containMatchedYetTrack[2] = {kFALSE, kFALSE};
   AliMUONVTrackStore *fakeTrackStore = AliMUONESDInterface::NewTrackStore();
   Int_t nTracks = (Int_t)esd->GetNumberOfMuonTracks() ;
-  for (Int_t iTrack = 0; iTrack <  nTracks;  iTrack++) {
+  for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
     
     AliESDMuonTrack* esdTrack = esd->GetMuonTrack(iTrack);
     
     // skip ghosts
     if (!esdTrack->ContainTrackerData()) continue;
-    nTrackerTracks++;
     containTrack[0] = kTRUE;
     
     // trigger condition
     Bool_t trigger = esdTrack->ContainTriggerData();
-    TString trig = "trig:";
-    if (trigger) {
-      trig += "yes";
-      containTrack[1] = kTRUE;
-    } else trig += "no";
+    TString trig = trigger ? "trig:yes" : "trig:no";
+    if (trigger) containTrack[1] = kTRUE;
     
     // acceptance condition
-    TString acc = "acc:";
-    Double_t thetaTrackAbsEnd = TMath::ATan(esdTrack->GetRAtAbsorberEnd()/505.) * TMath::RadToDeg();
+    Double_t rAbs = esdTrack->GetRAtAbsorberEnd();
+    Double_t thetaTrackAbsEnd = TMath::ATan(rAbs/505.) * TMath::RadToDeg();
     Double_t eta = esdTrack->Eta();
-    if (thetaTrackAbsEnd >= 2. && thetaTrackAbsEnd <= 9. && eta >= -4. && eta <= -2.5) acc += "in";
-    else acc += "out";
+    Bool_t inAcc = (thetaTrackAbsEnd >= 2. && thetaTrackAbsEnd <= 10. && eta >= -4. && eta <= -2.5);
+    TString acc = inAcc ? "acc:in" : "acc:out";
     
     // fill global counters
-    fTrackCounters->Count(Form("track:reconstructed/run:%d/%s/%s/%s", fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data()));
+    if ((!fMatchTrig || trigger) && (!fApplyAccCut || inAcc)) nTrackerTracks++;
+    fTrackCounters->Count(Form("track:reconstructed/run:%d/%s/%s/%s/%s", fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data(), centrality.Data()));
     
     // find the corresponding MUON track
     AliMUONTrack* muonTrack = static_cast<AliMUONTrack*>(muonTrackStore->FindObject(esdTrack->GetUniqueID()));
@@ -360,16 +458,19 @@ void AliAnalysisTaskMuonFakes::UserExec(Option_t *)
     Double_t dca = esdTrack->GetDCA();
     
     // fill global histograms
-    ((TH1F*)fList->UncheckedAt(kNumberOfClusters))->Fill(nClusters);
-    ((TH1F*)fList->UncheckedAt(kNumberOfChamberHit))->Fill(nChamberHit);
-    ((TH1F*)fList->UncheckedAt(kChi2PerDof))->Fill(normalizedChi2);
-    ((TH1F*)fList->UncheckedAt(kP))->Fill(p);
-    ((TH1F*)fList->UncheckedAt(kPt))->Fill(pT);
-    ((TH1F*)fList->UncheckedAt(kEta))->Fill(eta);
-    ((TH1F*)fList->UncheckedAt(kPhi))->Fill(phi);
-    ((TH1F*)fList->UncheckedAt(kDCA))->Fill(dca);
-    ((TH1F*)fList->UncheckedAt(kChi2PerDofVsNClusters))->Fill(nClusters,normalizedChi2);
-    ((TH1F*)fList->UncheckedAt(kChi2PerDofVsNChamberHit))->Fill(nChamberHit,normalizedChi2);
+    if ((!fMatchTrig || trigger) && (!fApplyAccCut || inAcc)) {
+      ((TH1F*)fList->UncheckedAt(kNumberOfClusters))->Fill(nClusters);
+      ((TH1F*)fList->UncheckedAt(kNumberOfChamberHit))->Fill(nChamberHit);
+      ((TH1F*)fList->UncheckedAt(kChi2PerDof))->Fill(normalizedChi2);
+      ((TH1F*)fList->UncheckedAt(kP))->Fill(p);
+      ((TH1F*)fList->UncheckedAt(kPt))->Fill(pT);
+      ((TH1F*)fList->UncheckedAt(kEta))->Fill(eta);
+      ((TH1F*)fList->UncheckedAt(kPhi))->Fill(phi);
+      ((TH1F*)fList->UncheckedAt(kDCA))->Fill(dca);
+      ((TH1F*)fList->UncheckedAt(kRAbs))->Fill(rAbs);
+      ((TH1F*)fList->UncheckedAt(kChi2PerDofVsNClusters))->Fill(nClusters,normalizedChi2);
+      ((TH1F*)fList->UncheckedAt(kChi2PerDofVsNChamberHit))->Fill(nChamberHit,normalizedChi2);
+    }
     
     // try to match, by position, the reconstructed track with a simulated one
     Int_t nMatchClustersByPosition = 0;
@@ -387,10 +488,10 @@ void AliAnalysisTaskMuonFakes::UserExec(Option_t *)
     if (MCLabelByLabel >= 0 && MCLabelByPosition >= 0 && MCLabelByLabel != MCLabelByPosition) labelCase += "match other";
     else if (MCLabelByLabel >= 0) labelCase += "match";
     else labelCase += "not match";
-    fMatchedTrackCounters->Count(Form("%s/%s/run:%d/%s/%s/%s", positionCase.Data(), labelCase.Data(), fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data()));
+    fMatchedTrackCounters->Count(Form("%s/%s/run:%d/%s/%s/%s/%s", positionCase.Data(), labelCase.Data(), fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data(), centrality.Data()));
     if ((MCLabelByLabel >= 0 && MCLabelByPosition < 0) || (MCLabelByLabel < 0 && MCLabelByPosition >= 0))
-      fFakeTrackCounters->Count(Form("track:fake?/run:%d/file:%s/event:%d/%s/%s/%s", fCurrentRunNumber, fCurrentFileName.Data(),
-                                    esd->GetEventNumberInFile(), trig.Data(), selected.Data(), acc.Data()));
+      fFakeTrackCounters->Count(Form("track:fake?/run:%d/file:%s/event:%d/%s/%s/%s/%s", fCurrentRunNumber, fCurrentFileName.Data(),
+                                    esd->GetEventNumberInFile(), trig.Data(), selected.Data(), acc.Data(), centrality.Data()));
 
     // take actions according to the matching result we are interested in
     Int_t nMatchClusters = (fUseLabel) ? nMatchClustersByLabel : nMatchClustersByPosition;
@@ -398,7 +499,7 @@ void AliAnalysisTaskMuonFakes::UserExec(Option_t *)
     if (matchedTrackRef) {
       
       // fill global counters
-      fTrackCounters->Count(Form("track:matched/run:%d/%s/%s/%s", fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data()));
+      fTrackCounters->Count(Form("track:matched/run:%d/%s/%s/%s/%s", fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data(), centrality.Data()));
       
       // track matched with a trackRef that is not reconstructible
       if (!matchedTrackRef->IsValid(fRequestedStationMask, fRequest2ChInSameSt45)) {
@@ -407,24 +508,30 @@ void AliAnalysisTaskMuonFakes::UserExec(Option_t *)
        if (trigger) containMatchedYetTrack[1] = kTRUE;
        
        // fill global counters
-       fTrackCounters->Count(Form("track:matchedyet/run:%d/%s/%s/%s", fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data()));
-       fFakeTrackCounters->Count(Form("track:matchedyet/run:%d/file:%s/event:%d/%s/%s/%s", fCurrentRunNumber, fCurrentFileName.Data(),
-                                      esd->GetEventNumberInFile(), trig.Data(), selected.Data(), acc.Data()));
+       fTrackCounters->Count(Form("track:matchedyet/run:%d/%s/%s/%s/%s", fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data(), centrality.Data()));
+       fFakeTrackCounters->Count(Form("track:matchedyet/run:%d/file:%s/event:%d/%s/%s/%s/%s", fCurrentRunNumber, fCurrentFileName.Data(),
+                                      esd->GetEventNumberInFile(), trig.Data(), selected.Data(), acc.Data(), centrality.Data()));
       }
       
       // fill histograms
-      ((TH1F*)fList->UncheckedAt(kFractionOfMatchedClusters))->Fill(((Float_t) nMatchClusters) / ((Float_t) nClusters));
-      ((TH1F*)fList->UncheckedAt(kNumberOfClustersMC))->Fill(matchedTrackRef->GetNClusters());
-      ((TH1F*)fList->UncheckedAt(kNumberOfClustersM))->Fill(nClusters);
-      ((TH1F*)fList->UncheckedAt(kNumberOfChamberHitM))->Fill(nChamberHit);
-      ((TH1F*)fList->UncheckedAt(kChi2PerDofM))->Fill(normalizedChi2);
-      ((TH1F*)fList->UncheckedAt(kPM))->Fill(p);
-      ((TH1F*)fList->UncheckedAt(kPtM))->Fill(pT);
-      ((TH1F*)fList->UncheckedAt(kEtaM))->Fill(eta);
-      ((TH1F*)fList->UncheckedAt(kPhiM))->Fill(phi);
-      ((TH1F*)fList->UncheckedAt(kDCAM))->Fill(dca);
-      ((TH1F*)fList->UncheckedAt(kChi2PerDofVsNClustersM))->Fill(nClusters,normalizedChi2);
-      ((TH1F*)fList->UncheckedAt(kChi2PerDofVsNChamberHitM))->Fill(nChamberHit,normalizedChi2);
+      if ((!fMatchTrig || trigger) && (!fApplyAccCut || inAcc)) {
+       ((TH1F*)fList->UncheckedAt(kFractionOfMatchedClusters))->Fill(((Float_t) nMatchClusters) / ((Float_t) nClusters));
+       ((TH1F*)fList->UncheckedAt(kNumberOfClustersMC))->Fill(matchedTrackRef->GetNClusters());
+       ((TH1F*)fList->UncheckedAt(kNumberOfClustersM))->Fill(nClusters);
+       ((TH1F*)fList->UncheckedAt(kNumberOfChamberHitM))->Fill(nChamberHit);
+       ((TH1F*)fList->UncheckedAt(kChi2PerDofM))->Fill(normalizedChi2);
+       ((TH1F*)fList->UncheckedAt(kPM))->Fill(p);
+       ((TH1F*)fList->UncheckedAt(kPtM))->Fill(pT);
+       ((TH1F*)fList->UncheckedAt(kEtaM))->Fill(eta);
+       ((TH1F*)fList->UncheckedAt(kPhiM))->Fill(phi);
+       ((TH1F*)fList->UncheckedAt(kDCAM))->Fill(dca);
+       ((TH1F*)fList->UncheckedAt(kRAbsM))->Fill(rAbs);
+       ((TH1F*)fList->UncheckedAt(kChi2PerDofVsNClustersM))->Fill(nClusters,normalizedChi2);
+       ((TH1F*)fList->UncheckedAt(kChi2PerDofVsNChamberHitM))->Fill(nChamberHit,normalizedChi2);
+      }
+      
+      // flag matched tracks
+      esdTrack->SetLabel(matchedTrackRef->GetUniqueID());
       
       // remove already matched trackRefs
       trackRefStore->Remove(*matchedTrackRef);
@@ -435,21 +542,27 @@ void AliAnalysisTaskMuonFakes::UserExec(Option_t *)
       if (trigger) containFakeTrack[1] = kTRUE;
       
       // fill global counters
-      fTrackCounters->Count(Form("track:fake/run:%d/%s/%s/%s", fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data()));
-      fFakeTrackCounters->Count(Form("track:fake/run:%d/file:%s/event:%d/%s/%s/%s", fCurrentRunNumber, fCurrentFileName.Data(),
-                                    esd->GetEventNumberInFile(), trig.Data(), selected.Data(), acc.Data()));
+      fTrackCounters->Count(Form("track:fake/run:%d/%s/%s/%s/%s", fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data(), centrality.Data()));
+      fFakeTrackCounters->Count(Form("track:fake/run:%d/file:%s/event:%d/%s/%s/%s/%s", fCurrentRunNumber, fCurrentFileName.Data(),
+                                    esd->GetEventNumberInFile(), trig.Data(), selected.Data(), acc.Data(), centrality.Data()));
       
       // fill histograms
-      ((TH1F*)fList->UncheckedAt(kNumberOfClustersF))->Fill(nClusters);
-      ((TH1F*)fList->UncheckedAt(kNumberOfChamberHitF))->Fill(nChamberHit);
-      ((TH1F*)fList->UncheckedAt(kChi2PerDofF))->Fill(normalizedChi2);
-      ((TH1F*)fList->UncheckedAt(kPF))->Fill(p);
-      ((TH1F*)fList->UncheckedAt(kPtF))->Fill(pT);
-      ((TH1F*)fList->UncheckedAt(kEtaF))->Fill(eta);
-      ((TH1F*)fList->UncheckedAt(kPhiF))->Fill(phi);
-      ((TH1F*)fList->UncheckedAt(kDCAF))->Fill(dca);
-      ((TH1F*)fList->UncheckedAt(kChi2PerDofVsNClustersF))->Fill(nClusters,normalizedChi2);
-      ((TH1F*)fList->UncheckedAt(kChi2PerDofVsNChamberHitF))->Fill(nChamberHit,normalizedChi2);
+      if ((!fMatchTrig || trigger) && (!fApplyAccCut || inAcc)) {
+       ((TH1F*)fList->UncheckedAt(kNumberOfClustersF))->Fill(nClusters);
+       ((TH1F*)fList->UncheckedAt(kNumberOfChamberHitF))->Fill(nChamberHit);
+       ((TH1F*)fList->UncheckedAt(kChi2PerDofF))->Fill(normalizedChi2);
+       ((TH1F*)fList->UncheckedAt(kPF))->Fill(p);
+       ((TH1F*)fList->UncheckedAt(kPtF))->Fill(pT);
+       ((TH1F*)fList->UncheckedAt(kEtaF))->Fill(eta);
+       ((TH1F*)fList->UncheckedAt(kPhiF))->Fill(phi);
+       ((TH1F*)fList->UncheckedAt(kDCAF))->Fill(dca);
+       ((TH1F*)fList->UncheckedAt(kRAbsF))->Fill(rAbs);
+       ((TH1F*)fList->UncheckedAt(kChi2PerDofVsNClustersF))->Fill(nClusters,normalizedChi2);
+       ((TH1F*)fList->UncheckedAt(kChi2PerDofVsNChamberHitF))->Fill(nChamberHit,normalizedChi2);
+      }
+      
+      // flag fake tracks
+      esdTrack->SetLabel(-1);
       
       // store fake tracks
       fakeTrackStore->Add(*muonTrack);
@@ -460,23 +573,23 @@ void AliAnalysisTaskMuonFakes::UserExec(Option_t *)
   
   // fill histogram and global counters
   ((TH1F*)fList->UncheckedAt(kNumberOfTracks))->Fill(nTrackerTracks);
-  if (containTrack[0]) fEventCounters->Count(Form("event:any/run:%d/trig:any/%s", fCurrentRunNumber, selected.Data()));
-  if (containTrack[1]) fEventCounters->Count(Form("event:any/run:%d/trig:yes/%s", fCurrentRunNumber, selected.Data()));
-  if (containFakeTrack[0]) fEventCounters->Count(Form("event:fake/run:%d/trig:any/%s", fCurrentRunNumber, selected.Data()));
-  if (containFakeTrack[1]) fEventCounters->Count(Form("event:fake/run:%d/trig:yes/%s", fCurrentRunNumber, selected.Data()));
-  if (containMatchedYetTrack[0]) fEventCounters->Count(Form("event:matchedyet/run:%d/trig:any/%s", fCurrentRunNumber, selected.Data()));
-  if (containMatchedYetTrack[1]) fEventCounters->Count(Form("event:matchedyet/run:%d/trig:yes/%s", fCurrentRunNumber, selected.Data()));
+  if (containTrack[0]) fEventCounters->Count(Form("event:any/run:%d/trig:any/%s/%s", fCurrentRunNumber, selected.Data(), centrality.Data()));
+  if (containTrack[1]) fEventCounters->Count(Form("event:any/run:%d/trig:yes/%s/%s", fCurrentRunNumber, selected.Data(), centrality.Data()));
+  if (containFakeTrack[0]) fEventCounters->Count(Form("event:fake/run:%d/trig:any/%s/%s", fCurrentRunNumber, selected.Data(), centrality.Data()));
+  if (containFakeTrack[1]) fEventCounters->Count(Form("event:fake/run:%d/trig:yes/%s/%s", fCurrentRunNumber, selected.Data(), centrality.Data()));
+  if (containMatchedYetTrack[0]) fEventCounters->Count(Form("event:matchedyet/run:%d/trig:any/%s/%s", fCurrentRunNumber, selected.Data(), centrality.Data()));
+  if (containMatchedYetTrack[1]) fEventCounters->Count(Form("event:matchedyet/run:%d/trig:yes/%s/%s", fCurrentRunNumber, selected.Data(), centrality.Data()));
   
   // count the number of not connected and additional fake tracks
   if (fakeTrackStore->GetSize() > 0) {
     
     // remove the most connected fake tracks
-    Int_t nFreeMissingTracks = RemoveConnectedFakes(*fakeTrackStore, *trackRefStore);
+    Int_t nFreeMissingTracks = RemoveConnectedFakes(*fakeTrackStore, *trackRefStore, selected, centrality);
     
     if (fakeTrackStore->GetSize() > 0) {
       
       // fill global counters
-      fEventCounters->Count(Form("event:notconnected/run:%d/trig:any/%s", fCurrentRunNumber, selected.Data()));
+      fEventCounters->Count(Form("event:notconnected/run:%d/trig:any/%s/%s", fCurrentRunNumber, selected.Data(), centrality.Data()));
       
       // check status of remaining fakes with respect to the matching with trigger
       Bool_t containMatchedFake = kFALSE;
@@ -489,7 +602,7 @@ void AliAnalysisTaskMuonFakes::UserExec(Option_t *)
       }
       
       // fill global counters
-      if (containMatchedFake) fEventCounters->Count(Form("event:notconnected/run:%d/trig:yes/%s", fCurrentRunNumber, selected.Data()));
+      if (containMatchedFake) fEventCounters->Count(Form("event:notconnected/run:%d/trig:yes/%s/%s", fCurrentRunNumber, selected.Data(), centrality.Data()));
       
       // remove the remaining free reconstructible tracks
       Int_t nAdditionalTracks = fakeTrackStore->GetSize() - nFreeMissingTracks;
@@ -498,18 +611,18 @@ void AliAnalysisTaskMuonFakes::UserExec(Option_t *)
        
        // fill histogram and global counters
        ((TH1F*)fList->UncheckedAt(kNumberOfAdditionalTracks))->Fill(nAdditionalTracks);
-       fEventCounters->Count(Form("event:additional/run:%d/trig:any/%s", fCurrentRunNumber, selected.Data()));
+       fEventCounters->Count(Form("event:additional/run:%d/trig:any/%s/%s", fCurrentRunNumber, selected.Data(), centrality.Data()));
        if (!containUnmatchedFake) { // all matched
-         fTrackCounters->Count(Form("track:additional/run:%d/trig:yes/%s/acc:unknown", fCurrentRunNumber, selected.Data()), nAdditionalTracks);
-         fFakeTrackCounters->Count(Form("track:additional/run:%d/file:%s/event:%d/trig:yes/%s/acc:unknown", fCurrentRunNumber, fCurrentFileName.Data(), esd->GetEventNumberInFile(), selected.Data()), nAdditionalTracks);
-         fEventCounters->Count(Form("event:additional/run:%d/trig:yes/%s", fCurrentRunNumber, selected.Data()));
+         fTrackCounters->Count(Form("track:additional/run:%d/trig:yes/%s/acc:unknown/%s", fCurrentRunNumber, selected.Data(), centrality.Data()), nAdditionalTracks);
+         fFakeTrackCounters->Count(Form("track:additional/run:%d/file:%s/event:%d/trig:yes/%s/acc:unknown/%s", fCurrentRunNumber, fCurrentFileName.Data(), esd->GetEventNumberInFile(), selected.Data(), centrality.Data()), nAdditionalTracks);
+         fEventCounters->Count(Form("event:additional/run:%d/trig:yes/%s/%s", fCurrentRunNumber, selected.Data(), centrality.Data()));
        } else if (!containMatchedFake) { // none matched
-         fTrackCounters->Count(Form("track:additional/run:%d/trig:no/%s/acc:unknown", fCurrentRunNumber, selected.Data()), nAdditionalTracks);
-         fFakeTrackCounters->Count(Form("track:additional/run:%d/file:%s/event:%d/trig:no/%s/acc:unknown", fCurrentRunNumber, fCurrentFileName.Data(), esd->GetEventNumberInFile(), selected.Data()), nAdditionalTracks);
+         fTrackCounters->Count(Form("track:additional/run:%d/trig:no/%s/acc:unknown/%s", fCurrentRunNumber, selected.Data(), centrality.Data()), nAdditionalTracks);
+         fFakeTrackCounters->Count(Form("track:additional/run:%d/file:%s/event:%d/trig:no/%s/acc:unknown/%s", fCurrentRunNumber, fCurrentFileName.Data(), esd->GetEventNumberInFile(), selected.Data(), centrality.Data()), nAdditionalTracks);
        } else { // mixed
-         fTrackCounters->Count(Form("track:additional/run:%d/trig:unknown/%s/acc:unknown", fCurrentRunNumber, selected.Data()), nAdditionalTracks);
-         fFakeTrackCounters->Count(Form("track:additional/run:%d/file:%s/event:%d/trig:unknown/%s/acc:unknown", fCurrentRunNumber, fCurrentFileName.Data(), esd->GetEventNumberInFile(), selected.Data()), nAdditionalTracks);
-         fEventCounters->Count(Form("event:additional/run:%d/trig:yes/%s", fCurrentRunNumber, selected.Data()));
+         fTrackCounters->Count(Form("track:additional/run:%d/trig:unknown/%s/acc:unknown/%s", fCurrentRunNumber, selected.Data(), centrality.Data()), nAdditionalTracks);
+         fFakeTrackCounters->Count(Form("track:additional/run:%d/file:%s/event:%d/trig:unknown/%s/acc:unknown/%s", fCurrentRunNumber, fCurrentFileName.Data(), esd->GetEventNumberInFile(), selected.Data(), centrality.Data()), nAdditionalTracks);
+         fEventCounters->Count(Form("event:additional/run:%d/trig:yes/%s/%s", fCurrentRunNumber, selected.Data(), centrality.Data()));
        }
        
       }
@@ -518,14 +631,134 @@ void AliAnalysisTaskMuonFakes::UserExec(Option_t *)
     
   }
   
+  // clean memory
   delete fakeTrackStore;
   
+  // double loop over ESD tracks, build pairs and fill histograms and counters according to their label
+  TLorentzVector vMu1, vMu2, vDiMu;
+  for (Int_t iTrack1 = 0; iTrack1 < nTracks; iTrack1++) {
+    AliESDMuonTrack* muonTrack1 = esd->GetMuonTrack(iTrack1);
+    
+    // skip ghosts
+    if (!muonTrack1->ContainTrackerData()) continue;
+    
+    // get track info
+    Bool_t trigger1 = muonTrack1->ContainTriggerData();
+    Double_t thetaAbs1 = TMath::ATan(muonTrack1->GetRAtAbsorberEnd()/505.) * TMath::RadToDeg();
+    Double_t eta1 = muonTrack1->Eta();
+    Bool_t acc1 = (thetaAbs1 >= 2. && thetaAbs1 <= 10. && eta1 >= -4. && eta1 <= -2.5);
+    Short_t charge1 = muonTrack1->Charge();
+    Int_t label1 = muonTrack1->GetLabel();
+    muonTrack1->LorentzP(vMu1);
+    
+    for (Int_t iTrack2 = iTrack1+1; iTrack2 < nTracks; iTrack2++) {
+      AliESDMuonTrack* muonTrack2 = esd->GetMuonTrack(iTrack2);
+      
+      // skip ghosts
+      if (!muonTrack2->ContainTrackerData()) continue;
+      
+      // keep only opposite sign pairs
+      Short_t charge2 = muonTrack2->Charge();
+      if (charge1*charge2 > 0) continue;
+      
+      // get track info
+      Bool_t trigger2 = muonTrack2->ContainTriggerData();
+      Double_t thetaAbs2 = TMath::ATan(muonTrack2->GetRAtAbsorberEnd()/505.) * TMath::RadToDeg();
+      Double_t eta2 = muonTrack2->Eta();
+      Bool_t acc2 = (thetaAbs2 >= 2. && thetaAbs2 <= 10. && eta2 >= -4. && eta2 <= -2.5);
+      Int_t label2 = muonTrack2->GetLabel();
+      muonTrack2->LorentzP(vMu2);
+      
+      // compute kinematics of the pair
+      vDiMu = vMu1 + vMu2;
+      Float_t mass = vDiMu.M();
+      Float_t p = vDiMu.P();
+      Float_t pt = vDiMu.Pt();
+      Float_t y = vDiMu.Rapidity();
+      Float_t eta = vDiMu.Eta();
+      Float_t phi = vDiMu.Phi();
+      if (phi < 0) phi += 2.*TMath::Pi();
+      
+      // trigger condition
+      TString trig = "trig:";
+      if (trigger1 && trigger2) trig += "2";
+      else if (trigger1 || trigger2) trig += "1";
+      else trig += "0";
+      
+      // acceptance condition
+      Bool_t inAcc = (acc1 && acc2 && eta >= -4. && eta <= -2.5);
+      TString acc = inAcc ? "acc:in" : "acc:out";
+      
+      // fill global histograms
+      if ((!fMatchTrig || (trigger1 && trigger2)) && (!fApplyAccCut || inAcc)) {
+       ((TH1F*)fList2->UncheckedAt(k2Mass))->Fill(mass);
+       ((TH1F*)fList2->UncheckedAt(k2P))->Fill(p);
+       ((TH1F*)fList2->UncheckedAt(k2Pt))->Fill(pt);
+       ((TH1F*)fList2->UncheckedAt(k2Y))->Fill(y);
+       ((TH1F*)fList2->UncheckedAt(k2Eta))->Fill(eta);
+       ((TH1F*)fList2->UncheckedAt(k2Phi))->Fill(phi);
+      }
+      
+      TString pair = "pair:";
+      
+      // fill histograms according to labels
+      if (label1 >= 0 && label2 >= 0) {
+       
+       pair += "matched";
+       
+       if ((!fMatchTrig || (trigger1 && trigger2)) && (!fApplyAccCut || inAcc)) {
+         ((TH1F*)fList2->UncheckedAt(k2MassM))->Fill(mass);
+         ((TH1F*)fList2->UncheckedAt(k2PM))->Fill(p);
+         ((TH1F*)fList2->UncheckedAt(k2PtM))->Fill(pt);
+         ((TH1F*)fList2->UncheckedAt(k2YM))->Fill(y);
+         ((TH1F*)fList2->UncheckedAt(k2EtaM))->Fill(eta);
+         ((TH1F*)fList2->UncheckedAt(k2PhiM))->Fill(phi);
+       }
+       
+      } else if (label1 >= 0 || label2 >= 0) {
+       
+       pair += "1fake";
+       
+       if ((!fMatchTrig || (trigger1 && trigger2)) && (!fApplyAccCut || inAcc)) {
+         ((TH1F*)fList2->UncheckedAt(k2MassF1))->Fill(mass);
+         ((TH1F*)fList2->UncheckedAt(k2PF1))->Fill(p);
+         ((TH1F*)fList2->UncheckedAt(k2PtF1))->Fill(pt);
+         ((TH1F*)fList2->UncheckedAt(k2YF1))->Fill(y);
+         ((TH1F*)fList2->UncheckedAt(k2EtaF1))->Fill(eta);
+         ((TH1F*)fList2->UncheckedAt(k2PhiF1))->Fill(phi);
+       }
+       
+      } else {
+       
+       pair += "2fakes";
+       
+       if ((!fMatchTrig || (trigger1 && trigger2)) && (!fApplyAccCut || inAcc)) {
+         ((TH1F*)fList2->UncheckedAt(k2MassF2))->Fill(mass);
+         ((TH1F*)fList2->UncheckedAt(k2PF2))->Fill(p);
+         ((TH1F*)fList2->UncheckedAt(k2PtF2))->Fill(pt);
+         ((TH1F*)fList2->UncheckedAt(k2YF2))->Fill(y);
+         ((TH1F*)fList2->UncheckedAt(k2EtaF2))->Fill(eta);
+         ((TH1F*)fList2->UncheckedAt(k2PhiF2))->Fill(phi);
+       }
+       
+      }
+      
+      // fill global counters
+      fPairCounters->Count(Form("pair:reconstructed/run:%d/%s/%s/%s/%s", fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data(), centrality.Data()));
+      fPairCounters->Count(Form("%s/run:%d/%s/%s/%s/%s", pair.Data(), fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data(), centrality.Data()));
+      
+    }
+    
+  }
+  
   // Post final data
   PostData(1, fList);
   PostData(2, fTrackCounters);
   PostData(3, fFakeTrackCounters);
   PostData(4, fMatchedTrackCounters);
   PostData(5, fEventCounters);
+  PostData(6, fList2);
+  PostData(7, fPairCounters);
 }
 
 //________________________________________________________________________
@@ -563,24 +796,28 @@ void AliAnalysisTaskMuonFakes::Terminate(Option_t *)
   
   // recover output objects
   fList = static_cast<TObjArray*> (GetOutputData(1));
-  if (!fList) return;
+  fList2 = static_cast<TObjArray*> (GetOutputData(6));
+  if (!fList || !fList2) return;
   fTrackCounters = static_cast<AliCounterCollection*> (GetOutputData(2));
   fFakeTrackCounters = static_cast<AliCounterCollection*> (GetOutputData(3));
   fMatchedTrackCounters = static_cast<AliCounterCollection*> (GetOutputData(4));
   fEventCounters = static_cast<AliCounterCollection*> (GetOutputData(5));
+  fPairCounters = static_cast<AliCounterCollection*> (GetOutputData(7));
   
   // add canvas to compare histograms
-  fCanvases = new TObjArray(1000);
+  fCanvases = new TObjArray(3);
   fCanvases->SetOwner();
-  TCanvas *cFakesSummary1 = new TCanvas("cFakesSummary1","cFakesSummary1",1200,600);
+  TCanvas *cFakesSummary1 = new TCanvas("cFakesSummary1","cFakesSummary1",900,900);
   fCanvases->AddAtAndExpand(cFakesSummary1, 0);
   TCanvas *cFakesSummary2 = new TCanvas("cFakesSummary2","cFakesSummary2",1200,600);
   fCanvases->AddAtAndExpand(cFakesSummary2, 1);
+  TCanvas *cFakesSummary3 = new TCanvas("cFakesSummary3","cFakesSummary3",900,600);
+  fCanvases->AddAtAndExpand(cFakesSummary3, 2);
   
   // display
-  Int_t iHist1[8] = {kNumberOfClusters, kChi2PerDof, kP, kEta, kNumberOfChamberHit, kDCA, kPt, kPhi};
-  cFakesSummary1->Divide(4,2);
-  for (Int_t i=0; i<8; i++) {
+  Int_t iHist1[9] = {kNumberOfClusters, kNumberOfChamberHit, kChi2PerDof, kDCA, kRAbs, kEta, kP, kPt, kPhi};
+  cFakesSummary1->Divide(3,3);
+  for (Int_t i=0; i<9; i++) {
     cFakesSummary1->cd(i+1);
     cFakesSummary1->GetPad(i+1)->SetLogy();
     ((TH1F*)fList->UncheckedAt(iHist1[i]))->SetMinimum(0.5);
@@ -604,6 +841,26 @@ void AliAnalysisTaskMuonFakes::Terminate(Option_t *)
     ((TH2F*)fList->UncheckedAt(iHist2[i]+2))->DrawCopy("sames");
   }
   
+  Int_t iHist3[6] = {k2Mass, k2P, k2Pt, k2Y, k2Eta, k2Phi};
+  cFakesSummary3->Divide(3,2);
+  for (Int_t i=0; i<6; i++) {
+    cFakesSummary3->cd(i+1);
+    cFakesSummary3->GetPad(i+1)->SetLogy();
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]))->SetMinimum(0.5);
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]))->DrawCopy();
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]+1))->SetLineColor(4);
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]+1))->DrawCopy("sames");
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]+2))->Add((TH1F*)fList2->UncheckedAt(iHist3[i]+3));
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]+2))->SetLineColor(2);
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]+2))->SetFillColor(2);
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]+2))->SetFillStyle(3017);
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]+2))->DrawCopy("sames");
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]+3))->SetLineColor(6);
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]+3))->SetFillColor(6);
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]+3))->SetFillStyle(3018);
+    ((TH1F*)fList2->UncheckedAt(iHist3[i]+3))->DrawCopy("sames");
+  }
+  
   // print
   if (fTrackCounters && fFakeTrackCounters && fMatchedTrackCounters && fEventCounters) {
     printf("\nGlobal statistics of reconstructed tracks matched or not with the trigger:\n");
@@ -616,11 +873,17 @@ void AliAnalysisTaskMuonFakes::Terminate(Option_t *)
     fEventCounters->Print("event/trig");
   }
   
-  printf("\nREMINDER: results are relevent provided that you use the same recoParams as for the reconstruction\n");
+  if (fPairCounters) {
+    printf("\nGlobal statistics of reconstructed track pairs matched or not with the trigger:\n");
+    fPairCounters->Print("pair/trig");
+  }
+  
+  printf("\nREMINDER: results are relevent provided that you use the same recoParams as for the reconstruction\n\n");
 }
 
 //________________________________________________________________________
-Int_t AliAnalysisTaskMuonFakes::RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStore &trackRefStore)
+Int_t AliAnalysisTaskMuonFakes::RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStore &trackRefStore,
+                                                    TString &selected, TString &centrality)
 {
   /// loop over reconstructible TrackRef not associated with reconstructed track:
   /// for each of them, find and remove the most connected the fake track, if any,
@@ -691,16 +954,13 @@ Int_t AliAnalysisTaskMuonFakes::RemoveConnectedFakes(AliMUONVTrackStore &fakeTra
       // acceptance condition
       Double_t thetaTrackAbsEnd = TMath::ATan(esdTrack->GetRAtAbsorberEnd()/505.) * TMath::RadToDeg();
       Double_t eta = esdTrack->Eta();
-      TString acc = (thetaTrackAbsEnd >= 2. && thetaTrackAbsEnd <= 9. && eta >= -4. && eta <= -2.5) ? "acc:in" : "acc:out";
-      
-      // check physics selection
-      TString selected = (fInputHandler && fInputHandler->IsEventSelected()) ? "selected:yes" : "selected:no";
+      TString acc = (thetaTrackAbsEnd >= 2. && thetaTrackAbsEnd <= 10. && eta >= -4. && eta <= -2.5) ? "acc:in" : "acc:out";
       
       // fill histogram and counters
       ((TH1F*)fList->UncheckedAt(kFractionOfConnectedClusters))->Fill(fractionOfConnectedClusters);
-      fTrackCounters->Count(Form("track:connected/run:%d/%s/%s/%s", fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data()));
-      fFakeTrackCounters->Count(Form("track:connected/run:%d/file:%s/event:%d/%s/%s/%s", fCurrentRunNumber, fCurrentFileName.Data(),
-                                    esd->GetEventNumberInFile(), trig.Data(), selected.Data(), acc.Data()));
+      fTrackCounters->Count(Form("track:connected/run:%d/%s/%s/%s/%s", fCurrentRunNumber, trig.Data(), selected.Data(), acc.Data(), centrality.Data()));
+      fFakeTrackCounters->Count(Form("track:connected/run:%d/file:%s/event:%d/%s/%s/%s/%s", fCurrentRunNumber, fCurrentFileName.Data(),
+                                    esd->GetEventNumberInFile(), trig.Data(), selected.Data(), acc.Data(), centrality.Data()));
       
       // remove the most connected fake track
       fakeTrackStore.Remove(*connectedFake);
index 47a16b0efe238bbfac4a98cb74fda6ab5af659e0..51c58d7e9e5f61d1512099e18237e651fc352cc8 100644 (file)
@@ -33,6 +33,12 @@ public:
   /// Set the flag to match reconstructed and simulated tracks by using the MC labels or by position
   void UseMCLabels(Bool_t flag = kTRUE) { fUseLabel = flag; }
   
+  /// set the flag to fill histograms only with tracks matched with trigger or not
+  void MatchTrigger(Bool_t flag = kTRUE) { fMatchTrig = flag; }
+  
+  /// set the flag to fill histograms only with tracks passing the acceptance cuts (Rabs, eta)
+  void ApplyAccCut(Bool_t flag = kTRUE) { fApplyAccCut = flag; }
+  
   /// Set the ocdb path toward the reconstruction parameters
   void RecoParamLocation(const char* ocdbPath) { fRecoParamLocation = ocdbPath; }
   
@@ -47,7 +53,8 @@ private:
   AliAnalysisTaskMuonFakes& operator = (const AliAnalysisTaskMuonFakes& rhs);
   
   // look for fake tracks still connected to a reconstructible simulated track
-  Int_t RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStore &trackRefStore);
+  Int_t RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStore &trackRefStore,
+                            TString &selected, TString &centrality);
   
 private:
   
@@ -85,39 +92,74 @@ private:
     kChi2PerDofVsNChamberHitF, ///< normalized chi2 of fake tracks versus number of fired chambers
     
     // physics quantities
-    kP,    ///< momentum of tracks
-    kPM,   ///< momentum of matched tracks
-    kPF,   ///< momentum of fake tracks
-    kPt,   ///< transverse momentum of tracks
-    kPtM,  ///< transverse momentum of matched tracks
-    kPtF,  ///< transverse momentum of fake tracks
-    kEta,  ///< pseudo-rapidity of tracks
-    kEtaM, ///< pseudo-rapidity of matched tracks
-    kEtaF, ///< pseudo-rapidity of fake tracks
-    kPhi,  ///< phi angle of tracks
-    kPhiM, ///< phi angle of matched tracks
-    kPhiF, ///< phi angle of fake tracks
-    kDCA,  ///< DCA of tracks
-    kDCAM, ///< DCA of matched tracks
-    kDCAF, ///< DCA of fake tracks
+    kP,     ///< momentum of tracks
+    kPM,    ///< momentum of matched tracks
+    kPF,    ///< momentum of fake tracks
+    kPt,    ///< transverse momentum of tracks
+    kPtM,   ///< transverse momentum of matched tracks
+    kPtF,   ///< transverse momentum of fake tracks
+    kEta,   ///< pseudo-rapidity of tracks
+    kEtaM,  ///< pseudo-rapidity of matched tracks
+    kEtaF,  ///< pseudo-rapidity of fake tracks
+    kPhi,   ///< phi angle of tracks
+    kPhiM,  ///< phi angle of matched tracks
+    kPhiF,  ///< phi angle of fake tracks
+    kDCA,   ///< DCA of tracks
+    kDCAM,  ///< DCA of matched tracks
+    kDCAF,  ///< DCA of fake tracks
+    kRAbs,  ///< R of tracks at the end of the absorber
+    kRAbsM, ///< R of matched tracks at the end of the absorber
+    kRAbsF, ///< R of fake tracks at the end of the absorber
+  };
+  
+  enum histo2Index {
+    // physics quantities
+    k2Mass,   ///< invariant mass of the pair
+    k2MassM,  ///< invariant mass of matched-matched pairs
+    k2MassF1, ///< invariant mass of matched-fake pairs
+    k2MassF2, ///< invariant mass of fake-fake pairs
+    k2P,      ///< momentum of the pair
+    k2PM,     ///< momentum of matched-matched pairs
+    k2PF1,    ///< momentum of matched-fake pairs
+    k2PF2,    ///< momentum of fake-fake pairs
+    k2Pt,     ///< transverse momentum of pair
+    k2PtM,    ///< transverse momentum of matched-matched pairs
+    k2PtF1,   ///< transverse momentum of matched-fake pairs
+    k2PtF2,   ///< transverse momentum of fake-fake pairs
+    k2Y,      ///< rapidity of pair
+    k2YM,     ///< rapidity of matched-matched pairs
+    k2YF1,    ///< rapidity of matched-fake pairs
+    k2YF2,    ///< rapidity of fake-fake pairs
+    k2Eta,    ///< pseudo-rapidity of pair
+    k2EtaM,   ///< pseudo-rapidity of matched-matched pairs
+    k2EtaF1,  ///< pseudo-rapidity of matched-fake pairs
+    k2EtaF2,  ///< pseudo-rapidity of fake-fake pairs
+    k2Phi,    ///< phi angle of pair
+    k2PhiM,   ///< phi angle of matched-matched pairs
+    k2PhiF1,  ///< phi angle of matched-fake pairs
+    k2PhiF2,  ///< phi angle of fake-fake pairs
   };
   
-  TObjArray* fList;     //!< list of output histograms
+  TObjArray* fList;     //!< list of output histograms about single tracks
+  TObjArray* fList2;     //!< list of output histograms about track pairs
   TObjArray* fCanvases; //!< List of canvases summarizing the results
   
   AliCounterCollection* fTrackCounters;        //!< global counters of tracks
   AliCounterCollection* fFakeTrackCounters;    //!< detailled counters of fake tracks
   AliCounterCollection* fMatchedTrackCounters; //!< detailled counters of matched tracks
   AliCounterCollection* fEventCounters;        //!< counters of events
+  AliCounterCollection* fPairCounters;         //!< global counters of track pairs
   
   TString  fCurrentFileName;      //!< current input file name
   UInt_t   fRequestedStationMask; //!< sigma cut to associate clusters with TrackRefs
   Bool_t   fRequest2ChInSameSt45; //!< 2 fired chambers requested in the same station (4 or 5) or not
   Double_t fSigmaCut;             //!< mask of requested stations
   Bool_t   fUseLabel;             ///< match reconstructed and simulated tracks by using the MC labels or by position
+  Bool_t   fMatchTrig;            ///< fill histograms with tracks matched with trigger only
+  Bool_t   fApplyAccCut;          ///< fill histograms with tracks passing the acceptance cuts (Rabs, eta) only
   TString  fRecoParamLocation;    ///< ocdb path toward the reconstruction parameters
   
-  ClassDef(AliAnalysisTaskMuonFakes, 1); // fake muon analysis
+  ClassDef(AliAnalysisTaskMuonFakes, 2); // fake muon analysis
 };
 
 #endif