bug fixed

[u/mrichter/AliRoot.git] / MUON / DecodeRecoCocktail.C

diff --git a/MUON/DecodeRecoCocktail.C b/MUON/DecodeRecoCocktail.C
index 6f00fe90d4469321f137670849d61a9e68dd739f..89f04146c94d29ed0f72d68e802762becccab002 100644 (file)
--- a/MUON/DecodeRecoCocktail.C
+++ b/MUON/DecodeRecoCocktail.C
@@ -45,38 +45,33 @@
 #include <TLorentzVector.h>
 #include <TParticle.h>
 #include <TSystem.h>
-#include <TGeoManager.h>
 #include "AliMUONRecoCheck.h"
 #include "AliMUONTrack.h"
-#include "AliMUONTrackParam.h"
 #include "AliMUONTrackLight.h"
 #include "AliMUONPairLight.h"
 #include "AliMUONVTrackStore.h"
-#include "AliMUONTrackExtrap.h"
 #include "AliESDEvent.h"
 #include "AliESDVertex.h"
+#include "AliESDMuonTrack.h"
 #include "AliMCEventHandler.h"
 #include "AliMCEvent.h"
-#include "AliMagF.h"
+#include "AliStack.h"
 /*TODO: need to update this with changes made to ITS
 #include "AliITSVertexerPPZ.h"
 #include "AliITSLoader.h"
 */
-#include "AliTracker.h"
 #endif
 
 
 void DecodeRecoCocktail(
     char* recodir=".",          // The directory containing galice.root for reconstructed data.
     char* simdir="generated/",  // The directory containing galice.root for simulated data.
-    char* outFileName = "MuonLight.root", // The output filename containing AliMUONTrackLight and AliMUONPairLight objects.
-    char* geoFilename = "geometry.root"  // The filename containing the geometry.
+    char* outFileName = "MuonLight.root" // The output filename containing AliMUONTrackLight and AliMUONPairLight objects.
   )
 {
 /// \param recodir      The directory containing galice.root for reconstructed data.
 /// \param simdir       The directory containing galice.root for simulated data.
 /// \param outFileName  The output filename containing AliMUONTrackLight and AliMUONPairLight objects.
-/// \param geoFilename  The filename containing the geometry.
 
   char startingDir[200]; 
   sprintf (startingDir,"%s",gSystem->pwd()); 
@@ -92,27 +87,7 @@ void DecodeRecoCocktail(
   treeOut->Branch("muons",&muonArray); 
   treeOut->Branch("dimuons",&dimuonArray); 
   
-  // Import TGeo geometry (needed by AliMUONTrackExtrap::ExtrapToVertex)
-  if (!gGeoManager) {
-    TGeoManager::Import(geoFilename);
-    if (!gGeoManager) {
-      Error("MUONmass_ESD", "getting geometry from file %s failed", geoFilename);
-      return;
-    }
-  }
-  
-  // set  mag field 
-  // waiting for mag field in CDB 
-  if (!TGeoGlobalMagField::Instance()->GetField()) {
-    printf("Loading field map...\n");
-    AliMagF* field = new AliMagF("Maps","Maps",2,1.,1., 10.,AliMagF::k5kG);
-    TGeoGlobalMagField::Instance()->SetField(field);
-  }
-  // set the magnetic field for track extrapolations
-  AliMUONTrackExtrap::SetField();
-
   AliMUONRecoCheck *rc = new AliMUONRecoCheck("AliESDs.root", simdir);
-  Int_t nev = rc->NumberOfEvents();
   
   /*TODO: need to update this with changes made to ITS
   AliITSLoader* ITSloader =  (AliITSLoader*) runLoaderSim->GetLoader("ITSLoader");
@@ -128,6 +103,7 @@ void DecodeRecoCocktail(
   
   TLorentzVector v; 
  
+  Int_t nev = rc->NumberOfEvents();
   for(Int_t ievent = 0; ievent < nev; ievent++){ // loop over events 
     
     /*TODO: need to update this with changes made to ITS
@@ -141,41 +117,50 @@ void DecodeRecoCocktail(
     muonArray->Clear();     // clean muon and dimuon arrays 
     dimuonArray->Clear(); 
     
-    AliMUONVTrackStore* recoTracks = rc->ReconstructedTracks(ievent);  // Use tracks from actual reconstruction, but this does not work anymore.
+    AliMUONVTrackStore* recoTracks = rc->ReconstructedTracks(ievent, kFALSE);  // Use tracks from actual reconstruction.
     //AliMUONVTrackStore* trackRefs = rc->ReconstructibleTracks(ievent);  // Only use reconstructible reference tracks.
     AliMUONVTrackStore* trackRefs = rc->TrackRefs(ievent);
     AliStack* pstack = (const_cast<AliMCEventHandler*>(rc->GetMCEventHandler()))->MCEvent()->Stack();
     
-    TIter next(recoTracks->CreateIterator());
-    AliMUONTrack* trackReco = NULL;
-    
-    Int_t nTrackReco = recoTracks->GetSize();
-    Int_t nTracksESD = rc->GetESDEvent()->GetNumberOfMuonTracks();
-    if (nTrackReco != nTracksESD) printf ("Tracks in recoTracks (%d) and in ESD (%d) do not match!\n", nTrackReco, nTracksESD);
+    // loop over ESD tracks
     Int_t nreftracks = 0;
-    Int_t itrRec = 0;
-    while ( (trackReco = static_cast<AliMUONTrack*>(next())) != NULL )
-    {
-      // assign parameters concerning the reconstructed tracks
-      AliMUONTrackLight muLight;
+    const AliESDEvent* esd = rc->GetESDEvent();
+    Int_t nTracks = (Int_t)esd->GetNumberOfMuonTracks() ;
+    for (Int_t iTrack = 0; iTrack <  nTracks;  iTrack++) {
+      
+      AliESDMuonTrack* esdTrack = esd->GetMuonTrack(iTrack);
       
-      muLight.FillFromESD(rc->GetESDEvent()->GetMuonTrack(itrRec));
-      // muLight.FillFromAliMUONTrack(trackReco);
+      // skip ghosts
+      if (!esdTrack->ContainTrackerData()) continue;
       
-      // find the reference track and store further information        
-      TParticle *part = muLight.FindRefTrack(trackReco, trackRefs, pstack);
-      if (part) { 
+      // find the corresponding MUON track
+      AliMUONTrack* trackReco = (AliMUONTrack*) recoTracks->FindObject(esdTrack->GetUniqueID());
+      
+      // try to match the reconstructed track with a simulated one
+      Int_t nMatchClusters = 0;
+      AliMUONTrack* matchedTrackRef = rc->FindCompatibleTrack(*trackReco, *trackRefs, nMatchClusters, kFALSE, 10.);
+      
+      if (matchedTrackRef) {
+       
+       //store new referenced track in the muonArray
+       AliMUONTrackLight* muLight = new ((*muonArray)[nreftracks++]) AliMUONTrackLight();
+       
+       // assign parameters concerning the reconstructed tracks
+       muLight->FillFromESD(esdTrack);
+       // muLight->FillFromAliMUONTrack(trackReco);
+       
+       // store further information related to the simulated track
+       muLight->SetTrackPythiaLine(matchedTrackRef->GetUniqueID());
+       TParticle *part = pstack->Particle(matchedTrackRef->GetUniqueID());
+       muLight->SetTrackPDGCode(part->GetPdgCode());
        v.SetPxPyPzE(part->Px(), part->Py(), part->Pz(), part->Energy());
-       muLight.SetPGen(v); 
-       muLight.FillMuonHistory(pstack, part);
+       muLight->SetPGen(v); 
+       muLight->FillMuonHistory(pstack, part);
        //        muLight.PrintInfo("A");
-       //store the referenced track in the muonArray:
-       TClonesArray &muons = *muonArray;
-       new (muons[nreftracks++]) AliMUONTrackLight(muLight);
-      } 
+       
+      }
       
-      itrRec++;
-    }  // end reco track
+    } // end esd track
     
     // now loop over muon pairs to build dimuons
     Int_t nmuons = muonArray->GetEntriesFast();