JETAN/AliJetParticlesReaderHLT.cxx

   1 // $Id$
   2
   3 //___________________________________________________________________________
   4 /////////////////////////////////////////////////////////////////////////////
   5 //                                                                         //
   6 // File reader for HLT tracks ESD                                          //
   7 //                                                                         //
   8 // loizides@ikf.uni-frankfurt.de                                           //
   9 /////////////////////////////////////////////////////////////////////////////
  10
  11 #include <Riostream.h>
  12 #include <TMath.h>
  13 #include <TString.h>
  14 #include <TObjString.h>
  15 #include <TTree.h>
  16 #include <TFile.h>
  17 #include <AliESD.h>
  18 #include <AliESDtrack.h>
  19 #include <AliESDHLTtrack.h>
  20 #include <AliL3Track.h>
  21 #include <AliL3Vertex.h>
  22 #include <AliKalmanTrack.h>
  23 #include <AliJetEventParticles.h>
  24 #include "AliJetParticlesReaderHLT.h"
  25
  26 ClassImp(AliJetParticlesReaderHLT)
  27
  28 AliJetParticlesReaderHLT::AliJetParticlesReaderHLT(Bool_t bMapper, const Char_t* esdfilename) :
  29   AliJetParticlesReaderESD(0,esdfilename),
  30   fTrackerType(bMapper),
  31   fMinHits(0),
  32   fMinWeight(0)
  33 {
  34   //constructor
  35 }
  36
  37 /********************************************************************/
  38
  39 AliJetParticlesReaderHLT::AliJetParticlesReaderHLT(
  40                                       Bool_t bMapper,
  41                                       TObjArray* dirs,
  42                                       const Char_t* esdfilename) :
  43   AliJetParticlesReaderESD(0,dirs,esdfilename),
  44   fTrackerType(bMapper),
  45   fMinHits(0),
  46   fMinWeight(0)
  47
  48 {
  49   //constructor
  50 }
  51
  52 /********************************************************************/
  53
  54 AliJetParticlesReaderHLT::~AliJetParticlesReaderHLT()
  55 {
  56   //desctructor
  57 }
  58
  59 Int_t AliJetParticlesReaderHLT::ReadESD(AliESD* esd)
  60 {
  61   //Reads one ESD
  62
  63   if (esd == 0)
  64    {
  65      Error("ReadESD","ESD is NULL");
  66      return kFALSE;
  67    }
  68
  69   Float_t mf = esd->GetMagneticField();
  70   if (mf <= 0.0)
  71   {
  72      Error("ReadESD","Magnetic Field is 0. Skipping to next event.");
  73      return kFALSE;
  74   }
  75   AliKalmanTrack::SetMagneticField(mf/10.);
  76
  77   Info("ReadESD","Reading Event %d",fCurrentDir*1000+fCurrentEvent);
  78   if((!fOwner) || (fEventParticles==0))
  79     fEventParticles = new AliJetEventParticles();
  80
  81   Int_t ntr=0;
  82   if(fTrackerType){
  83     ntr =esd->GetNumberOfHLTHoughTracks();
  84     Info("ReadESD","Found %d conformal tracks.",ntr);
  85   } else {
  86     ntr=esd->GetNumberOfHLTHoughTracks();
  87     Info("ReadESD","Found %d hough tracks.",ntr);
  88   }
  89   fEventParticles->Reset(ntr);
  90
  91   TString headdesc="";
  92   headdesc+="Run ";
  93   headdesc+=esd->GetRunNumber();
  94   headdesc+=": Ev ";
  95   headdesc+=esd->GetEventNumber();
  96   fEventParticles->SetHeader(headdesc);
  97
  98   Double_t vertexpos[3];//vertex position
  99   const AliESDVertex* kvertex = esd->GetVertex();
 100   if (kvertex == 0)
 101    {
 102      Info("ReadESD","ESD returned NULL pointer to vertex - assuming (0.0,0.0,0.0)");
 103      vertexpos[0] = 0.0;
 104      vertexpos[1] = 0.0;
 105      vertexpos[2] = 0.0;
 106    }
 107   else
 108    {
 109      kvertex->GetXYZ(vertexpos);
 110    }
 111
 112   fEventParticles->SetVertex(vertexpos[0],vertexpos[1],vertexpos[2]);
 113   //cout << vertexpos[0] << " " << vertexpos[1] << " " << vertexpos[2] << endl;
 114
 115   AliL3Track l3;
 116   AliL3Vertex v;
 117   v.SetX(vertexpos[0]);
 118   v.SetY(vertexpos[1]);
 119   v.SetZ(vertexpos[2]);
 120   Double_t xc=0.,yc=0.,zc=0.;
 121   for (Int_t i = 0;i<ntr; i++) {
 122     AliESDHLTtrack *kesdtrack;
 123     if(fTrackerType){
 124       kesdtrack=esd->GetHLTConfMapTrack(i);
 125     } else {
 126       kesdtrack=esd->GetHLTHoughTrack(i);
 127     }
 128
 129     if (kesdtrack == 0)
 130       {
 131         Error("ReadESD","Can not get track %d", i);
 132         continue;
 133       }
 134
 135     //const Float_t kpid=kesdtrack->GetPID();
 136     const Int_t knhits=kesdtrack->GetNHits();
 137     const Int_t kweight=kesdtrack->GetWeight();
 138     //cout << i << " " << kweight << " " << knhits << endl;
 139     if((fMinHits>0) && (knhits<fMinHits)) continue;
 140     if(kweight>1000) continue; //avoid ghosts
 141     if((fMinWeight>0) && (kweight<fMinWeight)) continue;
 142
 143     Float_t px=kesdtrack->GetPx();
 144     Float_t py=kesdtrack->GetPy();
 145     Float_t pz=kesdtrack->GetPz();
 146
 147     if(fTrackerType){
 148       //if(!kesdtrack->ComesFromMainVertex()) continue;
 149       //cout << kesdtrack->GetPx() << " " << kesdtrack->GetPy() << " " << kesdtrack->GetPz() << endl;
 150       l3.SetFirstPoint(kesdtrack->GetFirstPointX(),kesdtrack->GetFirstPointY(),kesdtrack->GetFirstPointZ());
 151       l3.SetLastPoint(kesdtrack->GetLastPointX(),kesdtrack->GetLastPointY(),kesdtrack->GetLastPointZ());
 152       l3.SetCharge(kesdtrack->GetCharge());
 153       l3.SetPt(kesdtrack->GetPt());
 154       l3.SetTgl(kesdtrack->GetTgl());
 155       l3.SetPsi(kesdtrack->GetPsi());
 156       l3.CalculateHelix();
 157       l3.GetClosestPoint(&v,xc,yc,zc);
 158       if(TMath::Abs(zc)>10.) continue;
 159       l3.SetFirstPoint(vertexpos[0],vertexpos[1],vertexpos[2]);
 160       //l3.CalculateHelix();
 161       l3.UpdateToFirstPoint();
 162       px=l3.GetPx();
 163       py=l3.GetPy();
 164       pz=l3.GetPz();
 165       //cout << px << " " << py << " " << pz << endl;
 166     }
 167
 168     const Float_t kpt=TMath::Sqrt(px*px+py*py);
 169     const Float_t kp=TMath::Sqrt(pz*pz+kpt*kpt);
 170     const Float_t keta=0.5*TMath::Log((kp+pz+1e-30)/(kp-pz+1e-30));
 171     const Float_t kphi=TMath::Pi()+TMath::ATan2(-py,-px);
 172
 173     if(IsAcceptedParticle(kpt,kphi,keta))
 174       fEventParticles->AddParticle(px,py,pz,kp,i,kesdtrack->GetMCid(),knhits,kpt,kphi,keta);
 175   } //loop over tracks
 176
 177   return kTRUE;
 178 }