]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
New class for applying cut in ESDMuonTracks. Example of how to use ESDMuonTrackCuts...
authormartinez <martinez@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 17 Jun 2008 16:39:51 +0000 (16:39 +0000)
committermartinez <martinez@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 17 Jun 2008 16:39:51 +0000 (16:39 +0000)
PWG3/PWG3muonLinkDef.h
PWG3/libPWG3muon.pkg
PWG3/muon/AliAnalysisTaskESDMuonFilter.cxx
PWG3/muon/AliESDMuonTrackCuts.cxx [new file with mode: 0644]
PWG3/muon/AliESDMuonTrackCuts.h [new file with mode: 0644]
PWG3/muon/AnalysisTrainMuonCAF.C
PWG3/muon/AnalysisTrainMuonLocal.C
PWG3/muon/RunESDMuonFilter.C [deleted file]

index 4d05f326e97e42da0f68acdf7e238a738dd4b3f2..4723719561054304ffb1b5d55f970ed128a0dc6f 100644 (file)
@@ -11,6 +11,7 @@
 #pragma link C++ class AliAnalysisTaskLUT+;
 #pragma link C++ class AliAnalysisTaskESDMuonFilter+;
 #pragma link C++ class AliAnalysisTaskTrigChEff+;
+#pragma link C++ class AliESDMuonTrackCuts+;
 #endif
 
 
index 7c8db04ac514aff25f73ddbc4c68cd091ce54f1c..ef4aeb0d525587b900c1c3bbbf4844f99072c291 100644 (file)
@@ -6,7 +6,8 @@ SRCS:=   muon/AliAnalysisTaskESDMuonFilter.cxx \
                  muon/AliAnalysisTaskLUT.cxx  \
                  muon/AliAnalysisTaskTrigChEff.cxx \
                  muon/AliAODDimuon.cxx \
-                 muon/AliAODEventInfo.cxx 
+                 muon/AliAODEventInfo.cxx \
+                 muon/AliESDMuonTrackCuts.cxx 
 
 HDRS:= $(SRCS:.cxx=.h)
 
index c0b377b9ee34f65925bbfe456b7248b2be67b471..0b580f869776e6223096aba2fb89c5c85d4bbae8 100644 (file)
@@ -15,7 +15,7 @@
 
 // Add the muon tracks to the generic AOD track branch during the 
 // filtering of the ESD - R. Arnaldi 5/5/08
-
 #include <TChain.h>
 #include <TFile.h>
 
@@ -37,15 +37,15 @@ ClassImp(AliAnalysisTaskESDMuonFilter)
 ////////////////////////////////////////////////////////////////////////
 
 AliAnalysisTaskESDMuonFilter::AliAnalysisTaskESDMuonFilter():
-  AliAnalysisTaskSE(),
-  fTrackFilter(0x0)
+    AliAnalysisTaskSE(),
+    fTrackFilter(0x0)
 {
   // Default constructor
 }
 
 AliAnalysisTaskESDMuonFilter::AliAnalysisTaskESDMuonFilter(const char* name):
-  AliAnalysisTaskSE(name),
-  fTrackFilter(0x0)
+    AliAnalysisTaskSE(name),
+    fTrackFilter(0x0)
 {
   // Constructor
 }
@@ -53,13 +53,13 @@ AliAnalysisTaskESDMuonFilter::AliAnalysisTaskESDMuonFilter(const char* name):
 void AliAnalysisTaskESDMuonFilter::UserCreateOutputObjects()
 {
   // Create the output container
-  if (fTrackFilter) OutputTree()->GetUserInfo()->Add(fTrackFilter);
+    if (fTrackFilter) OutputTree()->GetUserInfo()->Add(fTrackFilter);
 }
 
 void AliAnalysisTaskESDMuonFilter::Init()
 {
   // Initialization
-  if (fDebug > 1) AliInfo("Init() \n");
+    if (fDebug > 1) AliInfo("Init() \n");
 }
 
 
@@ -68,7 +68,7 @@ void AliAnalysisTaskESDMuonFilter::UserExec(Option_t */*option*/)
   // Execute analysis for current event                                            
   Long64_t ientry = Entry();
   printf("Muon Filter: Analysing event # %5d\n", (Int_t) ientry);
-  
+
   ConvertESDtoAOD();
 }
 
@@ -76,87 +76,93 @@ void AliAnalysisTaskESDMuonFilter::ConvertESDtoAOD()
 {
   // ESD Muon Filter analysis task executed for each event
   AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
-  
   // Define arrays for muons
   Double_t pos[3];
   Double_t p[3];
   Double_t pid[10];
-  
-  // has to be changed once the muon pid is provided by the ESD
-  for (Int_t i = 0; i < 10; pid[i++] = 0.);
-  pid[AliAODTrack::kMuon]=1.;
-  
   AliAODHeader* header = AODEvent()->GetHeader();
   AliAODTrack *aodTrack = 0x0;
-  AliESDMuonTrack *esdMuTrack = 0x0;
-  
   // Access to the AOD container of tracks
   TClonesArray &tracks = *(AODEvent()->GetTracks());
   Int_t jTracks = header->GetRefMultiplicity();
   
   // Read primary vertex from AOD event 
-  AliAODVertex *primary = AODEvent()->GetPrimaryVertex();
+  AliAODVertex * primary = AODEvent()->GetPrimaryVertex();
   primary->Print();
-  
   // Loop on muon tracks to fill the AOD track branch
   Int_t nMuTracks = esd->GetNumberOfMuonTracks();
   printf("Number of Muon Tracks=%d\n",nMuTracks);
-  
-  // Update number of positive and negative tracks from AOD event (M.G.)
+   
+   // Update number of positive and negative tracks from AOD event (M.G.)
   Int_t nPosTracks = header->GetRefMultiplicityPos();
   Int_t nNegTracks = header->GetRefMultiplicityNeg();
-  
+   
   for (Int_t nMuTrack = 0; nMuTrack < nMuTracks; ++nMuTrack) {
-    esdMuTrack = esd->GetMuonTrack(nMuTrack);
-    
-    if (!esdMuTrack->ContainTrackerData()) continue;
-    
-    p[0] = esdMuTrack->Px(); 
-    p[1] = esdMuTrack->Py(); 
-    p[2] = esdMuTrack->Pz();
-    
-    pos[0] = esdMuTrack->GetNonBendingCoor(); 
-    pos[1] = esdMuTrack->GetBendingCoor(); 
-    pos[2] = esdMuTrack->GetZ();
-    
-    aodTrack = new(tracks[jTracks++]) AliAODTrack(esdMuTrack->GetUniqueID(), // ID
-                                                 0, // label
-                                                 p, // momentum
-                                                 kTRUE, // cartesian coordinate system
-                                                 pos, // position
-                                                 kFALSE, // isDCA
-                                                 0x0, // covariance matrix
-                                                 esdMuTrack->Charge(), // charge
-                                                 0, // ITSClusterMap
-                                                 pid, // pid
-                                                 primary, // primary vertex
-                                                 kFALSE, // used for vertex fit?
-                                                 kFALSE, // used for primary vertex fit?
-                                                 AliAODTrack::kPrimary); // track type
-    
-    aodTrack->SetXYAtDCA(esdMuTrack->GetNonBendingCoorAtDCA(), esdMuTrack->GetBendingCoorAtDCA());
-    aodTrack->SetPxPyPzAtDCA(esdMuTrack->PxAtDCA(), esdMuTrack->PyAtDCA(), esdMuTrack->PzAtDCA());
-    aodTrack->ConvertAliPIDtoAODPID();
-    aodTrack->SetChi2perNDF(esdMuTrack->GetChi2() / (2.*esdMuTrack->GetNHit() - 5.));
-    aodTrack->SetChi2MatchTrigger(esdMuTrack->GetChi2MatchTrigger());
-    aodTrack->SetHitsPatternInTrigCh(esdMuTrack->GetHitsPatternInTrigCh());
-    aodTrack->SetMuonClusterMap(esdMuTrack->GetMuonClusterMap());
-    aodTrack->SetMatchTrigger(esdMuTrack->GetMatchTrigger());
-    
-    primary->AddDaughter(aodTrack);
-    
-    if (esdMuTrack->Charge() > 0) nPosTracks++;
-    else nNegTracks++;
-  }
-  
+
+       AliESDMuonTrack *esdMuTrack = esd->GetMuonTrack(nMuTrack); 
+       
+       UInt_t selectInfo = 0;
+       // Track selection
+       if (fTrackFilter) {
+          selectInfo = fTrackFilter->IsSelected(esdMuTrack);
+          if (!selectInfo) {
+            continue;
+          }  
+       }
+       p[0] = esdMuTrack->Px(); 
+       p[1] = esdMuTrack->Py(); 
+       p[2] = esdMuTrack->Pz();
+       pos[0] = primary->GetX(); 
+       pos[1] = primary->GetY(); 
+       pos[2] = primary->GetZ();
+        
+       // has to be changed once the muon pid is provided by the ESD
+       for (Int_t i = 0; i < 10; pid[i++] = 0.); pid[AliAODTrack::kMuon]=1.;
+       primary->AddDaughter(aodTrack =
+         new(tracks[jTracks++]) AliAODTrack(0 , // no ID provided
+                                            0, // no label provided
+                                            p,
+                                            kTRUE,
+                                            pos,
+                                            kFALSE,
+                                            NULL, // no covariance matrix provided
+                                            esdMuTrack->Charge(), 
+                                            0, // no ITSClusterMap
+                                            pid,
+                                            primary,
+                                            kTRUE,  // check if this is right
+                                            kTRUE,  // not used for vertex fit
+                                            AliAODTrack::kPrimary,
+                                            selectInfo)
+        );
+        
+     if (esdMuTrack->Charge() > 0) nPosTracks++;
+     else nNegTracks++;
+
+     aodTrack->ConvertAliPIDtoAODPID();  
+     aodTrack->SetHitsPatternInTrigCh(esdMuTrack->GetHitsPatternInTrigCh());
+     Int_t track2Trigger = esdMuTrack->GetMatchTrigger();
+     aodTrack->SetMatchTrigger(track2Trigger);
+     if (track2Trigger) 
+       aodTrack->SetChi2MatchTrigger(esdMuTrack->GetChi2MatchTrigger());
+     else 
+       aodTrack->SetChi2MatchTrigger(0.);
+ }
+
   header->SetRefMultiplicity(jTracks); 
   header->SetRefMultiplicityPos(nPosTracks);
   header->SetRefMultiplicityNeg(nNegTracks);
+  
+  return;
 }
 
 void AliAnalysisTaskESDMuonFilter::Terminate(Option_t */*option*/)
 {
-  // Terminate analysis
-  //
-  if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");
+// Terminate analysis
+//
+    if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");
 }
diff --git a/PWG3/muon/AliESDMuonTrackCuts.cxx b/PWG3/muon/AliESDMuonTrackCuts.cxx
new file mode 100644 (file)
index 0000000..1dfac25
--- /dev/null
@@ -0,0 +1,580 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+#include "AliESDMuonTrackCuts.h"
+
+#include <AliESDMuonTrack.h>
+#include <AliESD.h>
+#include <AliESDEvent.h>
+#include <AliLog.h>
+
+#include <TTree.h>
+#include <TStyle.h>
+#include <TCanvas.h>
+#include <TDirectory.h>
+
+//____________________________________________________________________
+ClassImp(AliESDMuonTrackCuts)
+
+// Cut names
+const Char_t* AliESDMuonTrackCuts::fgkCutNames[kNCuts] = {
+ "p",
+ "p_{T}",
+ "p_{x}",
+ "p_{y}",
+ "p_{z}",
+ "y",
+ "eta"
+};
+
+//____________________________________________________________________
+AliESDMuonTrackCuts::AliESDMuonTrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title),
+  fPMin(0),
+  fPMax(0),
+  fPtMin(0),
+  fPtMax(0),
+  fPxMin(0),
+  fPxMax(0),
+  fPyMin(0),
+  fPyMax(0),
+  fPzMin(0),
+  fPzMax(0),
+  fEtaMin(0),
+  fEtaMax(0),
+  fRapMin(0),
+  fRapMax(0),
+  fHistogramsOn(0),
+  fhCutStatistics(0),         
+  fhCutCorrelation(0)
+{
+  //
+  // constructor
+  //
+  Init();
+
+  //##############################################################################
+  // setting default cuts
+  SetPRange();
+  SetPtRange();
+  SetPxRange();
+  SetPyRange();
+  SetPzRange();
+  SetEtaRange();
+  SetRapRange();
+
+  SetHistogramsOn();
+}
+
+//_____________________________________________________________________________
+AliESDMuonTrackCuts::AliESDMuonTrackCuts(const AliESDMuonTrackCuts &c) : AliAnalysisCuts(c),
+  fPMin(0),
+  fPMax(0),
+  fPtMin(0),
+  fPtMax(0),
+  fPxMin(0),
+  fPxMax(0),
+  fPyMin(0),
+  fPyMax(0),
+  fPzMin(0),
+  fPzMax(0),
+  fEtaMin(0),
+  fEtaMax(0),
+  fRapMin(0),
+  fRapMax(0),
+  fHistogramsOn(0),
+  fhCutStatistics(0),         
+  fhCutCorrelation(0)
+{
+  //
+  // copy constructor
+  //
+  ((AliESDMuonTrackCuts &) c).Copy(*this);
+}
+
+AliESDMuonTrackCuts::~AliESDMuonTrackCuts()
+{
+  //
+  // destructor
+  //
+  for (Int_t i=0; i<2; i++) {    
+    if (fhPt[i])
+      delete fhPt[i];
+    if (fhEta[i])
+      delete fhEta[i];
+  }
+
+  if (fhCutStatistics)
+    delete fhCutStatistics;             
+  if (fhCutCorrelation)
+    delete fhCutCorrelation;            
+}
+
+void AliESDMuonTrackCuts::Init()
+{
+  //
+  // sets everything to zero
+  //
+  fPMin = 0;
+  fPMax = 0;
+  fPtMin = 0;
+  fPtMax = 0;
+  fPxMin = 0;
+  fPxMax = 0;
+  fPyMin = 0;
+  fPyMax = 0;
+  fPzMin = 0;
+  fPzMax = 0;
+  fEtaMin = 0;
+  fEtaMax = 0;
+  fRapMin = 0;
+  fRapMax = 0;
+
+  fHistogramsOn = kFALSE;
+
+  for (Int_t i=0; i<2; ++i)
+  {
+    fhPt[i] = 0;
+    fhEta[i] = 0;
+  }
+  fhCutStatistics = 0;
+  fhCutCorrelation = 0;
+}
+
+//_____________________________________________________________________________
+AliESDMuonTrackCuts &AliESDMuonTrackCuts::operator=(const AliESDMuonTrackCuts &c)
+{
+  //
+  // Assignment operator
+  //
+  if (this != &c) ((AliESDMuonTrackCuts &) c).Copy(*this);
+  return *this;
+}
+
+//_____________________________________________________________________________
+void AliESDMuonTrackCuts::Copy(TObject &c) const
+{
+  //
+  // Copy function
+  //
+
+   AliESDMuonTrackCuts& target = (AliESDMuonTrackCuts &) c;
+  target.Init();
+// 
+  target.fPMin = fPMin;
+  target.fPMax = fPMax;
+  target.fPtMin = fPtMin;
+  target.fPtMax = fPtMax;
+  target.fPxMin = fPxMin;
+  target.fPxMax = fPxMax;
+  target.fPyMin = fPyMin;
+  target.fPyMax = fPyMax;
+  target.fPzMin = fPzMin;
+  target.fPzMax = fPzMax;
+  target.fEtaMin = fEtaMin;
+  target.fEtaMax = fEtaMax;
+  target.fRapMin = fRapMin;
+  target.fRapMax = fRapMax;
+
+  target.fHistogramsOn = fHistogramsOn;
+
+  for (Int_t i=0; i<2; ++i)
+  {     
+    if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone();
+    if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone();
+  }
+
+  if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
+  if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
+
+  TNamed::Copy(c);
+}
+
+//_____________________________________________________________________________
+Long64_t AliESDMuonTrackCuts::Merge(TCollection* list) {
+  // Merge a list of AliESDMuonTrackCuts objects with this (needed for PROOF)
+  // Returns the number of merged objects (including this)
+
+  if (!list)
+    return 0;
+  
+  if (list->IsEmpty())
+    return 1;
+
+  if (!fHistogramsOn)
+    return 0;
+
+  TIterator* iter = list->MakeIterator();
+  TObject* obj;
+
+  // collection of measured and generated histograms
+  Int_t count = 0;
+  while ((obj = iter->Next())) {
+
+    AliESDMuonTrackCuts* entry = dynamic_cast<AliESDMuonTrackCuts*>(obj);
+    if (entry == 0)
+      continue;
+
+    if (!entry->fHistogramsOn)
+      continue;
+    
+    for (Int_t i=0; i<2; i++) {
+      fhPt[i]->Add(entry->fhPt[i]); 
+      fhEta[i]->Add(entry->fhEta[i]); 
+    }      
+
+    fhCutStatistics->Add(entry->fhCutStatistics);        
+    fhCutCorrelation ->Add(entry->fhCutCorrelation);      
+
+    count++;
+  }
+
+  return count+1;
+}
+
+void AliESDMuonTrackCuts::EnableNeededBranches(TTree* tree)
+{
+  // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack
+
+  tree->SetBranchStatus("fTracks.fFlags", 1);
+  tree->SetBranchStatus("fTracks.fP*", 1);
+  tree->SetBranchStatus("fTracks.fR*", 1);  //detector response probability
+}
+
+//____________________________________________________________________
+Bool_t
+AliESDMuonTrackCuts::AcceptTrack(AliESDMuonTrack* esdMuTrack) {
+  // 
+  // figure out if the tracks survives all the track cuts defined
+  //
+  // the different kinematic values are first
+  // retrieved from the track. then it is found out what cuts the
+  // track did not survive and finally the cuts are imposed.
+
+  // getting the kinematic variables of the track
+  // (assuming the mass is known)
+  Double_t p[3];
+  esdMuTrack->PxPyPz(p);
+  Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
+  Float_t pt       = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
+  Float_t energy   = TMath::Sqrt(TMath::Power(esdMuTrack->M(),2) + TMath::Power(momentum,2));
+
+
+  //y-eta related calculations
+  Float_t eta = -100.;
+  Float_t y   = -100.;
+  if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
+    eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
+  if((energy != TMath::Abs(p[2]))&&(momentum != 0))
+    y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
+     
+  //########################################################################
+  // cut the track?
+  
+  Bool_t cuts[kNCuts];
+  for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
+  
+  // track kinematics cut
+  if((momentum < fPMin) || (momentum > fPMax)) 
+    cuts[0]=kTRUE;
+  if((pt < fPtMin) || (pt > fPtMax)) 
+    cuts[1] = kTRUE;
+  if((p[0] < fPxMin) || (p[0] > fPxMax)) 
+    cuts[2] = kTRUE;
+  if((p[1] < fPyMin) || (p[1] > fPyMax)) 
+    cuts[3] = kTRUE;
+  if((p[2] < fPzMin) || (p[2] > fPzMax))
+    cuts[4] = kTRUE;
+  if((eta < fEtaMin) || (eta > fEtaMax)) 
+    cuts[5] = kTRUE;
+  if((y < fRapMin) || (y > fRapMax)) 
+    cuts[6] = kTRUE;
+
+  Bool_t cut=kFALSE;
+  for (Int_t i=0; i<kNCuts; i++) 
+    if (cuts[i]) cut = kTRUE;
+  
+  //########################################################################
+  // filling histograms
+  if (fHistogramsOn) {
+    fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
+    
+    if (cut)
+      fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
+    
+    for (Int_t i=0; i<kNCuts; i++) {
+      if (cuts[i])
+       fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
+      
+      for (Int_t j=i; j<kNCuts; j++) {
+       if (cuts[i] && cuts[j]) {
+         Float_t xC = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
+         Float_t yC = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
+         fhCutCorrelation->Fill(xC, yC);
+       }
+      }
+    }
+    
+    fhPt[0]->Fill(pt);
+    fhEta[0]->Fill(eta);
+
+  }
+
+  //########################################################################
+  // cut the track!
+  if (cut) return kFALSE;
+
+  //########################################################################
+  // filling histograms after cut
+  if (fHistogramsOn) {
+// 
+    fhPt[1]->Fill(pt);
+    fhEta[1]->Fill(eta);
+    
+  }
+
+  return kTRUE;
+}
+
+//____________________________________________________________________
+TObjArray* AliESDMuonTrackCuts::GetAcceptedTracks(AliESD* esd)
+{
+  //
+  // returns an array of all tracks that pass the cuts
+  //
+
+  TObjArray* acceptedTracks = new TObjArray();
+
+  // loop over esd tracks
+  for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
+    AliESDMuonTrack* track = esd->GetMuonTrack(iTrack);
+
+    if (AcceptTrack(track))
+      acceptedTracks->Add(track);
+  }
+
+  return acceptedTracks;
+}
+
+//____________________________________________________________________
+Int_t AliESDMuonTrackCuts::CountAcceptedTracks(AliESD* esd)
+{
+  //
+  // returns an the number of tracks that pass the cuts
+  //
+
+  Int_t count = 0;
+
+  // loop over esd tracks
+  for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
+    AliESDMuonTrack* track = esd->GetMuonTrack(iTrack);
+
+    if (AcceptTrack(track))
+      count++;
+  }
+
+  return count;
+}
+
+//____________________________________________________________________
+TObjArray* AliESDMuonTrackCuts::GetAcceptedTracks(AliESDEvent* esd)
+{
+  //
+  // returns an array of all tracks that pass the cuts
+  //
+
+  TObjArray* acceptedTracks = new TObjArray();
+
+  // loop over esd tracks
+  for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
+    AliESDMuonTrack* track = esd->GetMuonTrack(iTrack);
+
+    if (AcceptTrack(track))
+      acceptedTracks->Add(track);
+  }
+
+  return acceptedTracks;
+}
+
+//____________________________________________________________________
+Int_t AliESDMuonTrackCuts::CountAcceptedTracks(AliESDEvent* esd)
+{
+  //
+  // returns an the number of tracks that pass the cuts
+  //
+
+  Int_t count = 0;
+
+  // loop over esd tracks
+  for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
+    AliESDMuonTrack* track = esd->GetMuonTrack(iTrack);
+
+    if (AcceptTrack(track))
+      count++;
+  }
+
+  return count;
+}
+
+//____________________________________________________________________
+ void AliESDMuonTrackCuts::DefineHistograms(Int_t color) {
+   // 
+   // diagnostics histograms are defined
+   // 
+
+   fHistogramsOn=kTRUE;
+   
+   Bool_t oldStatus = TH1::AddDirectoryStatus();
+   TH1::AddDirectory(kFALSE);
+   
+   //###################################################################################
+   // defining histograms
+
+   fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
+
+   fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
+   fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
+
+   fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
+  
+   for (Int_t i=0; i<kNCuts; i++) {
+     fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
+     fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
+     fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
+   } 
+
+  fhCutStatistics  ->SetLineColor(color);
+  fhCutCorrelation ->SetLineColor(color);
+  fhCutStatistics  ->SetLineWidth(2);
+  fhCutCorrelation ->SetLineWidth(2);
+
+  Char_t str[256];
+  for (Int_t i=0; i<2; i++) {
+    if (i==0) sprintf(str," ");
+    else sprintf(str,"_cut");
+
+    fhPt[i]                  = new TH1F(Form("pt%s",str)     ,"p_{T} distribution;p_{T} (GeV/c)",500,0.0,100.0);
+    fhEta[i]                 = new TH1F(Form("eta%s",str)     ,"#eta distribution;#eta",40,-2.0,2.0);
+  }
+
+  TH1::AddDirectory(oldStatus);
+}
+
+//____________________________________________________________________
+Bool_t AliESDMuonTrackCuts::LoadHistograms(const Char_t* dir)
+{
+  //
+  // loads the histograms from a file
+  // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
+  //
+
+  if (!dir)
+    dir = GetName();
+
+  if (!gDirectory->cd(dir))
+    return kFALSE;
+
+  fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
+  fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
+
+  Char_t str[5];
+  for (Int_t i=0; i<2; i++) {
+    if (i==0)
+    {
+      gDirectory->cd("before_cuts");
+      str[0] = 0;
+    }
+    else
+    {
+      gDirectory->cd("after_cuts");
+      sprintf(str,"_cut");
+    }
+
+    fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("pt%s",str)));
+    fhEta[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("eta%s",str)));
+
+    gDirectory->cd("../");
+  }
+
+  gDirectory->cd("..");
+
+  return kTRUE;
+}
+
+//____________________________________________________________________
+void AliESDMuonTrackCuts::SaveHistograms(const Char_t* dir) {
+  //
+  // saves the histograms in a directory (dir)
+  //
+
+  if (!fHistogramsOn) {
+    AliDebug(0, "Histograms not on - cannot save histograms!!!");
+    return;
+  }
+
+  if (!dir)
+    dir = GetName();
+
+  gDirectory->mkdir(dir);
+  gDirectory->cd(dir);
+
+  gDirectory->mkdir("before_cuts");
+  gDirectory->mkdir("after_cuts");
+
+  fhCutStatistics->Write();
+  fhCutCorrelation->Write();
+
+  for (Int_t i=0; i<2; i++) {
+    if (i==0)
+      gDirectory->cd("before_cuts");
+    else
+      gDirectory->cd("after_cuts");
+
+    fhPt[i]                  ->Write();
+    fhEta[i]                 ->Write();
+    
+    gDirectory->cd("../");
+  }
+
+  gDirectory->cd("../");
+}
+
+//____________________________________________________________________
+void AliESDMuonTrackCuts::DrawHistograms()
+{
+  gStyle->SetPalette(1);
+  gStyle->SetFrameFillColor(10);
+  gStyle->SetCanvasColor(10);
+  
+  TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Cut Results", 800, 500);
+  canvas1->Divide(2, 1);
+
+  canvas1->cd(1);
+  fhCutStatistics->SetStats(kFALSE);
+  fhCutStatistics->LabelsOption("v");
+  gPad->SetBottomMargin(0.3);
+  fhCutStatistics->Draw();
+
+  canvas1->cd(2);
+  fhCutCorrelation->SetStats(kFALSE);
+  fhCutCorrelation->LabelsOption("v");
+  gPad->SetBottomMargin(0.3);
+  gPad->SetLeftMargin(0.3);
+  fhCutCorrelation->Draw("COLZ");
+  canvas1->Update();
+  canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));
+  
+}
+
diff --git a/PWG3/muon/AliESDMuonTrackCuts.h b/PWG3/muon/AliESDMuonTrackCuts.h
new file mode 100644 (file)
index 0000000..a19196a
--- /dev/null
@@ -0,0 +1,97 @@
+//
+//  Class for handling of ESD Muon track cuts 
+//  (based on ANALYSIS/AliESDtrackCuts).
+//
+//  The class manages some kinematic cuts. Two methods
+//  can be used to figure out if an ESD Muon track survives the cuts:
+//  AcceptTrack which takes a single AliESDMuonTrack as argument and
+//  returns kTRUE/kFALSE or GetAcceptedTracks which takes an AliESD
+//  object and returns an TObjArray (of AliESDMuonTracks) with the tracks
+//  in the ESD that survived the cuts.
+//
+
+#ifndef ALIESDMUONTRACKCUTS_H
+#define ALIESDMUONTRACKCUTS_H
+
+#include <TF1.h>
+#include <TH2.h>
+#include "AliAnalysisCuts.h"
+
+class AliESD;
+class AliESDEvent;
+class AliESDMuonTrack;
+class AliLog;
+class TTree;
+
+class AliESDMuonTrackCuts : public AliAnalysisCuts
+{
+public:
+  AliESDMuonTrackCuts(const Char_t* name = "AliESDMuonTrackCuts", const Char_t* title = "");
+  virtual ~AliESDMuonTrackCuts();
+  Bool_t IsSelected(TObject* obj)
+       {return AcceptTrack((AliESDMuonTrack*)obj);}
+  Bool_t IsSelected(TList* /*list*/) {return kTRUE;}
+  Bool_t AcceptTrack(AliESDMuonTrack* esdMuTrack);
+  TObjArray* GetAcceptedTracks(AliESD* esd);
+  Int_t CountAcceptedTracks(AliESD* esd);
+  TObjArray* GetAcceptedTracks(AliESDEvent* esd);
+  Int_t CountAcceptedTracks(AliESDEvent* esd);
+
+  virtual Long64_t Merge(TCollection* list);
+  virtual void Copy(TObject &c) const;
+  AliESDMuonTrackCuts(const AliESDMuonTrackCuts& pd);  // Copy Constructor
+  AliESDMuonTrackCuts &operator=(const AliESDMuonTrackCuts &c);
+
+  //######################################################
+
+  // track kinematic cut setters
+  void SetPRange(Float_t r1=0, Float_t r2=1e10)       {fPMin=r1;   fPMax=r2;}
+  void SetPtRange(Float_t r1=0, Float_t r2=1e10)      {fPtMin=r1;  fPtMax=r2;}
+  void SetPxRange(Float_t r1=-1e10, Float_t r2=1e10)  {fPxMin=r1;  fPxMax=r2;}
+  void SetPyRange(Float_t r1=-1e10, Float_t r2=1e10)  {fPyMin=r1;  fPyMax=r2;}
+  void SetPzRange(Float_t r1=-1e10, Float_t r2=1e10)  {fPzMin=r1;  fPzMax=r2;}
+  void SetEtaRange(Float_t r1=-1e10, Float_t r2=1e10) {fEtaMin=r1; fEtaMax=r2;}
+  void SetRapRange(Float_t r1=-1e10, Float_t r2=1e10) {fRapMin=r1; fRapMax=r2;}
+
+  //######################################################
+  
+  void SetHistogramsOn(Bool_t b=kFALSE) {fHistogramsOn = b;}
+  void DefineHistograms(Int_t color=1);
+  virtual Bool_t LoadHistograms(const Char_t* dir = 0);
+  void SaveHistograms(const Char_t* dir = 0);
+  void DrawHistograms();
+
+  static void EnableNeededBranches(TTree* tree);
+
+protected:
+  void Init(); // sets everything to 0
+
+  enum { kNCuts = 7 };
+
+  //######################################################
+  static const Char_t* fgkCutNames[kNCuts]; //! names of cuts (for internal use)
+
+  // esd kinematics cuts
+  Float_t fPMin,   fPMax;             // definition of the range of the P
+  Float_t fPtMin,  fPtMax;            // definition of the range of the Pt
+  Float_t fPxMin,  fPxMax;            // definition of the range of the Px
+  Float_t fPyMin,  fPyMax;            // definition of the range of the Py
+  Float_t fPzMin,  fPzMax;            // definition of the range of the Pz
+  Float_t fEtaMin, fEtaMax;           // definition of the range of the eta
+  Float_t fRapMin, fRapMax;           // definition of the range of the y
+
+  //######################################################
+  // diagnostics histograms
+  Bool_t fHistogramsOn;               // histograms on/off
+  
+  TH1F* fhPt[2];                      //-> pt of esd tracks
+  TH1F* fhEta[2];                     //-> eta of esd tracks
+
+  TH1F*  fhCutStatistics;             //-> statistics of what cuts the tracks did not survive
+  TH2F*  fhCutCorrelation;            //-> 2d statistics plot
+
+  ClassDef(AliESDMuonTrackCuts, 1)
+};
+
+
+#endif
index 01db2fde2f6ee8af59883028b4bbb1498ec05c04..8e27c339993637ee54479e2d6f4dbc755419a434 100644 (file)
-void AnalysisTrainMuonCAF(char* fileout = "AliAOD.root", char *datasetname = "myDataSet", Int_t nev=1234567890)\r
-{\r
-// Macro to produce a generic AOD starting from an ESD file. \r
-// The AOD is filled with two tasks: \r
-// 1- with the first one (AliAnalysisTaskESDfilter), \r
-//    all the branches of the AOD are filled apart from the muons. \r
-// 2- with the second task (AliAnalysisTaskESDMuonFilter) \r
-//    muons tracks are added to the tracks branch \r
-// This macro works on the CAF\r
-// R. Arnaldi 4/5/08\r
-\r
-  gSystem->Load("libTree.so");\r
-  gSystem->Load("libGeom.so");\r
-  gSystem->Load("libVMC.so");\r
-  gSystem->Load("libPhysics.so");\r
-    \r
-  // Reset user processes if CAF if not responding anymore\r
-  // TProof::Reset("lxb6046");\r
-\r
-  // Connect to proof\r
-  TProof::Open("lxb6046"); // may be username@lxb6046 if user not the same as on local\r
-\r
-  // Clear packages if changing ROOT version on CAF or local\r
-  // gProof->ClearPackages();\r
-  // Enable proof debugging if needed\r
-  // gProof->SetLogLevel(5);\r
-\r
-  // Common packages\r
-  gProof->UploadPackage("STEERBase.par");\r
-  gProof->EnablePackage("STEERBase");\r
-  gProof->UploadPackage("ESD.par");\r
-  gProof->EnablePackage("ESD");\r
-  gProof->UploadPackage("AOD.par");\r
-  gProof->EnablePackage("AOD");\r
-  gProof->UploadPackage("ANALYSIS.par");\r
-  gProof->EnablePackage("ANALYSIS");\r
-  gProof->UploadPackage("ANALYSISalice.par");\r
-  gProof->EnablePackage("ANALYSISalice");\r
-  // Analysis-specific\r
-  // --- Enable the PWG3base Package\r
-  gProof->UploadPackage("PWG3muon.par");\r
-  gProof->EnablePackage("PWG3muon");\r
-\r
-  // Chain from files staged on CAF\r
-  // gROOT->LoadMacro("CreateESDChain.C");\r
-  // TChain* chain = CreateESDChain("ESD1503X_v1.txt",3);\r
-  // TChain* chain = CreateESDChain("ESD82XX_30Kshort.txt", 10);\r
-  \r
-  // Chain from datasets\r
-  gROOT->LoadMacro("CreateChainFromDataSet.C");\r
-  ds = gProof->GetDataSet(datasetname)->GetStagedSubset();\r
-  chain = CreateChainFromDataSet(ds, "esdTree");   \r
-  \r
-  // Make the analysis manager\r
-  AliAnalysisManager *mgr  = new AliAnalysisManager("Analysis Train", "Analysis train");\r
-  \r
-  // ESD input handler\r
-  AliESDInputHandler *esdHandler = new AliESDInputHandler();\r
-  esdHandler->SetInactiveBranches("FMD CaloCluster");\r
-  // AOD output handler\r
-  AliAODHandler* aodHandler   = new AliAODHandler();\r
-  //aodHandler->SetOutputFileName(fileout);\r
-  aodHandler->SetOutputFileName("AOD.root");\r
-\r
-  mgr->SetInputEventHandler(esdHandler);\r
-  mgr->SetOutputEventHandler(aodHandler);\r
-  \r
-  // Set of cuts plugged into the ESD filter\r
-  // \r
-  // standard\r
-  AliESDtrackCuts* esdTrackCutsL = new AliESDtrackCuts("AliESDtrackCuts", "Loose");\r
-  esdTrackCutsL->SetMinNClustersTPC(50);\r
-  esdTrackCutsL->SetMaxChi2PerClusterTPC(3.5);\r
-  esdTrackCutsL->SetMaxCovDiagonalElements(2,2,0.5,0.5,2);\r
-  esdTrackCutsL->SetRequireTPCRefit(kTRUE);\r
-  esdTrackCutsL->SetMinNsigmaToVertex(3);\r
-  esdTrackCutsL->SetRequireSigmaToVertex(kTRUE);\r
-  esdTrackCutsL->SetAcceptKingDaughters(kFALSE);\r
-  //\r
-  // hard\r
-  AliESDtrackCuts* esdTrackCutsH = new AliESDtrackCuts("AliESDtrackCuts", "Hard");\r
-  esdTrackCutsH->SetMinNClustersTPC(100);\r
-  esdTrackCutsH->SetMaxChi2PerClusterTPC(2.0);\r
-  esdTrackCutsH->SetMaxCovDiagonalElements(2,2,0.5,0.5,2);\r
-  esdTrackCutsH->SetRequireTPCRefit(kTRUE);\r
-  esdTrackCutsH->SetMinNsigmaToVertex(2);\r
-  esdTrackCutsH->SetRequireSigmaToVertex(kTRUE);\r
-  esdTrackCutsH->SetAcceptKingDaughters(kFALSE);\r
-  //\r
-  AliAnalysisFilter* trackFilter = new AliAnalysisFilter("trackFilter");\r
-  trackFilter->AddCuts(esdTrackCutsL);\r
-  trackFilter->AddCuts(esdTrackCutsH);\r
-\r
-  // ESD filter task putting standard info to output generic AOD \r
-  AliAnalysisTaskESDfilter *esdfilter = new AliAnalysisTaskESDfilter("ESD Filter");\r
-  //esdfilter->SetTrackFilter(trackFilter);\r
-  esdfilter->SetDebugLevel(10);\r
-  mgr->AddTask(esdfilter);\r
-  \r
-  // ESD filter task putting muon info to output generic AOD \r
-  AliAnalysisTaskESDMuonFilter *esdmuonfilter = new AliAnalysisTaskESDMuonFilter("ESD Muon Filter");\r
-  mgr->AddTask(esdmuonfilter);\r
-\r
-  // Containers for input/output\r
-  AliAnalysisDataContainer *cin_esd = mgr->CreateContainer("cESD",TChain::Class(), AliAnalysisManager::kInputContainer);\r
-  // Output AOD container. \r
-  AliAnalysisDataContainer *cout_aod = mgr->CreateContainer("cAOD", TTree::Class(), AliAnalysisManager::kOutputContainer, "default");\r
-                                                           \r
-  // Connect containers to tasks slots\r
-  mgr->ConnectInput  (esdfilter,  0, cin_esd  );\r
-  mgr->ConnectOutput (esdfilter,  0, cout_aod );\r
-\r
-  mgr->ConnectInput  (esdmuonfilter,  0, cin_esd);\r
-  mgr->ConnectOutput (esdmuonfilter,  0, cout_aod );\r
-\r
-  //\r
-  // Run the analysis\r
-  //   \r
-  if (mgr->InitAnalysis()) {\r
-      mgr->PrintStatus();\r
-      mgr->StartAnalysis("proof",chain,nev);\r
-  }   \r
-}\r
-\r
+void AnalysisTrainMuonCAF(char* fileout = "AliAOD.root", char *datasetname = "myDataSet", Int_t nev=1234567890)
+{
+// Macro to produce a generic AOD starting from an ESD file. 
+// The AOD is filled with two tasks: 
+// 1- with the first one (AliAnalysisTaskESDfilter), 
+//    all the branches of the AOD are filled apart from the muons. 
+// 2- with the second task (AliAnalysisTaskESDMuonFilter) 
+//    muons tracks are added to the tracks branch 
+// This macro works on the CAF
+// R. Arnaldi 4/5/08
+
+  gSystem->Load("libTree.so");
+  gSystem->Load("libGeom.so");
+  gSystem->Load("libVMC.so");
+  gSystem->Load("libPhysics.so");
+    
+  // Reset user processes if CAF if not responding anymore
+  // TProof::Reset("lxb6046");
+
+  // Connect to proof
+  TProof::Open("lxb6046"); // may be username@lxb6046 if user not the same as on local
+
+  // Clear packages if changing ROOT version on CAF or local
+  // gProof->ClearPackages();
+  // Enable proof debugging if needed
+  // gProof->SetLogLevel(5);
+
+  // Common packages
+  gProof->UploadPackage("STEERBase.par");
+  gProof->EnablePackage("STEERBase");
+  gProof->UploadPackage("ESD.par");
+  gProof->EnablePackage("ESD");
+  gProof->UploadPackage("AOD.par");
+  gProof->EnablePackage("AOD");
+  gProof->UploadPackage("ANALYSIS.par");
+  gProof->EnablePackage("ANALYSIS");
+  gProof->UploadPackage("ANALYSISalice.par");
+  gProof->EnablePackage("ANALYSISalice");
+  // Analysis-specific
+  // --- Enable the PWG3base Package
+  gProof->UploadPackage("PWG3muon.par");
+  gProof->EnablePackage("PWG3muon");
+
+  // Chain from files staged on CAF
+  // gROOT->LoadMacro("CreateESDChain.C");
+  // TChain* chain = CreateESDChain("ESD1503X_v1.txt",3);
+  // TChain* chain = CreateESDChain("ESD82XX_30Kshort.txt", 10);
+  
+  // Chain from datasets
+  gROOT->LoadMacro("CreateChainFromDataSet.C");
+  ds = gProof->GetDataSet(datasetname)->GetStagedSubset();
+  chain = CreateChainFromDataSet(ds, "esdTree");   
+  
+  // Make the analysis manager
+  AliAnalysisManager *mgr  = new AliAnalysisManager("Analysis Train", "Analysis train");
+  
+  // ESD input handler
+  AliESDInputHandler *esdHandler = new AliESDInputHandler();
+  esdHandler->SetInactiveBranches("FMD CaloCluster");
+  
+  // AOD output handler
+  AliAODHandler* aodHandler   = new AliAODHandler();
+  aodHandler->SetOutputFileName(fileout);
+  //aodHandler->SetOutputFileName("AOD.root");
+
+  mgr->SetInputEventHandler(esdHandler);
+  mgr->SetOutputEventHandler(aodHandler);
+  
+  // Set of cuts plugged into the ESD filter
+  // 
+  // standard
+  AliESDtrackCuts* esdTrackCutsL = new AliESDtrackCuts("AliESDtrackCuts", "Loose");
+  esdTrackCutsL->SetMinNClustersTPC(50);
+  esdTrackCutsL->SetMaxChi2PerClusterTPC(3.5);
+  esdTrackCutsL->SetMaxCovDiagonalElements(2,2,0.5,0.5,2);
+  esdTrackCutsL->SetRequireTPCRefit(kTRUE);
+  esdTrackCutsL->SetMinNsigmaToVertex(3);
+  esdTrackCutsL->SetRequireSigmaToVertex(kTRUE);
+  esdTrackCutsL->SetAcceptKingDaughters(kFALSE);
+  //
+  // hard cuts
+  AliESDtrackCuts* esdTrackCutsH = new AliESDtrackCuts("AliESDtrackCuts", "Hard");
+  esdTrackCutsH->SetMinNClustersTPC(100);
+  esdTrackCutsH->SetMaxChi2PerClusterTPC(2.0);
+  esdTrackCutsH->SetMaxCovDiagonalElements(2,2,0.5,0.5,2);
+  esdTrackCutsH->SetRequireTPCRefit(kTRUE);
+  esdTrackCutsH->SetMinNsigmaToVertex(2);
+  esdTrackCutsH->SetRequireSigmaToVertex(kTRUE);
+  esdTrackCutsH->SetAcceptKingDaughters(kFALSE);
+  esdTrackCutsH->SetPRange(0.,2.);
+  //
+  //  muon cuts
+  AliESDMuonTrackCuts* esdMuonTrackCuts = new AliESDMuonTrackCuts("AliESDMuonTrackCuts", "test");
+  esdMuonTrackCuts->SetPRange(0.,20.);
+  //esdMuonTrackCuts->SetPtRange(0.,0.5);   // example of kinematic cuts that can be applied
+  
+  // track filter (to reject tracks not surviving the cuts - refers to all particles apart from muons)
+  AliAnalysisFilter* trackFilter = new AliAnalysisFilter("trackFilter");
+  trackFilter->AddCuts(esdTrackCutsH);
+  
+  // muon track filter  (to reject muon tracks not surviving the cuts)
+  AliAnalysisFilter* trackMuonFilter = new AliAnalysisFilter("trackMuonFilter");
+  trackMuonFilter->AddCuts(esdMuonTrackCuts);
+
+  // ESD filter task putting standard info to output generic AOD 
+  AliAnalysisTaskESDfilter *esdfilter = new AliAnalysisTaskESDfilter("ESD Filter");
+  //esdfilter->SetTrackFilter(trackFilter);
+  esdfilter->SetDebugLevel(10);
+  mgr->AddTask(esdfilter);
+  
+  // ESD filter task putting muon info to output generic AOD 
+  AliAnalysisTaskESDMuonFilter *esdmuonfilter = new AliAnalysisTaskESDMuonFilter("ESD Muon Filter");
+  esdmuonfilter->SetTrackFilter(trackMuonFilter);
+  mgr->AddTask(esdmuonfilter);
+
+  // Containers for input/output
+  AliAnalysisDataContainer *cin_esd = mgr->CreateContainer("cESD",TChain::Class(), AliAnalysisManager::kInputContainer);
+  // Output AOD container. 
+  AliAnalysisDataContainer *cout_aod = mgr->CreateContainer("cAOD", TTree::Class(), AliAnalysisManager::kOutputContainer, "default");
+                                                           
+  // Connect containers to tasks slots
+  mgr->ConnectInput  (esdfilter,  0, cin_esd  );
+  mgr->ConnectOutput (esdfilter,  0, cout_aod );
+
+  mgr->ConnectInput  (esdmuonfilter,  0, cin_esd);
+  mgr->ConnectOutput (esdmuonfilter,  0, cout_aod );
+
+  //
+  // Run the analysis
+  //   
+  if (mgr->InitAnalysis()) {
+      mgr->PrintStatus();
+      mgr->StartAnalysis("proof",chain,nev);
+  }   
+}
+
index 72b9f28a04e9905155889a0696d59850219574b6..402dae3c33baa476bec182f08e12b5d20cc34ad8 100644 (file)
-void AnalysisTrainMuonLocal(char* filein = "AliESDs.root", char* fileout = "AliAOD.root")\r
-\r
-// Macro to produce a generic AOD starting from an ESD file. \r
-// The AOD is filled with two tasks: \r
-// 1- with the first one (AliAnalysisTaskESDfilter), \r
-//    all the branches of the AOD are filled apart from the muons. \r
-// 2- with the second task (AliAnalysisTaskESDMuonFilter) \r
-//    muons tracks are added to the tracks branch \r
-// This macro works locally\r
-// R. Arnaldi 5/5/08\r
-\r
-{\r
-    gSystem->Load("libTree.so");\r
-    gSystem->Load("libGeom.so");\r
-    gSystem->Load("libVMC.so");\r
-    gSystem->Load("libPhysics.so");\r
-\r
-    // If analysis is .par based:\r
-\r
-    // Common packages\r
-    SetupPar("STEERBase");\r
-    SetupPar("ESD");\r
-    SetupPar("AOD");\r
-    SetupPar("ANALYSIS");\r
-    SetupPar("ANALYSISalice");\r
-    // Analysis-specific packages\r
-    SetupPar("PWG3muon");      \r
-  \r
-    // Input ESD file\r
-    TChain* chain = new TChain("esdTree");  \r
-    chain->Add(filein);\r
-   \r
-    // Define the analysis manager\r
-    AliAnalysisManager *mgr  = new AliAnalysisManager("Analysis Train", "Analysis train");\r
-    \r
-    // ESD input handler\r
-    AliESDInputHandler *esdHandler = new AliESDInputHandler();\r
-    esdHandler->SetInactiveBranches("FMD CaloCluster");\r
-    \r
-    // AOD output handler\r
-    AliAODHandler* aodHandler   = new AliAODHandler();\r
-    aodHandler->SetOutputFileName(fileout);\r
-\r
-    mgr->SetInputEventHandler(esdHandler);\r
-    mgr->SetOutputEventHandler(aodHandler);\r
-    \r
-    // Set of cuts for the ESD filter\r
-    // \r
-    // standard cuts\r
-    AliESDtrackCuts* esdTrackCutsL = new AliESDtrackCuts("AliESDtrackCuts", "Loose");\r
-    esdTrackCutsL->SetMinNClustersTPC(50);\r
-    esdTrackCutsL->SetMaxChi2PerClusterTPC(3.5);\r
-    esdTrackCutsL->SetMaxCovDiagonalElements(2,2,0.5,0.5,2);\r
-    esdTrackCutsL->SetRequireTPCRefit(kTRUE);\r
-    esdTrackCutsL->SetMinNsigmaToVertex(3);\r
-    esdTrackCutsL->SetRequireSigmaToVertex(kTRUE);\r
-    esdTrackCutsL->SetAcceptKingDaughters(kFALSE);\r
-    //\r
-    // hard cuts\r
-    AliESDtrackCuts* esdTrackCutsH = new AliESDtrackCuts("AliESDtrackCuts", "Hard");\r
-    esdTrackCutsH->SetMinNClustersTPC(100);\r
-    esdTrackCutsH->SetMaxChi2PerClusterTPC(2.0);\r
-    esdTrackCutsH->SetMaxCovDiagonalElements(2,2,0.5,0.5,2);\r
-    esdTrackCutsH->SetRequireTPCRefit(kTRUE);\r
-    esdTrackCutsH->SetMinNsigmaToVertex(2);\r
-    esdTrackCutsH->SetRequireSigmaToVertex(kTRUE);\r
-    esdTrackCutsH->SetAcceptKingDaughters(kFALSE);\r
-    //\r
-    AliAnalysisFilter* trackFilter = new AliAnalysisFilter("trackFilter");\r
-    trackFilter->AddCuts(esdTrackCutsL);\r
-    trackFilter->AddCuts(esdTrackCutsH);\r
-\r
-    // ESD filter task putting standard info in the output generic AOD \r
-    AliAnalysisTaskESDfilter *esdfilter = new AliAnalysisTaskESDfilter("ESD Filter");\r
-    //esdfilter->SetTrackFilter(trackFilter);\r
-    esdfilter->SetDebugLevel(10);\r
-    mgr->AddTask(esdfilter);\r
-    \r
-    // ESD filter task putting muon info in the output generic AOD \r
-    AliAnalysisTaskESDMuonFilter *esdmuonfilter = new AliAnalysisTaskESDMuonFilter("ESD Muon Filter");\r
-    mgr->AddTask(esdmuonfilter);\r
-\r
-    // Containers for input/output\r
-    AliAnalysisDataContainer *cin_esd = mgr->CreateContainer("cESD",TChain::Class(), \r
-                                                            AliAnalysisManager::kInputContainer);\r
-    // Output AOD container. \r
-    AliAnalysisDataContainer *cout_aod = mgr->CreateContainer("cAOD", TTree::Class(),\r
-                                                             AliAnalysisManager::kOutputContainer, "default");\r
-                                                             \r
-    // Connect containers to tasks slots\r
-    mgr->ConnectInput  (esdfilter,  0, cin_esd  );\r
-    mgr->ConnectOutput (esdfilter,  0, cout_aod );\r
-\r
-    mgr->ConnectInput  (esdmuonfilter,  0, cin_esd);\r
-    mgr->ConnectOutput (esdmuonfilter,  0, cout_aod );\r
\r
-    //\r
-    // Run the analysis\r
-    //    \r
-    if (mgr->InitAnalysis()) {\r
-        mgr->PrintStatus();\r
-        mgr->StartAnalysis("local",chain);\r
-    }   \r
-}\r
-\r
-//______________________________________________________________________________\r
-void SetupPar(char* pararchivename)\r
-{\r
-    if (pararchivename) {\r
-       char processline[1024];\r
-       sprintf(processline,".! tar xvzf %s.par",pararchivename);\r
-       gROOT->ProcessLine(processline);\r
-       TString ocwd = gSystem->WorkingDirectory();\r
-       gSystem->ChangeDirectory(pararchivename);\r
-       \r
-       // check for BUILD.sh and execute\r
-       if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) {\r
-           printf("*******************************\n");\r
-           printf("*** Building PAR archive    ***\n");\r
-           printf("*******************************\n");\r
-           \r
-           if (gSystem->Exec("PROOF-INF/BUILD.sh")) {\r
-               Error("runProcess","Cannot Build the PAR Archive! - Abort!");\r
-               return -1;\r
-           }\r
-       }\r
-       // check for SETUP.C and execute\r
-       if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) {\r
-           printf("*******************************\n");\r
-           printf("*** Setup PAR archive       ***\n");\r
-           printf("*******************************\n");\r
-           gROOT->Macro("PROOF-INF/SETUP.C");\r
-       }\r
-       \r
-       gSystem->ChangeDirectory(ocwd.Data());\r
-   printf("Current dir: %s\n", ocwd.Data());\r
-    } \r
-}\r
+void AnalysisTrainMuonLocal(char* filein = "AliESDs.root", char* fileout = "AliAOD.root")
+
+// Macro to produce a generic AOD starting from an ESD file. 
+// The AOD is filled with two tasks: 
+// 1- with the first one (AliAnalysisTaskESDfilter), 
+//    all the branches of the AOD are filled apart from the muons. 
+// 2- with the second task (AliAnalysisTaskESDMuonFilter) 
+//    muons tracks are added to the tracks branch 
+// There is the possibility to apply cuts on the muon tracks in order 
+// to reject muons before filling the AOD
+// This macro works locally
+// R. Arnaldi 5/5/08
+
+{
+    gSystem->Load("libTree.so");
+    gSystem->Load("libGeom.so");
+    gSystem->Load("libVMC.so");
+    gSystem->Load("libPhysics.so");
+
+    // If analysis is .par based:
+
+    // Common packages
+    SetupPar("STEERBase");
+    SetupPar("ESD");
+    SetupPar("AOD");
+    SetupPar("ANALYSIS");
+    SetupPar("ANALYSISalice");
+    // Analysis-specific packages
+    SetupPar("PWG3muon");      
+  
+    // Input ESD file
+    TChain* chain = new TChain("esdTree");  
+    chain->Add(filein);
+   
+    // Define the analysis manager
+    AliAnalysisManager *mgr  = new AliAnalysisManager("Analysis Train", "Analysis train");
+    
+    // ESD input handler
+    AliESDInputHandler *esdHandler = new AliESDInputHandler();
+    esdHandler->SetInactiveBranches("FMD CaloCluster");
+    
+    // AOD output handler
+    AliAODHandler* aodHandler   = new AliAODHandler();
+    aodHandler->SetOutputFileName(fileout);
+
+    mgr->SetInputEventHandler(esdHandler);
+    mgr->SetOutputEventHandler(aodHandler);
+    
+    // Set of cuts for the ESD filter
+    // 
+    // standard cut
+    AliESDtrackCuts* esdTrackCutsL = new AliESDtrackCuts("AliESDtrackCuts", "Loose");
+    esdTrackCutsL->SetMinNClustersTPC(50);
+    esdTrackCutsL->SetMaxChi2PerClusterTPC(3.5);
+    esdTrackCutsL->SetMaxCovDiagonalElements(2,2,0.5,0.5,2);
+    esdTrackCutsL->SetRequireTPCRefit(kTRUE);
+    esdTrackCutsL->SetMinNsigmaToVertex(3);
+    esdTrackCutsL->SetRequireSigmaToVertex(kTRUE);
+    esdTrackCutsL->SetAcceptKingDaughters(kFALSE);
+    //
+    // hard cuts
+    AliESDtrackCuts* esdTrackCutsH = new AliESDtrackCuts("AliESDtrackCuts", "Hard");
+    esdTrackCutsH->SetMinNClustersTPC(100);
+    esdTrackCutsH->SetMaxChi2PerClusterTPC(2.0);
+    esdTrackCutsH->SetMaxCovDiagonalElements(2,2,0.5,0.5,2);
+    esdTrackCutsH->SetRequireTPCRefit(kTRUE);
+    esdTrackCutsH->SetMinNsigmaToVertex(2);
+    esdTrackCutsH->SetRequireSigmaToVertex(kTRUE);
+    esdTrackCutsH->SetAcceptKingDaughters(kFALSE);
+    esdTrackCutsH->SetPRange(0.,2.);
+    //
+    //  muon cuts
+    AliESDMuonTrackCuts* esdMuonTrackCuts = new AliESDMuonTrackCuts("AliESDMuonTrackCuts", "test");
+    esdMuonTrackCuts->SetPRange(0.,20.);
+    //esdMuonTrackCuts->SetPtRange(0.,0.5);   // examples of kinematic cuts that can be applied
+    esdMuonTrackCuts->SetHistogramsOn(kTRUE);  // methods to draw control histos
+    esdMuonTrackCuts->DefineHistograms();
+    esdMuonTrackCuts->DrawHistograms();
+    
+    // track filter (to reject tracks not surviving the cuts - refers to all particles apart from muons)
+    AliAnalysisFilter* trackFilter = new AliAnalysisFilter("trackFilter");
+    trackFilter->AddCuts(esdTrackCutsH);
+    
+    // muon track filter  (to reject muon tracks not surviving the cuts)
+    AliAnalysisFilter* trackMuonFilter = new AliAnalysisFilter("trackMuonFilter");
+    trackMuonFilter->AddCuts(esdMuonTrackCuts);
+
+    // ESD filter task putting standard info in the output generic AOD 
+    AliAnalysisTaskESDfilter *esdfilter = new AliAnalysisTaskESDfilter("ESD Filter");
+    //esdfilter->SetTrackFilter(trackFilter);
+    esdfilter->SetDebugLevel(10);
+    mgr->AddTask(esdfilter);
+    
+    // ESD filter task putting muon info in the output generic AOD 
+    AliAnalysisTaskESDMuonFilter *esdmuonfilter = new AliAnalysisTaskESDMuonFilter("ESD Muon Filter");
+    esdmuonfilter->SetTrackFilter(trackMuonFilter);
+    mgr->AddTask(esdmuonfilter);
+
+    // Containers for input/output
+    AliAnalysisDataContainer *cin_esd = mgr->CreateContainer("cESD",TChain::Class(), 
+                                                            AliAnalysisManager::kInputContainer);
+    // Output AOD container. 
+    AliAnalysisDataContainer *cout_aod = mgr->CreateContainer("cAOD", TTree::Class(),
+                                                             AliAnalysisManager::kOutputContainer, "default");
+                                                             
+    // Connect containers to tasks slots
+    mgr->ConnectInput  (esdfilter,  0, cin_esd  );
+    mgr->ConnectOutput (esdfilter,  0, cout_aod );
+
+    mgr->ConnectInput  (esdmuonfilter,  0, cin_esd);
+    mgr->ConnectOutput (esdmuonfilter,  0, cout_aod );
+    //
+    // Run the analysis
+    //    
+    if (mgr->InitAnalysis()) {
+        mgr->PrintStatus();
+        mgr->StartAnalysis("local",chain);
+    }   
+}
+
+//______________________________________________________________________________
+void SetupPar(char* pararchivename)
+{
+    if (pararchivename) {
+       char processline[1024];
+       sprintf(processline,".! tar xvzf %s.par",pararchivename);
+       gROOT->ProcessLine(processline);
+       TString ocwd = gSystem->WorkingDirectory();
+       gSystem->ChangeDirectory(pararchivename);
+       
+       // check for BUILD.sh and execute
+       if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) {
+           printf("*******************************\n");
+           printf("*** Building PAR archive    ***\n");
+           printf("*******************************\n");
+           
+           if (gSystem->Exec("PROOF-INF/BUILD.sh")) {
+               Error("runProcess","Cannot Build the PAR Archive! - Abort!");
+               return -1;
+           }
+       }
+       // check for SETUP.C and execute
+       if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) {
+           printf("*******************************\n");
+           printf("*** Setup PAR archive       ***\n");
+           printf("*******************************\n");
+           gROOT->Macro("PROOF-INF/SETUP.C");
+       }
+       
+       gSystem->ChangeDirectory(ocwd.Data());
+   printf("Current dir: %s\n", ocwd.Data());
+    } 
+}
diff --git a/PWG3/muon/RunESDMuonFilter.C b/PWG3/muon/RunESDMuonFilter.C
deleted file mode 100644 (file)
index 4a34ab3..0000000
+++ /dev/null
@@ -1,89 +0,0 @@
-// Macro to run AliAnalysisTaskESDMuonFilter
-//
-
-void RunESDMuonFilter(char* filein = "AliESDs.root", char* fileout = "AliMuonAOD.root" ){
-     
-    gSystem->Load("libTree.so"); 
-    gSystem->Load("libGeom.so");
-    gSystem->Load("libVMC.so");
-    gSystem->Load("libPhysics");
-    
-    // for analysis .par file based
-    
-    setupPar("STEERBase");
-    setupPar("ESD");
-    setupPar("AOD");
-    setupPar("ANALYSIS");
-    setupPar("ANALYSISalice");
-    setupPar("PWG3base");
-    
-    // Input ESD file
-    TChain* chain = new TChain("esdTree");  
-    chain->Add(filein);
-    AliESDInputHandler* esdHandler = new AliESDInputHandler();
-
-    // Make aod output handler
-    AliAODHandler* aodHandler = new AliAODHandler();
-    aodHandler->SetOutputFileName(fileout);
-    
-    // Make the analysis manager
-    AliAnalysisManager *mgr  = new AliAnalysisManager("Muon AOD Manager", "Muon AOD Manager");
-    mgr->SetInputEventHandler(esdHandler);
-    mgr->SetOutputEventHandler(aodHandler);
-    mgr-> SetDebugLevel(10);
-    
-    // Task for MUON AOD generation 
-    AliAnalysisTaskESDMuonFilter *esdfilter = new AliAnalysisTaskESDMuonFilter("ESD Muon Filter");
-    mgr->AddTask(esdfilter);
-  
-    // Create containers for input/output
-    AliAnalysisDataContainer *cinput1 = mgr->CreateContainer("cESD",TChain::Class(), 
-                                                            AliAnalysisManager::kInputContainer);
-    AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("cAOD", TTree::Class(),
-                                                             AliAnalysisManager::kOutputContainer, "default");
-
-    mgr->ConnectInput  (esdfilter,  0, cinput1  );
-    mgr->ConnectOutput (esdfilter,  0, coutput1 );
-    
-    // Run the analysis    
-    mgr->InitAnalysis();
-    mgr->PrintStatus();
-    mgr->StartAnalysis("local",chain);
-}
-
-
-Int_t setupPar(const char* pararchivename) {
-  ///////////////////
-  // Setup PAR File//
-  ///////////////////
-  if (pararchivename) {
-    char processline[1024];
-    sprintf(processline,".! tar xvzf %s.par",pararchivename);
-    gROOT->ProcessLine(processline);
-    const char* ocwd = gSystem->WorkingDirectory();
-    gSystem->ChangeDirectory(pararchivename);
-
-    // check for BUILD.sh and execute
-    if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) {
-      printf("*******************************\n");
-      printf("*** Building PAR archive    ***\n");
-      printf("*******************************\n");
-
-      if (gSystem->Exec("PROOF-INF/BUILD.sh")) {
-        Error("runAnalysis","Cannot Build the PAR Archive! - Abort!");
-        return -1;
-      }
-    }
-    // check for SETUP.C and execute
-    if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) {
-      printf("*******************************\n");
-      printf("*** Setup PAR archive       ***\n");
-      printf("*******************************\n");
-      gROOT->Macro("PROOF-INF/SETUP.C");
-    }
-
-    gSystem->ChangeDirectory("../");
-  }
-
-  return 1;
-}