New code to work both with the ESD and MC (Mercedes)
authorhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 21 Jul 2006 15:49:44 +0000 (15:49 +0000)
committerhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 21 Jul 2006 15:49:44 +0000 (15:49 +0000)
JETAN/AliJetESDReader.cxx
JETAN/AliJetESDmcReader.cxx [new file with mode: 0644]
JETAN/AliJetESDmcReader.h [new file with mode: 0644]
JETAN/AliJetReaderHeader.h
JETAN/JetAnalysisLinkDef.h
JETAN/libJETAN.pkg

index fbc81fa..2086bba 100755 (executable)
  * about the suitability of this software for any purpose. It is          *
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
+
+//------------------------------------------------------------------------- 
 // Jet ESD Reader 
 // ESD reader for jet analysis
 // Author: Mercedes Lopez Noriega (mercedes.lopez.noriega@cern.ch)
+//------------------------------------------------------------------------- 
+
 
 #include <Riostream.h>
 #include <TSystem.h>
 
 ClassImp(AliJetESDReader)
 
-AliJetESDReader::AliJetESDReader()
+  AliJetESDReader::AliJetESDReader()
 {
   // Constructor    
-  printf("\nIn reader constructor\n");
   fReaderHeader = 0x0;
   fMass = 0;
   fSign = 0;
@@ -47,12 +49,10 @@ AliJetESDReader::~AliJetESDReader()
 //____________________________________________________________________________
 
 void AliJetESDReader::OpenInputFiles()
+
 {
-  // Open input files
-  printf("\nOpening files\n");
   // chain for the ESDs
   fChain   = new TChain("esdTree");
-  fChainMC = new TChain("mcStackTree");
 
   // get directory and pattern name from the header
   const char* dirName=fReaderHeader->GetDirectory();
@@ -61,24 +61,26 @@ void AliJetESDReader::OpenInputFiles()
   // Add files matching patters to the chain
   void *dir  = gSystem->OpenDirectory(dirName);
   const char *name = 0x0;
+  int nesd = fReaderHeader->GetNesd();
+  int a = 0;
   while ((name = gSystem->GetDirEntry(dir))){
+    if (a>=nesd) continue;
     if (strstr(name,pattern)){
       printf("Adding %s\n",name);
       char path[256];
       sprintf(path,"%s/%s",dirName,name);
       fChain->AddFile(path,-1);
-      fChainMC->AddFile(path,-1);
-     }
+      a++;
+    }
   }
-
-  gSystem ->FreeDirectory(dir);
-  fChain  ->SetBranchAddress("ESD",    &fESD);
-  fChainMC->SetBranchAddress("Header", &fAliHeader);
-  fChainMC->SetBranchAddress("Stack",  &fArrayMC);
-
+  printf("%d ESDs added\n",a);
+  
+  gSystem->FreeDirectory(dir);
+  fChain->SetBranchAddress("ESD",    &fESD);
+  
   int nMax = fChain->GetEntries(); 
   printf("\nTotal number of events in chain= %d",nMax);
-
+  
   // set number of events in header
   if (fReaderHeader->GetLastEvent() == -1)
     fReaderHeader->SetLastEvent(nMax);
@@ -92,62 +94,59 @@ void AliJetESDReader::OpenInputFiles()
 
 void AliJetESDReader::FillMomentumArray(Int_t event)
 {
-  // Fill the momentum array for each track
+  // Fill momentum array
+
   Int_t goodTrack = 0;
   Int_t nt = 0;
   Float_t pt, eta;
   TVector3 p3;
-
+  
   // clear momentum array
   ClearArray();
-
+  
   // get event from chain
   fChain->GetEntry(event);
-  fChainMC->GetEntry(event);
-
+  
   // get number of tracks in event (for the loop)
   nt = fESD->GetNumberOfTracks();
- // temporary storage of signal and pt cut flag
+  
+  // temporary storage of signal and pt cut flag
   Int_t* sflag  = new Int_t[nt];
   Int_t* cflag  = new Int_t[nt];
-
+  
   // get cuts set by user
   Float_t ptMin = fReaderHeader->GetPtCut();
   Float_t etaMin = fReaderHeader->GetFiducialEtaMin();
   Float_t etaMax = fReaderHeader->GetFiducialEtaMax();  
-
+  
   //loop over tracks
   for (Int_t it = 0; it < nt; it++) {
-      AliESDtrack *track = fESD->GetTrack(it);
-      UInt_t status = track->GetStatus();
-
-      Double_t mom[3];
-      track->GetPxPyPz(mom);
-
-      p3.SetXYZ(mom[0],mom[1],mom[2]);
-      pt = p3.Pt();
-      if (((status & AliESDtrack::kITSrefit) == 0) ||
-         ((status & AliESDtrack::kTPCrefit) == 0)) continue;    // quality check
-      if (((AliJetESDReaderHeader*) fReaderHeader)->ReadSignalOnly() 
-         && TMath::Abs(track->GetLabel()) > 10000)  continue;   // quality check
-      if (((AliJetESDReaderHeader*) fReaderHeader)->ReadBkgdOnly() 
-         && TMath::Abs(track->GetLabel()) < 10000)  continue;   // quality check
-      eta = p3.Eta();
-      if ( (eta > etaMax) || (eta < etaMin)) continue;           // checking eta cut
-
-      new ((*fMomentumArray)[goodTrack]) TLorentzVector(p3,p3.Mag());
-      sflag[goodTrack]=0;
-      if (TMath::Abs(track->GetLabel()) < 10000) sflag[goodTrack]=1;
-      cflag[goodTrack]=0;
-      if (pt > ptMin) cflag[goodTrack]=1;                       // pt cut
-      goodTrack++;
+    AliESDtrack *track = fESD->GetTrack(it);
+    UInt_t status = track->GetStatus();
+
+    Double_t mom[3];
+    track->GetPxPyPz(mom);
+    p3.SetXYZ(mom[0],mom[1],mom[2]);
+    pt = p3.Pt();
+    if (((status & AliESDtrack::kITSrefit) == 0) ||
+       ((status & AliESDtrack::kTPCrefit) == 0)) continue;    // quality check
+    if (((AliJetESDReaderHeader*) fReaderHeader)->ReadSignalOnly() 
+       && TMath::Abs(track->GetLabel()) > 10000)  continue;   // quality check
+    if (((AliJetESDReaderHeader*) fReaderHeader)->ReadBkgdOnly() 
+       && TMath::Abs(track->GetLabel()) < 10000)  continue;   // quality check
+    eta = p3.Eta();
+    if ( (eta > etaMax) || (eta < etaMin)) continue;           // checking eta cut
+    
+    new ((*fMomentumArray)[goodTrack]) TLorentzVector(p3,p3.Mag());
+    sflag[goodTrack]=0;
+    if (TMath::Abs(track->GetLabel()) < 10000) sflag[goodTrack]=1;
+    cflag[goodTrack]=0;
+    if (pt > ptMin) cflag[goodTrack]=1;                       // pt cut
+    goodTrack++;
   }
   // set the signal flags
   fSignalFlag.Set(goodTrack,sflag);
   fCutFlag.Set(goodTrack,cflag);
-
-  //printf("\nIn event %d, number of good tracks %d \n", event, goodTrack);
 }
 
 
diff --git a/JETAN/AliJetESDmcReader.cxx b/JETAN/AliJetESDmcReader.cxx
new file mode 100644 (file)
index 0000000..bad1a3e
--- /dev/null
@@ -0,0 +1,172 @@
+/**************************************************************************
+ * 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.                  *
+ **************************************************************************/
+// Jet ESD Reader 
+// ESD reader for jet analysis (it reads the esd and the MC trees)
+// Author: Mercedes Lopez Noriega (mercedes.lopez.noriega@cern.ch)
+
+#include <Riostream.h>
+#include <TSystem.h>
+#include <TLorentzVector.h>
+#include <TVector3.h>
+#include <TPDGCode.h>
+#include <TParticle.h>
+#include <TParticlePDG.h>
+#include <TMath.h>
+#include "AliJetESDmcReader.h"
+#include "AliJetESDReaderHeader.h"
+#include "AliESD.h"
+#include "AliESDtrack.h"
+
+ClassImp(AliJetESDmcReader)
+
+AliJetESDmcReader::AliJetESDmcReader()
+{
+  // Constructor    
+  fReaderHeader = 0x0;
+  fMass = 0;
+  fSign = 0;
+}
+
+//____________________________________________________________________________
+
+AliJetESDmcReader::~AliJetESDmcReader()
+{
+  // Destructor
+}
+
+//____________________________________________________________________________
+
+void AliJetESDmcReader::OpenInputFiles()
+
+{
+  // chain for the ESDs
+  fChain   = new TChain("esdTree");
+  fChainMC = new TChain("mcStackTree");
+
+  // get directory and pattern name from the header
+  const char* dirName=fReaderHeader->GetDirectory();
+  const char* pattern=fReaderHeader->GetPattern();
+    
+  // Add files matching patters to the chain
+  void *dir  = gSystem->OpenDirectory(dirName);
+  const char *name = 0x0;
+  int nesd = fReaderHeader->GetNesd();
+  int a = 0;
+  while ((name = gSystem->GetDirEntry(dir))){
+    if (a>=nesd) continue;
+    if (strstr(name,pattern)){
+      printf("Adding %s\n",name);
+      char path[256];
+      sprintf(path,"%s/%s",dirName,name);
+      fChain->AddFile(path,-1);
+      fChainMC->AddFile(path,-1);
+      a++;
+    }
+  }
+  printf("%d ESDs added\n",a);
+  
+  gSystem->FreeDirectory(dir);
+  fChain->SetBranchAddress("ESD",    &fESD);
+  fChainMC->SetBranchAddress("Header", &fAliHeader);
+  fChainMC->SetBranchAddress("Stack",  &fArrayMC);
+
+  int nMax = fChain->GetEntries(); 
+  printf("\nTotal number of events in chain= %d\n",nMax);
+
+  // set number of events in header
+  if (fReaderHeader->GetLastEvent() == -1)
+    fReaderHeader->SetLastEvent(nMax);
+  else {
+    Int_t nUsr = fReaderHeader->GetLastEvent();
+    fReaderHeader->SetLastEvent(TMath::Min(nMax,nUsr));
+  }
+}
+
+//____________________________________________________________________________
+
+void AliJetESDmcReader::FillMomentumArray(Int_t event)
+{
+  // Fill momentum array
+
+  Int_t goodTrack = 0;
+  Int_t nt = 0; 
+  Int_t pdgCode = 0;
+  Float_t pt, eta;
+  Float_t energyMC, pxMC, pyMC, pzMC, ptMC; // Monte Carlo
+  TVector3 p3;
+
+  // clear momentum array
+  ClearArray();
+
+  // get event from chain
+  fChain->GetEntry(event);
+  fChainMC->GetEntry(event);
+  // get number of tracks in event (for the loop)
+  nt = fESD->GetNumberOfTracks();
+ // temporary storage of signal and pt cut flag
+  Int_t* sflag  = new Int_t[nt];
+  Int_t* cflag  = new Int_t[nt];
+
+  // get cuts set by user
+  Float_t ptMin = fReaderHeader->GetPtCut();
+  Float_t etaMin = fReaderHeader->GetFiducialEtaMin();
+  Float_t etaMax = fReaderHeader->GetFiducialEtaMax();  
+
+  //loop over tracks
+  for (Int_t it = 0; it < nt; it++) {
+      AliESDtrack *track = fESD->GetTrack(it);
+      UInt_t status = track->GetStatus();
+
+      Double_t mom[3];
+      track->GetPxPyPz(mom);
+      p3.SetXYZ(mom[0],mom[1],mom[2]);
+      pt = p3.Pt();
+      if (((status & AliESDtrack::kITSrefit) == 0) ||
+         ((status & AliESDtrack::kTPCrefit) == 0)) continue;    // quality check
+      if (((AliJetESDReaderHeader*) fReaderHeader)->ReadSignalOnly() 
+         && TMath::Abs(track->GetLabel()) > 10000)  continue;   // quality check
+      if (((AliJetESDReaderHeader*) fReaderHeader)->ReadBkgdOnly() 
+         && TMath::Abs(track->GetLabel()) < 10000)  continue;   // quality check
+      eta = p3.Eta();
+      if ( (eta > etaMax) || (eta < etaMin)) continue;           // checking eta cut
+      
+      new ((*fMomentumArray)[goodTrack]) TLorentzVector(p3,p3.Mag());
+      
+      sflag[goodTrack]=0;
+      if (TMath::Abs(track->GetLabel()) < 10000) sflag[goodTrack]=1;
+      cflag[goodTrack]=0;
+      if (pt > ptMin) cflag[goodTrack]=1;                       // pt cut
+      goodTrack++;
+
+      // Monte Carlo information
+      Int_t label = TMath::Abs(track->GetLabel());
+      TClonesArray &arrayMC = *fArrayMC;
+      TParticle *part = (TParticle*)arrayMC[label]; //particle
+      ptMC = part->Pt();
+      pdgCode = part->GetPdgCode();
+      energyMC = part->Energy();
+      pxMC = part->Px(); pyMC = part->Py(); pzMC = part->Pz();
+  }
+  // set the signal flags
+  fSignalFlag.Set(goodTrack,sflag);
+  fCutFlag.Set(goodTrack,cflag);
+
+}
+
+
diff --git a/JETAN/AliJetESDmcReader.h b/JETAN/AliJetESDmcReader.h
new file mode 100644 (file)
index 0000000..5fc1f49
--- /dev/null
@@ -0,0 +1,36 @@
+#ifndef ALIJETESDMCREADER_H
+#define ALIJETESDMCREADER_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+// Jet ESD Reader 
+// ESD reader for jet analysis (it reads the esd and the MC trees)
+// Author: Mercedes Lopez Noriega (mercedes.lopez.noriega@cern.ch)
+
+#include "AliJetReader.h"
+class AliJetESDReaderHeader;
+
+
+class AliJetESDmcReader : public AliJetReader
+{
+ public: 
+  AliJetESDmcReader();
+  virtual ~AliJetESDmcReader();
+
+  // Getters
+  Float_t GetTrackMass() const {return fMass;}  // returns mass of the track
+  Int_t   GetTrackSign() const {return fSign;}  // returns sign of the track
+
+  // Setters
+  void FillMomentumArray(Int_t event); 
+  void OpenInputFiles();
+   
+ protected:
+  Float_t fMass;    // Particle mass
+  Int_t   fSign;    // Particle sign
+
+  ClassDef(AliJetESDmcReader,1)
+};
+#endif
index 7466b42..e323233 100755 (executable)
@@ -4,10 +4,11 @@
 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  * See cxx source for full Copyright notice                               */
  
-//-------------------------------------------------------------------------
+//---------------------------------------------------------------------
 // base class for Jet Reader Header 
+//
 // Author: jgcn@mda.cinvestav.mx
-//-------------------------------------------------------------------------
+//---------------------------------------------------------------------
   
 #include <TNamed.h>
 #include <TString.h>
@@ -21,22 +22,24 @@ class AliJetReaderHeader : public TNamed
   virtual ~AliJetReaderHeader();
   
   // Getters
-  virtual TString GetComment() const {return fComment;}
+  virtual const TString GetComment() {return fComment;}
   virtual const char* GetDirectory() {return fDir.Data();}
   virtual const char* GetPattern() {return fPattern.Data();}
   virtual Float_t     GetFiducialEtaMin() const {return fFiducialEtaMin;}
   virtual Float_t     GetFiducialEtaMax() const {return fFiducialEtaMax;}  
   virtual Float_t GetPtCut()       const  {return fPtCut;}
-  Float_t GetDCA()         const  {return fDCA;}       // not working so far..(always 0)
-  Float_t GetTLength()     const  {return fTLength;}   // not working so far.. (always 0)
+  Float_t GetDCA() const  {return fDCA;}       // not working so far..(always 0)
+  Float_t GetTLength() const  {return fTLength;}   // not working so far.. (always 0)
+  Int_t   GetNesd() const {return fNesd;}
   Int_t   GetNEvents() const {return fLast-fFirst;}
-  Int_t   GetLastEvent() const {return fLast;}
   Int_t   GetFirstEvent() const {return fFirst;}
+  Int_t   GetLastEvent() const {return fLast;}
 
   // Setters
   virtual void SetComment(const char* s) {fComment=TString(s);}
   virtual void SetPattern(const char* s) {fPattern=TString(s);}
   virtual void SetDirectory(const char* s) {fDir=TString(s);}
+  virtual void SetNumberOfESD(Int_t i=1) {fNesd = i;}
   virtual void SetFirstEvent(Int_t i=0) {fFirst=i;}
   virtual void SetLastEvent(Int_t i=-1) {fLast=i;}
   virtual void SetFiducialEta(Float_t etamin, Float_t etamax) 
@@ -47,6 +50,7 @@ class AliJetReaderHeader : public TNamed
  
  protected:
 
+  Int_t fNesd;             // Number of ESDs to read
   Int_t fFirst;            // First and last events analyzed
   Int_t fLast;             // in current set of files
   Float_t fFiducialEtaMin; // Fiducial minimum eta
index ad26b6c..c1d6d62 100644 (file)
@@ -15,6 +15,7 @@
 #pragma link C++ class AliJetESDReaderHeader+;
 #pragma link C++ class AliJetReader+;
 #pragma link C++ class AliJetESDReader+;
+#pragma link C++ class AliJetESDmcReader+;
 #pragma link C++ class AliJetKineReader+;
 #pragma link C++ class AliJetKineReaderHeader+;
 #pragma link C++ class AliJetControlPlots+;
index befd8f3..07cb99a 100644 (file)
@@ -3,8 +3,8 @@
 SRCS =         AliJet.cxx AliJetHeader.cxx AliPxconeJetHeader.cxx \
        AliJetFinder.cxx AliPxconeJetFinder.cxx AliJetReaderHeader.cxx \
        AliJetESDReaderHeader.cxx AliJetReader.cxx AliJetESDReader.cxx \
-       AliJetControlPlots.cxx AliUA1JetHeader.cxx AliUA1JetFinder.cxx \
-       AliLeading.cxx \
+       AliJetESDmcReader.cxx AliJetControlPlots.cxx AliUA1JetHeader.cxx \
+       AliUA1JetFinder.cxx AliLeading.cxx \
        AliJetKineReader.cxx AliJetKineReaderHeader.cxx \
        AliJetProductionData.cxx AliJetProductionDataPDC2004.cxx \
        AliJetAnalysis.cxx AliJetDistributions.cxx \