]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - MUON/MUONmassPlot_ESD.C
Simplifications following change of library loading and $ALICE_ROOT
[u/mrichter/AliRoot.git] / MUON / MUONmassPlot_ESD.C
index 7237e48b3e306f26052ab2a1c8e58518b509f365..1ebae5bcbd57ce774eaef6186aedfc12737039b7 100644 (file)
 
 // STEER includes
 #include "AliLog.h"
-#include "AliMagFMaps.h"
+#include "AliCDBManager.h"
 #include "AliESDEvent.h"
 #include "AliESDVertex.h"
-#include "AliTracker.h"
 #include "AliESDMuonTrack.h"
 
 // MUON includes
+#include "AliMUONCDB.h"
 #include "AliMUONTrackParam.h"
 #include "AliMUONTrackExtrap.h"
 #include "AliMUONESDInterface.h"
+#include "AliMUONConstants.h"
 #endif
 
 /// \ingroup macros
 /// \file MUONmassPlot_ESD.C
-/// \brief Macro MUONefficiency.C for ESD
+/// \brief Macro MUONmassPlot_ESD.C for ESD
 ///
 /// \author Ch. Finck, Subatech, April. 2004
 ///
 ///
 /// Add parameters and histograms for analysis 
 
-Bool_t MUONmassPlot(Int_t ExtrapToVertex = -1, char* geoFilename = "geometry.root", 
-                 Int_t FirstEvent = 0, Int_t LastEvent = 10000, char* esdFileName = "AliESDs.root", Int_t ResType = 553, 
-                  Float_t Chi2Cut = 100., Float_t PtCutMin = 1., Float_t PtCutMax = 10000.,
-                  Float_t massMin = 9.17,Float_t massMax = 9.77)
+Bool_t MUONmassPlot(const char* esdFileName = "AliESDs.root", const char* geoFilename = "geometry.root",
+                   const char* ocdbPath = "local://$ALICE_ROOT/OCDB",
+                   Int_t FirstEvent = 0, Int_t LastEvent = -1, Int_t ExtrapToVertex = -1,
+                   Int_t ResType = 553, Float_t Chi2Cut = 100., Float_t PtCutMin = 1.,
+                   Float_t PtCutMax = 10000., Float_t massMin = 9.17,Float_t massMax = 9.77)
 {
-/// \param ExtrapToVertex (default -1)
-///   -        <0: no extrapolation;
-///   -        =0: extrapolation to (0,0,0);
-///   -        >0: extrapolation to ESDVertex if available, else to (0,0,0)
-/// \param FirstEvent (default 0)
-/// \param LastEvent (default 0)
-/// \param ResType    553 for Upsilon, anything else for J/Psi (default 553)
-/// \param Chi2Cut    to keep only tracks with chi2 per d.o.f. < Chi2Cut (default 100)
-/// \param PtCutMin   to keep only tracks with transverse momentum > PtCutMin (default 1)
-/// \param PtCutMax   to keep only tracks with transverse momentum < PtCutMax (default 10000)
-/// \param massMin   (default 9.17 for Upsilon) 
-/// \param massMax   (default 9.77 for Upsilon);  
-///         to calculate the reconstruction efficiency for resonances with invariant mass
-///         massMin < mass < massMax.
+  /// \param FirstEvent (default 0)
+  /// \param LastEvent (default 10000)
+  /// \param ExtrapToVertex (default -1)
+  ///   -      <0: no extrapolation;
+  ///   -      =0: extrapolation to (0,0,0);
+  ///   -      >0: extrapolation to ESDVertex if available, else to (0,0,0)
+  /// \param ResType    553 for Upsilon, anything else for J/Psi (default 553)
+  /// \param Chi2Cut    to keep only tracks with chi2 per d.o.f. < Chi2Cut (default 100)
+  /// \param PtCutMin   to keep only tracks with transverse momentum > PtCutMin (default 1)
+  /// \param PtCutMax   to keep only tracks with transverse momentum < PtCutMax (default 10000)
+  /// \param massMin   (default 9.17 for Upsilon) 
+  /// \param massMax   (default 9.77 for Upsilon);  
+  ///         to calculate the reconstruction efficiency for resonances with invariant mass
+  ///         massMin < mass < massMax.
 
   cout << "MUONmassPlot " << endl;
   cout << "FirstEvent " << FirstEvent << endl;
-  cout << "LastEvent " << LastEvent << endl;
+  cout << "LastEvent " << ((LastEvent>=0) ? Form("%d",LastEvent) : "all") << endl;
   cout << "ResType " << ResType << endl;
   cout << "Chi2Cut " << Chi2Cut << endl;
   cout << "PtCutMin " << PtCutMin << endl;
@@ -118,7 +120,7 @@ Bool_t MUONmassPlot(Int_t ExtrapToVertex = -1, char* geoFilename = "geometry.roo
   Double_t fPxRec1, fPyRec1, fPzRec1, fE1;
   Double_t fPxRec2, fPyRec2, fPzRec2, fE2;
 
-  Int_t ntrackhits, nevents;
+  Int_t ntrackhits;
   Double_t fitfmin;
   Double_t fZVertex=0;
   Double_t fYVertex=0;
@@ -127,7 +129,7 @@ Bool_t MUONmassPlot(Int_t ExtrapToVertex = -1, char* geoFilename = "geometry.roo
   Double_t errYVtx=0;
  
   TLorentzVector fV1, fV2, fVtot;
-
+  
   // Import TGeo geometry (needed by AliMUONTrackExtrap::ExtrapToVertex)
   if (!gGeoManager) {
     TGeoManager::Import(geoFilename);
@@ -137,48 +139,55 @@ Bool_t MUONmassPlot(Int_t ExtrapToVertex = -1, char* geoFilename = "geometry.roo
     }
   }
   
-  // set mag field
-  // waiting for mag field in CDB 
-  printf("Loading field map...\n");
-  AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 1, 1., 10., AliMagFMaps::k5kG);
-  AliTracker::SetFieldMap(field, kFALSE);
-
   // open the ESD file
   TFile* esdFile = TFile::Open(esdFileName);
   if (!esdFile || !esdFile->IsOpen()) {
     Error("MUONmass_ESD", "opening ESD file %s failed", esdFileName);
     return kFALSE;
   }
-  
   AliESDEvent* esd = new AliESDEvent();
   TTree* tree = (TTree*) esdFile->Get("esdTree");
   if (!tree) {
-    Error("CheckESD", "no ESD tree found");
+    Error("MUONmass_ESD", "no ESD tree found");
     return kFALSE;
   }
-//  tree->SetBranchAddress("ESD", &esd);
   esd->ReadFromTree(tree);
   
-  nevents = (Int_t)tree->GetEntries();
+  // get run number
+  if (tree->GetEvent(0) <= 0) {
+    Error("MUONmass_ESD", "no ESD object found for event 0");
+    return kFALSE;
+  }
+  Int_t runNumber = esd->GetRunNumber();
+  
+  // load necessary data from OCDB
+  AliCDBManager::Instance()->SetDefaultStorage(ocdbPath);
+  AliCDBManager::Instance()->SetSpecificStorage("GRP/GRP/Data","local://.");
+  AliCDBManager::Instance()->SetRun(runNumber);
+  if (!AliMUONCDB::LoadField()) return kFALSE;
+  
+  // set the magnetic field for track extrapolations
+  AliMUONTrackExtrap::SetField();
   
   AliMUONTrackParam trackParam;
-
+  AliMUONTrackParam trackParamAtAbsEnd;
+  
   // Loop over events
-  for (Int_t iEvent = FirstEvent; iEvent <= TMath::Min(LastEvent, nevents - 1); iEvent++) {
+  Int_t nevents = (LastEvent >= 0) ? TMath::Min(LastEvent, (Int_t)tree->GetEntries()-1) : (Int_t)tree->GetEntries()-1;
+  for (Int_t iEvent = FirstEvent; iEvent <= nevents; iEvent++) {
 
     // get the event summary data
-    tree->GetEvent(iEvent);
-    if (!esd) {
-      Error("CheckESD", "no ESD object found for event %d", iEvent);
+    if (tree->GetEvent(iEvent) <= 0) {
+      Error("MUONmass_ESD", "no ESD object found for event %d", iEvent);
       return kFALSE;
     }
-
+    
     // get the SPD reconstructed vertex (vertexer) and fill the histogram
     AliESDVertex* Vertex = (AliESDVertex*) esd->GetVertex();
     if (Vertex->GetNContributors()) {
-      fZVertex = Vertex->GetZv();
-      fYVertex = Vertex->GetYv();
-      fXVertex = Vertex->GetXv();
+      fZVertex = Vertex->GetZ();
+      fYVertex = Vertex->GetY();
+      fXVertex = Vertex->GetX();
       errXVtx = Vertex->GetXRes();
       errYVtx = Vertex->GetYRes();
     }
@@ -189,8 +198,6 @@ Bool_t MUONmassPlot(Int_t ExtrapToVertex = -1, char* geoFilename = "geometry.roo
     //    printf("\n Nb of events analysed: %d\r",iEvent);
     //      cout << " number of tracks: " << nTracks  <<endl;
   
-    // set the magnetic field for track extrapolations
-    AliMUONTrackExtrap::SetField(AliTracker::GetFieldMap());
     // loop over all reconstructed tracks (also first track of combination)
     for (Int_t iTrack = 0; iTrack <  nTracks;  iTrack++) {
 
@@ -210,6 +217,14 @@ Bool_t MUONmassPlot(Int_t ExtrapToVertex = -1, char* geoFilename = "geometry.roo
        AliMUONESDInterface::SetParamAtVertex(trackParam, *muonTrack); // put the new parameters in this copy of AliESDMuonTrack
       }
 
+      // compute track position at the end of the absorber
+      AliMUONESDInterface::GetParamAtFirstCluster(*muonTrack, trackParamAtAbsEnd);
+      AliMUONTrackExtrap::ExtrapToZ(&trackParamAtAbsEnd, AliMUONConstants::AbsZEnd());
+      Double_t xAbs = trackParamAtAbsEnd.GetNonBendingCoor();
+      Double_t yAbs = trackParamAtAbsEnd.GetBendingCoor();
+      Double_t dAbs1 = TMath::Sqrt(xAbs*xAbs + yAbs*yAbs);
+      Double_t aAbs1 = TMath::ATan(-dAbs1/AliMUONConstants::AbsZEnd()) * TMath::RadToDeg();
+      
       fCharge1 = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
       
       muonTrack->LorentzP(fV1);
@@ -233,7 +248,7 @@ Bool_t MUONmassPlot(Int_t ExtrapToVertex = -1, char* geoFilename = "geometry.roo
 
       // condition for good track (Chi2Cut and PtCut)
 
-      if ((ch1 < Chi2Cut) && (pt1 > PtCutMin) && (pt1 < PtCutMax)) {
+//      if ((ch1 < Chi2Cut) && (pt1 > PtCutMin) && (pt1 < PtCutMax)) {
 
        // fill histos hPtMuon and hChi2PerDof
        hPtMuon->Fill(pt1);
@@ -266,6 +281,14 @@ Bool_t MUONmassPlot(Int_t ExtrapToVertex = -1, char* geoFilename = "geometry.roo
            AliMUONESDInterface::SetParamAtVertex(trackParam, *muonTrack2); // put the new parameters in this copy of AliESDMuonTrack
          }
          
+         // compute track position at the end of the absorber
+         AliMUONESDInterface::GetParamAtFirstCluster(*muonTrack2, trackParamAtAbsEnd);
+         AliMUONTrackExtrap::ExtrapToZ(&trackParamAtAbsEnd, AliMUONConstants::AbsZEnd());
+         xAbs = trackParamAtAbsEnd.GetNonBendingCoor();
+         yAbs = trackParamAtAbsEnd.GetBendingCoor();
+         Double_t dAbs2 = TMath::Sqrt(xAbs*xAbs + yAbs*yAbs);
+         Double_t aAbs2 = TMath::ATan(-dAbs2/AliMUONConstants::AbsZEnd()) * TMath::RadToDeg();
+         
          fCharge2 = Int_t(TMath::Sign(1.,muonTrack2->GetInverseBendingMomentum()));
 
          muonTrack2->LorentzP(fV2);
@@ -280,7 +303,16 @@ Bool_t MUONmassPlot(Int_t ExtrapToVertex = -1, char* geoFilename = "geometry.roo
          Float_t ch2 = fitfmin  / (2.0 * ntrackhits - 5);
 
          // condition for good track (Chi2Cut and PtCut)
-         if ((ch2 < Chi2Cut) && (pt2 > PtCutMin)  && (pt2 < PtCutMax)) {
+//       if ((ch2 < Chi2Cut) && (pt2 > PtCutMin)  && (pt2 < PtCutMax)) {
+//       if (aAbs1 < 2. || aAbs2 < 2.) {
+//       if (aAbs1 > 2. && aAbs2 > 2. && (aAbs1 < 2.1 || aAbs2 < 2.1)) {
+//       if (aAbs1 > 2. && aAbs2 > 2. && (dAbs1 < 17.8 || dAbs2 < 17.8)) {
+//       if (dAbs1 > 17.8 && dAbs2 > 17.8 && (dAbs1 < 18. || dAbs2 < 18.)) {
+//       if (dAbs1 > 18. && dAbs2 > 18. && (dAbs1 < 18.2 || dAbs2 < 18.2)) {
+//       if (dAbs1 > 18.2 && dAbs2 > 18.2 && (aAbs1 < 2.1 || aAbs2 < 2.1)) {
+//       if (dAbs1 > 18. && dAbs2 > 18.) {
+//       if (aAbs1 > 2.1 && aAbs2 > 2.1) {
+//       if (muonTrack->GetMatchTrigger() && muonTrack2->GetMatchTrigger()) {
 
            // condition for opposite charges
            if ((fCharge1 * fCharge2) == -1) {
@@ -310,10 +342,10 @@ Bool_t MUONmassPlot(Int_t ExtrapToVertex = -1, char* geoFilename = "geometry.roo
              }
 
            } // if (fCharge1 * fCharge2) == -1)
-         } // if ((ch2 < Chi2Cut) && (pt2 > PtCutMin) && (pt2 < PtCutMax))
+//       } // if ((ch2 < Chi2Cut) && (pt2 > PtCutMin) && (pt2 < PtCutMax))
          delete muonTrack2;
        } //  for (Int_t iTrack2 = iTrack + 1; iTrack2 < iTrack; iTrack2++)
-      } // if (ch1 < Chi2Cut) && (pt1 > PtCutMin)&& (pt1 < PtCutMax) )
+//      } // if (ch1 < Chi2Cut) && (pt1 > PtCutMin)&& (pt1 < PtCutMax) )
       delete muonTrack;
     } // for (Int_t iTrack = 0; iTrack < nrectracks; iTrack++)