Adding HLT tracks to ESD
authorhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 11 May 2004 11:58:55 +0000 (11:58 +0000)
committerhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 11 May 2004 11:58:55 +0000 (11:58 +0000)
STEER/AliESD.cxx
STEER/AliESD.h
STEER/AliESDHLTtrack.cxx [new file with mode: 0644]
STEER/AliESDHLTtrack.h [new file with mode: 0644]
STEER/AliReconstruction.cxx
STEER/AliReconstruction.h
STEER/STEERLinkDef.h
STEER/libSTEER.pkg

index fffee719f52c961c43516a6a1b71aba9457c72c8..a00b62f88f72350d9d826ea1f36e0aa5c746c4a6 100644 (file)
@@ -36,6 +36,8 @@ AliESD::AliESD():
   fT0zVertex(0),
   fPrimaryVertex(),
   fTracks("AliESDtrack",15000),
+  fHLTConfMapTracks("AliESDHLTtrack",25000),
+  fHLTHoughTracks("AliESDHLTtrack",15000),
   fMuonTracks("AliESDMuonTrack",30),
   fPmdTracks("AliESDPmdTrack",3000),
   fV0s("AliESDv0",200),
@@ -79,6 +81,8 @@ void AliESD::Print(Option_t *) const
   printf("Event from reconstruction version %d \n",fRecoVersion);
   printf("Number of tracks: \n");
   printf("                 charged   %d\n",GetNumberOfTracks()-GetNumberOfPHOSParticles()-GetNumberOfEMCALParticles());
+  printf("                 hlt CF    %d\n", GetNumberOfHLTConfMapTracks());
+  printf("                 hlt HT    %d\n", GetNumberOfHLTHoughTracks());
   printf("                 phos      %d\n", GetNumberOfPHOSParticles());
   printf("                 emcal     %d\n", GetNumberOfEMCALParticles());
   printf("                 muon      %d\n", GetNumberOfMuonTracks());
index c80dda9be9d553ff8dc863d61421a2528982a88e..4693582a4cf91bb494bdc37ed516801f72d8e11a 100644 (file)
@@ -21,6 +21,7 @@
 #include "AliESDVertex.h"
 #include "AliESDcascade.h"
 #include "AliESDtrack.h"
+#include "AliESDHLTtrack.h"
 #include "AliESDv0.h"
 
 class AliESD : public TObject {
@@ -37,6 +38,12 @@ public:
   AliESDtrack *GetTrack(Int_t i) const {
     return (AliESDtrack *)fTracks.UncheckedAt(i);
   }
+  AliESDHLTtrack *GetHLTConfMapTrack(Int_t i) const {
+    return (AliESDHLTtrack *)fHLTConfMapTracks.UncheckedAt(i);
+  }
+  AliESDHLTtrack *GetHLTHoughTrack(Int_t i) const {
+    return (AliESDHLTtrack *)fHLTHoughTracks.UncheckedAt(i);
+  }
   AliESDMuonTrack *GetMuonTrack(Int_t i) const {
     return (AliESDMuonTrack *)fMuonTracks.UncheckedAt(i);
   }
@@ -47,7 +54,12 @@ public:
   void AddTrack(const AliESDtrack *t) {
     new(fTracks[fTracks.GetEntriesFast()]) AliESDtrack(*t);
   }
-
+  void AddHLTConfMapTrack(const AliESDHLTtrack *t) {
+    new(fHLTConfMapTracks[fHLTConfMapTracks.GetEntriesFast()]) AliESDHLTtrack(*t);
+  }
+  void AddHLTHoughTrack(const AliESDHLTtrack *t) {
+    new(fHLTHoughTracks[fHLTHoughTracks.GetEntriesFast()]) AliESDHLTtrack(*t);
+  }
   void AddMuonTrack(const AliESDMuonTrack *t) {
     new(fMuonTracks[fMuonTracks.GetEntriesFast()]) AliESDMuonTrack(*t);
   }
@@ -79,6 +91,8 @@ public:
   Long_t GetTrigger() const {return fTrigger;}
   
   Int_t GetNumberOfTracks()     const {return fTracks.GetEntriesFast();}
+  Int_t GetNumberOfHLTConfMapTracks()     const {return fHLTConfMapTracks.GetEntriesFast();}
+  Int_t GetNumberOfHLTHoughTracks()     const {return fHLTHoughTracks.GetEntriesFast();}
   Int_t GetNumberOfMuonTracks() const {return fMuonTracks.GetEntriesFast();}
   Int_t GetNumberOfPmdTracks() const {return fPmdTracks.GetEntriesFast();}
   Int_t GetNumberOfV0s()      const {return fV0s.GetEntriesFast();}
@@ -109,6 +123,8 @@ protected:
   AliESDVertex  fPrimaryVertex;  // Primary vertex estimated by the ITS
 
   TClonesArray  fTracks;         // ESD tracks
+  TClonesArray  fHLTConfMapTracks; // HLT ESD tracks from Conformal Mapper method
+  TClonesArray  fHLTHoughTracks; // HLT ESD tracks from Hough Transform method
   TClonesArray  fMuonTracks;     // MUON ESD tracks
   TClonesArray  fPmdTracks;      // PMD ESD tracks
   TClonesArray  fV0s;            // V0 vertices
@@ -118,8 +134,8 @@ protected:
   Int_t         fFirstPHOSParticle; // First PHOS particle in the fTracks list 
   Int_t         fFirstEMCALParticle;// First EMCAL particle in the fTracks list 
  
-  ClassDef(AliESD,6)  //ESD class 
-                      //ver. 2: Magnetic Field Added; skowron
+  ClassDef(AliESD,7)  //ESD class 
+
 };
 
 #endif 
diff --git a/STEER/AliESDHLTtrack.cxx b/STEER/AliESDHLTtrack.cxx
new file mode 100644 (file)
index 0000000..14a9d47
--- /dev/null
@@ -0,0 +1,61 @@
+/**************************************************************************
+ * 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.                  *
+ **************************************************************************/
+
+//-----------------------------------------------------------------
+//           Implementation of the ESD HLT track class
+//   ESD = Event Summary Data
+//   HLT = High Level Trigger
+//   This is the class to deal with during the phisical analysis of data
+//-----------------------------------------------------------------
+
+#include "TMath.h"
+#include "AliESDHLTtrack.h"
+
+ClassImp(AliESDHLTtrack)
+
+AliESDHLTtrack::AliESDHLTtrack() : TObject()
+{
+  fNHits = 0;
+  fMCid = 0;
+  fWeight = 0;
+  fFromMainVertex = kFALSE;
+  fRowRange[0] = fRowRange[1] = 0;
+  fSector = 0;
+  fFirstPoint[0] = fFirstPoint[1] = fFirstPoint[2] = 0;
+  fLastPoint[0] = fLastPoint[1] = fLastPoint[2] = 0;
+  fQ = 0;
+  fTanl = 0;
+  fPsi = 0;
+  fPt = 0;
+  fPterr = 0;
+  fPsierr = 0;
+  fTanlerr = 0;
+  fBinX = 0;
+  fBinY = 0;
+  fSizeX = 0;
+  fSizeY = 0;
+  fPID =0;
+}
+
+Double_t AliESDHLTtrack::GetP() const
+{
+  // Returns total momentum.  
+  return TMath::Abs(GetPt())*sqrt(1. + GetTgl()*GetTgl());
+}
+
+Double_t AliESDHLTtrack::GetPseudoRapidity() const
+{
+  return 0.5 * TMath::Log((GetP() + GetPz()) / (GetP() - GetPz()));
+}
diff --git a/STEER/AliESDHLTtrack.h b/STEER/AliESDHLTtrack.h
new file mode 100644 (file)
index 0000000..7e73495
--- /dev/null
@@ -0,0 +1,121 @@
+#ifndef ALIESDHLTTRACK_H
+#define ALIESDHLTTRACK_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+//-------------------------------------------------------------------------
+//                          Class AliESDHLTtrack
+//   This is the class to handle HLT reconstruted TPC tracks
+//-------------------------------------------------------------------------
+#include "TObject.h"
+
+class AliESDHLTtrack : public TObject {
+public:
+  AliESDHLTtrack();
+  virtual ~AliESDHLTtrack() {}
+
+  // getters
+  Int_t GetNHits() const {return fNHits;}
+
+  Int_t GetMCid() const {return fMCid;}
+
+  Int_t GetWeight() const {return fWeight;}
+
+  Bool_t ComesFromMainVertex() const {return fFromMainVertex;}
+
+  Int_t GetFirstRow() const {return fRowRange[0];}
+  Int_t GetLastRow()  const {return fRowRange[1];}
+  Int_t GetSector()   const {return fSector;}
+
+  Double_t GetFirstPointX() {return fFirstPoint[0];}
+  Double_t GetFirstPointY() {return fFirstPoint[1];}
+  Double_t GetFirstPointZ() {return fFirstPoint[2];}
+  Double_t GetLastPointX() {return fLastPoint[0];}
+  Double_t GetLastPointY() {return fLastPoint[1];}
+  Double_t GetLastPointZ() {return fLastPoint[2];}
+
+  Int_t GetCharge() const {return fQ;}
+  Double_t GetPt() const {return fPt;}
+  Double_t GetTgl() const {return fTanl;}
+  Double_t GetPsi() const {return fPsi;}
+
+  Double_t GetPterr()  const {return fPterr;}
+  Double_t GetPsierr() const {return fPsierr;}
+  Double_t GetTglerr() const {return fTanlerr;}
+
+  Float_t    GetBinX()   const {return fBinX;}
+  Float_t    GetBinY()   const {return fBinY;}
+  Float_t    GetSizeX()  const {return fSizeX;}
+  Float_t    GetSizeY()  const {return fSizeY;}
+
+  Double_t GetPx() const {return fPt*cos(fPsi);}
+  Double_t GetPy() const {return fPt*sin(fPsi);}
+  Double_t GetPz() const {return fPt*fTanl;}
+
+  Double_t GetP() const;
+  Double_t GetPseudoRapidity() const;
+
+  Float_t GetPID() {return fPID;}
+
+  // setters
+  void SetNHits(Int_t f) {fNHits = f;}
+
+  void SetMCid(Int_t f) {fMCid = f;}
+
+  void SetWeight(Int_t f) {fWeight = f;}
+  
+  void ComesFromMainVertex(Bool_t f) {fFromMainVertex = f;}
+  
+  void SetRowRange(Int_t f,Int_t g) {fRowRange[0]=f; fRowRange[1]=g;}
+  void SetSector(Int_t f) {fSector = f;}
+
+  void SetFirstPoint(Double_t f,Double_t g,Double_t h) {fFirstPoint[0]=f; fFirstPoint[1]=g; fFirstPoint[2]=h;}
+  void SetLastPoint(Double_t f,Double_t g,Double_t h) {fLastPoint[0]=f; fLastPoint[1]=g; fLastPoint[2]=h;}
+
+  void SetCharge(Int_t f) {fQ = f;}
+  void SetTgl(Double_t f) {fTanl =f;}
+  void SetPsi(Double_t f) {fPsi = f;}
+  void SetPt(Double_t f) {fPt = f;}
+
+  void SetPterr(Double_t f) {fPterr = f;}
+  void SetPsierr(Double_t f) {fPsierr = f;}
+  void SetTglerr(Double_t f) {fTanlerr = f;}
+
+  void SetBinXY(Float_t binx,Float_t biny,Float_t sizex,Float_t sizey) {fBinX = binx; fBinY = biny; fSizeX = sizex; fSizeY = sizey;}
+
+  void SetPID(Float_t pid) {fPID = pid;}
+
+protected:
+  Int_t fNHits;  // Number of assigned clusters
+
+  Int_t fMCid;  //Assigned id from MC data.
+
+  Int_t fWeight; //Weight associated to Hough Transform
+
+  Bool_t   fFromMainVertex; // true if tracks origin is the main vertex, otherwise false
+  
+  Int_t fRowRange[2]; //Subsector where this track was build
+  Int_t fSector;      //Sector # where  this track was build
+
+  Double_t fFirstPoint[3]; //First and last track point in TPC
+  Double_t fLastPoint[3];
+
+  Int_t    fQ;    //track charge
+  Double_t fTanl; //tan of dipangle
+  Double_t fPsi;  //azimuthal angle of the momentum 
+  Double_t fPt;   //transverse momentum
+
+  Double_t fPterr;
+  Double_t fPsierr;
+  Double_t fTanlerr;
+
+  Float_t fBinX;
+  Float_t fBinY;
+  Float_t fSizeX;
+  Float_t fSizeY;
+  
+  Float_t fPID; //so far filled only for conformal mapper tracks
+
+  ClassDef(AliESDHLTtrack,1) //ESD HLT track class
+};
+
+#endif
index 7d45f9a7f0ad40ce2d36877538977200fbb6a58a..a1635a2e3ab2c5e431d31a48875585f0eb0da617 100644 (file)
 #include "AliGenEventHeader.h"
 #include "AliESDpid.h"
 #include "AliMagF.h"
+#include "../TPC/AliTPCReconstructor.h"
 
 ClassImp(AliReconstruction)
 
 
 //_____________________________________________________________________________
-const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "RICH", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "START", "VZERO", "CRT"};
+const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "RICH", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "START", "VZERO", "CRT", "HLT"};
 
 //_____________________________________________________________________________
 AliReconstruction::AliReconstruction(const char* gAliceFilename,
@@ -204,7 +205,16 @@ Bool_t AliReconstruction::Run()
   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
     TString detName = fgkDetectorName[iDet];
     TString recName = "Ali" + detName + "Reconstructor";
-    if (!gAlice->GetDetector(detName)) continue;
+    if (!gAlice->GetDetector(detName) && detName != "HLT") continue;
+
+    if(detName == "HLT") {
+      if (!gROOT->GetClass("AliLevel3")) {
+       gSystem->Load("libAliL3Src.so");
+       gSystem->Load("libAliL3Misc.so");
+       gSystem->Load("libAliL3Hough.so");
+       gSystem->Load("libAliL3Comp.so");
+      }
+    }
 
     AliReconstructor* reconstructor = NULL;
     // first check if a plugin is defined for the reconstructor
index 98b6aabb20a0946ee2d3fe2d1998d14dc08b61dc..59374566549d6d08317f162a120ef57a39512d05 100644 (file)
@@ -114,7 +114,7 @@ private:
   AliLoader*     fTOFLoader;          //! loader for TOF
   AliTracker*    fTOFTracker;         //! tracker for TOF
 
-  static const Int_t fgkNDetectors = 14;   //! number of detectors
+  static const Int_t fgkNDetectors = 15;   //! number of detectors
   static const char* fgkDetectorName[fgkNDetectors]; //! names of detectors
   TObjArray      fReconstructors;     //! array of reconstructor objects
 
index d165b5114042d7d61af172e9da1512abb0ba1597..96b90b9ad6769c9c12dccaa8572d9626c361328f 100644 (file)
@@ -68,6 +68,7 @@
 #pragma link C++ class  AliESDtrack+;
 #pragma link C++ class  AliESDMuonTrack+;
 #pragma link C++ class  AliESDPmdTrack+;
+#pragma link C++ class  AliESDHLTtrack+;
 #pragma link C++ class  AliReconstructor+;
 #pragma link C++ class  AliESDv0+;
 #pragma link C++ class  AliESDcascade+;
index b69f93b1b20643b970af0d914df864b44367b347..b57e6120af05ce3b4f0f908793b1902ba6c19abb 100644 (file)
@@ -16,7 +16,7 @@ AliMergeCombi.cxx AliMagFMaps.cxx AliFieldMap.cxx \
 AliGausCorr.cxx AliTrackReference.cxx AliESD.cxx \
 AliTrackMap.cxx AliTrackMapper.cxx AliCollisionGeometry.cxx \
 AliMemoryWatcher.cxx AliBarrelTrack.cxx \
-AliESDtrack.cxx AliESDMuonTrack.cxx AliESDPmdTrack.cxx AliESDv0.cxx AliESDcascade.cxx AliESDVertex.cxx AliESDpid.cxx \
+AliESDtrack.cxx AliESDMuonTrack.cxx AliESDPmdTrack.cxx AliESDHLTtrack.cxx AliESDv0.cxx AliESDcascade.cxx AliESDVertex.cxx AliESDpid.cxx \
 AliVertexer.cxx \
 AliMC.cxx AliSimulation.cxx AliReconstruction.cxx AliVertexGenFile.cxx \
 AliReconstructor.cxx