]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - MUON/DecodeRecoCocktail.C
Study jet p_T vs. total mult
[u/mrichter/AliRoot.git] / MUON / DecodeRecoCocktail.C
index 0b9f095cdd628e5bee324386b49fb995e8886cc1..89f04146c94d29ed0f72d68e802762becccab002 100644 (file)
 
 /* $Id$ */
 
-// A. De Falco, H. Woehri, INFN Cagliari, July 2006
-// This macro reads the generation/reconstruction files in the 
-// input directory (dirname), 
-// builds a tree that contains, for each event, an array of muons and dimuons
-// (TClonesArrays of AliMUONTrackLight and AliMUONPairLight objects) 
-// and writes the tree in an output file (outFileName) 
-// Note that if the path for the output file is not explicitly specified, 
-// it will be written in the directory containing the generation/reconstruction
+/// \ingroup macros
+/// \file DecodeRecoCocktail.C
+/// \brief add brief description
+///
+/// \author A. De Falco, H. Woehri, INFN Cagliari, July 2006
+///
+/// This macro reads the generation/reconstruction files in the 
+/// input directories (recodir and simdir), 
+/// builds a tree that contains, for each event, an array of muons and dimuons
+/// (TClonesArrays of AliMUONTrackLight and AliMUONPairLight objects) 
+/// and writes the tree in an output file (outFileName) 
+/// Note that if the path for the output file is not explicitly specified, 
+/// it will be written in the directory containing the generation/reconstruction
+/// 27-Nov-2006: modified by in order to loop on files
+///
+/// 13 Nov 2007:
+/// Updated this macro to work with new version of AliMUONRecoCheck. Also, we are
+/// now fetching reconstructed track data from ESD and not the track tree, because
+/// the AliMUONTrack objects for reconstructed tracks are no longer stored on disk.
+///  - Artur Szostak <artursz@iafrica.com>
 
-#include <iostream>
+#if !defined(__CINT__) || defined(__MAKECINT__)
+#include <Riostream.h>
 #include <TClonesArray.h>
 #include <TFile.h>
 #include <TTree.h>
 #include <TLorentzVector.h>
 #include <TParticle.h>
 #include <TSystem.h>
-#include <AliMUONRecoCheck.h>
-#include <AliMUONTrack.h>
-#include <AliMUONTrackParam.h>
-#include <AliRunLoader.h>
+#include "AliMUONRecoCheck.h"
+#include "AliMUONTrack.h"
 #include "AliMUONTrackLight.h"
 #include "AliMUONPairLight.h"
+#include "AliMUONVTrackStore.h"
+#include "AliESDEvent.h"
+#include "AliESDVertex.h"
+#include "AliESDMuonTrack.h"
+#include "AliMCEventHandler.h"
+#include "AliMCEvent.h"
+#include "AliStack.h"
+/*TODO: need to update this with changes made to ITS
+#include "AliITSVertexerPPZ.h"
+#include "AliITSLoader.h"
+*/
+#endif
 
-void DecodeRecoCocktail(char* dirname=".", char* outFileName = "MuonLight.root"){ 
-  const char *startingDir = gSystem->pwd(); 
-  gSystem->cd(dirname); 
-  TFile *fout = new TFile(outFileName,"recreate"); 
 
-  AliMUONRecoCheck *rc = new AliMUONRecoCheck("galice.root");
-  AliRunLoader *runLoader = rc->GetRunLoader();
-  
+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.
+  )
+{
+/// \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.
+
+  char startingDir[200]; 
+  sprintf (startingDir,"%s",gSystem->pwd()); 
+  gSystem->cd(recodir); 
+
   TClonesArray *muonArray   = new TClonesArray("AliMUONTrackLight",10); 
   TClonesArray *dimuonArray = new TClonesArray("AliMUONPairLight",10); 
   TTree *treeOut = new TTree("tree","tree"); 
 
+  TFile *fout = new TFile(outFileName,"recreate"); 
   fout->cd(); 
+
   treeOut->Branch("muons",&muonArray); 
   treeOut->Branch("dimuons",&dimuonArray); 
   
-  runLoader->LoadKinematics("READ");;
-  Int_t nev = runLoader->GetNumberOfEvents(); 
+  AliMUONRecoCheck *rc = new AliMUONRecoCheck("AliESDs.root", simdir);
+  
+  /*TODO: need to update this with changes made to ITS
+  AliITSLoader* ITSloader =  (AliITSLoader*) runLoaderSim->GetLoader("ITSLoader");
+  AliITSVertexerPPZ *dovert = 0; 
+  if (ITSloader) { 
+    dovert = new AliITSVertexerPPZ("default",0,0);
+    // dovert->SetDebug(0);
+    dovert->SetDiffPhiMax(0.05);
+    dovert->SetWindow(3.);
+  }
+  AliESDVertex *vert = 0;
+  */
   
   TLorentzVector v; 
+  Int_t nev = rc->NumberOfEvents();
   for(Int_t ievent = 0; ievent < nev; ievent++){ // loop over events 
-    printf ("Event %d of %d\n",ievent,nev);
+    
+    /*TODO: need to update this with changes made to ITS
+    runLoaderSim->GetHeader();
+    if (ITSloader) { 
+      vert = dovert->FindVertexForCurrentEvent(ievent);
+    }
+    */
+    
+    //printf ("Event %d of %d\n",ievent,nev);
     muonArray->Clear();     // clean muon and dimuon arrays 
     dimuonArray->Clear(); 
-    runLoader->GetEvent(ievent);
-    rc->ResetTracks();
-    rc->MakeTrackRef(); // make reconstructable tracks
     
-    TClonesArray * trackRecoArray = rc->GetTrackReco();
-    TClonesArray * trackRefArray = rc->GetMuonTrackRef();
-    Int_t nTrackReco = trackRecoArray->GetEntriesFast();
-    Int_t nreftracks = 0; 
-    for (Int_t itrRec = 0; itrRec<nTrackReco; itrRec++) { 
-      // assign parameters concerning the reconstructed tracks
-      AliMUONTrack *trackReco = (AliMUONTrack *)trackRecoArray->At(itrRec);
-      AliMUONTrackLight *muLight = new AliMUONTrackLight(); 
-      AliMUONTrackParam *trPar = trackReco->GetTrackParamAtVertex(); 
-      muLight->SetCharge(Int_t(TMath::Sign(1.,trPar->GetInverseBendingMomentum())));
-      muLight->SetPxPyPz(trPar->Px(),trPar->Py(), trPar->Pz()); 
-      muLight->SetTriggered(trackReco->GetMatchTrigger()); 
-      Double_t xyz[3] = { trPar->GetNonBendingCoor(), 
-                         trPar->GetBendingCoor(), 
-                         trPar->GetZ()};
-      muLight->SetVertex(xyz); 
-      // find the reference track and store further information
-      TParticle *part = muLight->FindRefTrack(trackReco,trackRefArray,runLoader); 
-      if (part) { 
+    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();
+    
+    // loop over ESD tracks
+    Int_t nreftracks = 0;
+    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);
+      
+      // skip ghosts
+      if (!esdTrack->ContainTrackerData()) continue;
+      
+      // 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(runLoader, part);
-       muLight->PrintInfo("A");
-       //store the referenced track in the muonArray:
-       TClonesArray &muons = *muonArray;
-       new (muons[nreftracks++]) AliMUONTrackLight(*muLight);
+       muLight->FillMuonHistory(pstack, part);
+       //        muLight.PrintInfo("A");
+       
       }
-    }
+      
+    } // end esd track
     
     // now loop over muon pairs to build dimuons
     Int_t nmuons = muonArray->GetEntriesFast(); 
@@ -103,27 +169,24 @@ void DecodeRecoCocktail(char* dirname=".", char* outFileName = "MuonLight.root")
       AliMUONTrackLight* mu1 = (AliMUONTrackLight*) muonArray->At(itrRec1); 
       for(Int_t itrRec2 = itrRec1+1; itrRec2 < nmuons; itrRec2++){
        AliMUONTrackLight* mu2 = (AliMUONTrackLight*) muonArray->At(itrRec2); 
-       AliMUONPairLight *dimuLight = new AliMUONPairLight(); 
-       dimuLight->SetMuons(*mu1, *mu2);
-       //      if(dimuLight->GetCreationProcess() == 2){
-       //        printf ("#dimuon = %d (%d, %d) \n", ndimuons, itrRec1, itrRec2);
-       // dimuLight->PrintInfo("A");
-       //      }
-       //store the referenced track in the dimuonArray:
+       AliMUONPairLight dimuLight;
+       dimuLight.SetMuons(*mu1, *mu2);
        TClonesArray &dimuons = *dimuonArray;
-       new (dimuons[ndimuons++]) AliMUONPairLight(*dimuLight);
+       new (dimuons[ndimuons++]) AliMUONPairLight(dimuLight);
       }
     }
-    Int_t ndimu2 = dimuonArray->GetEntriesFast(); 
-    //     printf ("dimuonArray has %d entries\n",ndimu2); 
-    //    dimuonArray->Dump(); 
-    //     printf ("filling tree\n"); 
-    //     // fill the tree
-    treeOut->Fill(); 
-    //     printf ("done\n"); 
     
-  }
+    treeOut->Fill(); 
+  } 
+  
   fout->cd(); 
   treeOut->Write(); 
   gSystem->cd(startingDir); 
+  fout->Close();
+  delete muonArray;
+  delete dimuonArray;
+  delete treeOut; 
+  delete fout;
+  delete rc;
 }
+