PHOS/AliPHOSRecParticle.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15 /* $Id$ */
  16 //_________________________________________________________________________
  17 //  A Reconstructed Particle in PHOS
  18 //  To become a general class of AliRoot ?
  19 //  Why should I put meaningless comments
  20 //  just to satisfy
  21 //  the code checker
  22 //
  23 //*-- Author: Yves Schutz (SUBATECH)
  24
  25
  26 // --- ROOT system ---
  27
  28 // --- Standard library ---
  29
  30
  31 // --- AliRoot header files ---
  32 #include "AliPHOSRecParticle.h"
  33 #include "AliPHOSGetter.h"
  34 #include "AliPHOSGeometry.h"
  35
  36 //____________________________________________________________________________
  37   AliPHOSRecParticle::AliPHOSRecParticle(): fPHOSTrackSegment(0)  ,  fDebug( kFALSE )
  38 {
  39   // ctor
  40   const Int_t nSPECIES = AliPID::kSPECIESN;
  41   for(Int_t i = 0; i<nSPECIES ; i++)
  42     fPID[i]=0.;
  43 }
  44
  45
  46 //____________________________________________________________________________
  47  AliPHOSRecParticle::AliPHOSRecParticle(const AliPHOSRecParticle & rp)
  48    : AliPHOSFastRecParticle(rp)
  49 {
  50   // copy ctor
  51
  52   fPHOSTrackSegment = rp.fPHOSTrackSegment ;
  53   fDebug            = kFALSE ;
  54   fType             = rp.fType ;
  55   fIndexInList      = rp.fIndexInList ;
  56
  57   fPdgCode     = rp.fPdgCode;
  58   fStatusCode  = rp.fStatusCode;
  59   fMother[0]   = rp.fMother[0];
  60   fMother[1]   = rp.fMother[1];
  61   fDaughter[0] = rp.fDaughter[0];
  62   fDaughter[1] = rp.fDaughter[1];
  63   fWeight      = rp.fWeight;
  64   fCalcMass    = rp.fCalcMass;
  65   fPx          = rp.fPx;
  66   fPy          = rp.fPy;
  67   fPz          = rp.fPz;
  68   fE           = rp.fE;
  69   fVx          = rp.fVx;
  70   fVy          = rp.fVy;
  71   fVz          = rp.fVz;
  72   fVt          = rp.fVt;
  73   fPolarTheta  = rp.fPolarTheta;
  74   fPolarPhi    = rp.fPolarPhi;
  75   fParticlePDG = rp.fParticlePDG;
  76   const Int_t nSPECIES = AliPID::kSPECIESN;
  77   for(Int_t i = 0; i<nSPECIES ; i++)
  78     fPID[i]=rp.fPID[i];
  79 }
  80
  81 //____________________________________________________________________________
  82 Int_t AliPHOSRecParticle::GetNPrimaries() const
  83 {
  84   return -1;
  85 }
  86
  87 //____________________________________________________________________________
  88 Int_t AliPHOSRecParticle::GetNPrimariesToRecParticles() const
  89 {
  90   // Get the number of primaries at the origine of the RecParticle
  91   Int_t rv = 0 ;
  92   AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
  93   Int_t emcRPindex = dynamic_cast<AliPHOSTrackSegment*>(gime->TrackSegments()->At(GetPHOSTSIndex()))->GetEmcIndex();
  94   dynamic_cast<AliPHOSEmcRecPoint*>(gime->EmcRecPoints()->At(emcRPindex))->GetPrimaries(rv) ;
  95   return rv ;
  96 }
  97
  98 //____________________________________________________________________________
  99 const TParticle * AliPHOSRecParticle::GetPrimary() const
 100 {
 101   // Get the primary particle at the origine of the RecParticle and
 102   // which has deposited the largest energy in SDigits
 103   AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
 104   if (!gime)
 105     Error("GetPrimary", "Getter not yet instantiated") ;
 106   gime->Event(gime->EventNumber(), "SRTPX") ;
 107   if(GetNPrimaries() == 0)
 108     return 0 ;
 109   if(GetNPrimaries() == 1)
 110     return GetPrimary(0) ;
 111   Int_t AbsId = 0;
 112   Int_t module ;
 113   const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
 114    Double_t x,z ;
 115   geom->ImpactOnEmc(static_cast<double>(Theta()),static_cast<double>(Phi()), module,z,x);
 116   Int_t amp = 0 ;
 117   Int_t iPrim=-1 ;
 118   if(module != 0){
 119     geom->RelPosToAbsId(module,x,z,AbsId) ;
 120    TClonesArray * sdigits = gime->SDigits() ;
 121    AliPHOSDigit * sdig ;
 122
 123    for(Int_t i = 0 ; i < sdigits->GetEntriesFast() ; i++){
 124      sdig = static_cast<AliPHOSDigit *>(sdigits->At(i)) ;
 125      if((sdig->GetId() == AbsId)){
 126        if((sdig->GetAmp() > amp) && (sdig->GetNprimary())){
 127          amp = sdig->GetAmp() ;
 128          iPrim = sdig->GetPrimary(1) ;
 129        }
 130        // do not scan rest of list
 131        if(sdig->GetId() > AbsId)
 132          continue ;
 133      }
 134    }
 135   }
 136   if(iPrim >= 0)
 137     return gime->Primary(iPrim) ;
 138   else
 139     return 0 ;
 140 }
 141
 142 //____________________________________________________________________________
 143 const Int_t AliPHOSRecParticle::GetPrimaryIndex() const
 144 {
 145   // Get the primary track index in TreeK which deposits the most energy
 146   // in Digits which forms EmcRecPoint, which produces TrackSegment,
 147   // which the RecParticle is created from
 148
 149   AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
 150   if (!gime)
 151     Error("GetPrimary", "Getter not yet instantiated") ;
 152   gime->Event(gime->EventNumber(), "SRTPX") ;
 153
 154   // Get TrackSegment corresponding to this RecParticle
 155   AliPHOSTrackSegment *ts          = gime->TrackSegment(fPHOSTrackSegment);
 156
 157   // Get EmcRecPoint corresponding to this TrackSegment
 158   Int_t emcRecPointIndex = ts->GetEmcIndex();
 159   AliPHOSEmcRecPoint  *emcRecPoint = gime->EmcRecPoint(emcRecPointIndex);
 160
 161   // Get the list of digits forming this EmcRecParticle
 162   Int_t  nDigits   = emcRecPoint->GetDigitsMultiplicity();
 163   Int_t *digitList = emcRecPoint->GetDigitsList();
 164
 165   // Find the digit with maximum amplitude
 166   AliPHOSDigit *digit = 0;
 167   Int_t maxAmp = 0;
 168   Int_t bestDigitIndex = -1;
 169   for (Int_t iDigit=0; iDigit<nDigits; iDigit++) {
 170     digit = gime->Digit(digitList[nDigits]);
 171     if (digit->GetAmp() > maxAmp) {
 172       maxAmp = digit->GetAmp();
 173       bestDigitIndex = iDigit;
 174     }
 175   }
 176   digit = gime->Digit(digitList[bestDigitIndex]);
 177
 178   // Get the list of primary tracks producing this digit
 179   // and find which track has more track energy.
 180   Int_t nPrimary = digit->GetNprimary();
 181   TParticle *track = 0;
 182   Double_t energyEM     = 0;
 183   Double_t energyHadron = 0;
 184   Int_t    trackEM      = 0;
 185   Int_t    trackHadron  = 0;
 186   for (Int_t iPrim=0; iPrim<nPrimary; iPrim++) {
 187     Int_t iPrimary = digit->GetPrimary(iPrim);
 188     track = gime->Primary(iPrimary);
 189     Int_t pdgCode   = track->GetPdgCode();
 190     Double_t energy = track->Energy();
 191     if (pdgCode==22 || TMath::Abs(pdgCode)==11) {
 192       if (energy > energyEM) {
 193         energyEM = energy;
 194         trackEM = iPrimary;
 195       }
 196     }
 197     else {
 198       if (energy > energyHadron) {
 199         energyHadron = energy;
 200         trackHadron = iPrimary;
 201       }
 202     }
 203   }
 204   // Preferences are given to electromagnetic tracks
 205   if (trackEM     != 0) return trackEM;     // track is gamma or e+-
 206   if (trackHadron != 0) return trackHadron; // track is hadron
 207   return -12345;                              // no track found :(
 208 }
 209
 210 //____________________________________________________________________________
 211 const TParticle * AliPHOSRecParticle::GetPrimary(Int_t index) const
 212 {
 213   // Get one of the primary particles at the origine of the RecParticle
 214   if ( index > GetNPrimariesToRecParticles() ) {
 215     if (fDebug)
 216       Warning("GetPrimary", "AliPHOSRecParticle::GetPrimary -> %d is larger that the number of primaries %d",
 217               index, GetNPrimaries()) ;
 218     return 0 ;
 219   }
 220   else {
 221     Int_t dummy ;
 222     AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
 223
 224     Int_t emcRPindex = dynamic_cast<AliPHOSTrackSegment*>(gime->TrackSegments()->At(GetPHOSTSIndex()))->GetEmcIndex();
 225     Int_t primaryindex = dynamic_cast<AliPHOSEmcRecPoint*>(gime->EmcRecPoints()->At(emcRPindex))->GetPrimaries(dummy)[index] ;
 226     return gime->Primary(primaryindex) ;
 227    }
 228   //  return 0 ;
 229 }
 230
 231 //____________________________________________________________________________
 232 void AliPHOSRecParticle::SetPID(Int_t type, Double_t weight)
 233 {
 234   // Set the probability densities that this reconstructed particle
 235   // has a type of i:
 236   // i       particle types
 237   // ----------------------
 238   // 0       electron
 239   // 1       muon
 240   // 2       pi+-
 241   // 3       K+-
 242   // 4       p/pbar
 243   // 5       photon
 244   // 6       pi0 at high pt
 245   // 7       neutron
 246   // 8       K0L
 247   // 9       Conversion electron
 248
 249   fPID[type] = weight ;
 250 }