PHOS/AliPHOSJet.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
  16 /* $Id$ */
  17
  18 //_________________________________________________________________________
  19 // Class to calculate jet chararacteristics
  20 //
  21 // This class implements for PHOS a jet finder for PHOS. It depends on a
  22 // energy seed
  23 // minimum energy, cone radius and movement of the cone.
  24 //*-- Author :  D.Peressounko
  25 //////////////////////////////////////////////////////////////////////////////
  26
  27 // --- ROOT system ---
  28 #include "TParticle.h"
  29
  30 // --- Standard library ---
  31
  32 // --- AliRoot header files ---
  33 #include "AliPHOSJet.h"
  34 //  #include "AliPHOSGetter.h"
  35
  36 ClassImp(AliPHOSJet)
  37
  38 //____________________________________________________________________________
  39 AliPHOSJet::AliPHOSJet():
  40   fNpart(0),
  41   fList(0),
  42   fConeRad(0),
  43   fMaxConeMove(0),
  44   fMinConeMove(0),
  45   fSumEnergy(0),
  46   fSumEta(0),
  47   fSumPhi(0),
  48   fEnergy(0),
  49   fEta(0),
  50   fPhi(0),
  51   fLEnergy(0),
  52   fLEta(0),
  53   fLPhi(0)
  54 {
  55   //Initialize members
  56 }
  57
  58 //____________________________________________________________________________
  59 AliPHOSJet::AliPHOSJet(const AliPHOSJet & jet) :
  60   TObject(jet),
  61   fNpart(0),
  62   fList(0),
  63   fConeRad(0),
  64   fMaxConeMove(0),
  65   fMinConeMove(0),
  66   fSumEnergy(0),
  67   fSumEta(0),
  68   fSumPhi(0),
  69   fEnergy(0),
  70   fEta(0),
  71   fPhi(0),
  72   fLEnergy(0),
  73   fLEta(0),
  74   fLPhi(0)
  75 {
  76   // copy ctor: no implementation yet
  77   Fatal("cpy ctor", "not implemented") ;
  78 }
  79
  80
  81 //____________________________________________________________________________
  82 AliPHOSJet::~AliPHOSJet(){
  83   //dtor
  84   if(fList){
  85     delete fList ;
  86     fList = 0 ;
  87   }
  88
  89 }
  90 //____________________________________________________________________________
  91 void AliPHOSJet::AddParticle(const TParticle * p, Int_t index){
  92   //adds particle to jet. Calculates change in jet direction,
  93   //due to addition of this particle and if it is smaller, than fMaxDev,
  94   //add particle, axcept new direction and return true.
  95
  96   fSumEnergy+=p->Energy() ;
  97   if(p->Pt()/p->Energy() > 0.001)
  98     fSumEta+=p->Eta()*p->Energy() ;
  99   else
 100     fSumEta+=100.*p->Pz() ;
 101   fSumPhi+=p->Phi()*p->Energy() ;
 102
 103   //check if this a leading particle?
 104   if(fLEnergy < p->Energy()){
 105     fLEnergy = p->Energy() ;
 106     if(p->Pt()/p->Energy() > 0.001)
 107       fLEta = p->Eta() ;
 108     else
 109       fLEta = 100.*p->Pz()/p->Energy() ;
 110     fLPhi = p->Phi() ;
 111   }
 112
 113   if(index >=0){ //add index to list of indexes
 114     if(!fList){
 115       fList = new TArrayI(100) ;
 116     }
 117     if(fList->GetSize()<=fNpart+1){
 118       TArrayI * tmp = new TArrayI(fNpart*2) ;
 119       tmp->Adopt(fList->GetSize(),fList->GetArray()) ; //note, we should not delete old array!
 120       fList=tmp ;
 121     }
 122     fList->AddAt(index,fNpart++) ;
 123   }
 124 }
 125 //____________________________________________________________________________
 126 void AliPHOSJet::AddDigit(Double_t e, Double_t eta, Double_t phi, Int_t index){
 127   //adds particle to jet. Calculates change in jet direction,
 128   //due to addition of this particle and if it is smaller, than fMaxDev,
 129   //add particle, axcept new direction and return true.
 130
 131   fSumEnergy+=e ;
 132   fSumEta+=eta*e ;
 133   fSumPhi+=phi*e ;
 134
 135   //check if this a leading particle?
 136   if(fLEnergy < e){
 137     fLEnergy = e;
 138     fLEta = eta ;
 139     fLPhi = phi ;
 140   }
 141
 142   if(index >=0){ //add index to list of indexes
 143     if(!fList){
 144       fList = new TArrayI(100) ;
 145     }
 146     if(fList->GetSize()<=fNpart+1){
 147       TArrayI * tmp = new TArrayI(fNpart*2) ;
 148       tmp->Adopt(fList->GetSize(),fList->GetArray()) ; //note, we should not delete old array!
 149       fList=tmp ;
 150     }
 151     fList->AddAt(index,fNpart++) ;
 152   }
 153 }
 154 // //____________________________________________________________________________
 155 // Bool_t AliPHOSJet::IsInJet(TParticle * p,Int_t mode)
 156 // {
 157 //   Double_t dEta ;
 158 //   Double_t dPhi ;
 159 //   Double_t energy ;
 160 //   if(!fEnergy){ //Final values not calculated yet, use intermediate
 161 //     if(fSumEnergy==0)
 162 //       return kTRUE ; //First particle
 163 // note p->Eta() causes fpe!
 164 //     dEta=(p->Eta() - fSumEta/fSumEnergy) ;
 165 //     dPhi=(p->Phi() - fSumPhi/fSumEnergy) ;
 166 //     energy = fSumEnergy ;
 167 //   }
 168 //   else{ //calculated with respect to final
 169 //     dEta=(p->Eta() - fEta) ;
 170 //     dPhi=(p->Phi() - fPhi) ;
 171 //     energy = fEnergy ;
 172 //   }
 173
 174 //   switch (mode){
 175 //   case 0:
 176 //     return IsInCone(eta,phi) ; //pure geometrical comparison
 177 //   case 1:
 178 //     return AcceptConeDeviation(dEta,dPhi,p->Energy() );
 179 //   default:
 180 //     AliError(Form("Unknown mode of cone calculation %d \n",mode ));
 181 //   }
 182 //   return kFALSE ;
 183 //}
 184 //____________________________________________________________________________
 185 Bool_t AliPHOSJet::AcceptConeDeviation(const TParticle * p)const
 186 { //Calculate cone deviation in case of inclusion of the given
 187   //particle to jet.
 188
 189   Double_t tmpEnergy = fSumEnergy + p->Energy() ;
 190   Double_t tmpEta = fSumEta ;
 191   if(p->Pt()/p->Energy() >0.001)
 192     tmpEta+= p->Eta()*p->Energy() ;
 193   else
 194     tmpEta+= 100.*p->Pz() ;
 195   tmpEta = tmpEta/tmpEnergy / (fSumEta/fSumEnergy) - 1 ;
 196   Double_t tmpPhi = fSumPhi + p->Phi()*p->Energy() ;
 197   tmpPhi = tmpPhi/tmpEnergy /(fSumPhi/fSumEnergy) - 1 ;
 198   Double_t dev = TMath::Sqrt(tmpEta*tmpEta + tmpPhi*tmpPhi) ;
 199   if(dev < fMaxConeMove)
 200     return kTRUE ;
 201   else
 202     return kFALSE ;
 203 }
 204 //____________________________________________________________________________
 205 Bool_t AliPHOSJet::AcceptConeDeviation(Double_t e, Double_t eta, Double_t phi)const
 206 { //Calculate cone deviation in case of inclusion of the given
 207   //particle to jet.
 208
 209   Double_t tmpEnergy = fSumEnergy + e ;
 210   Double_t tmpEta = fSumEta + eta*e ;
 211   tmpEta = tmpEta/tmpEnergy / (fSumEta/fSumEnergy) ;
 212   Double_t tmpPhi = fSumPhi + phi*e ;
 213   tmpPhi = tmpPhi/tmpEnergy /(fSumPhi/fSumEnergy) ;
 214   Double_t dev = TMath::Sqrt(tmpEta*tmpEta + tmpPhi*tmpPhi) ;
 215   if(dev<fMaxConeMove && dev > fMinConeMove)
 216     return kTRUE ;
 217   else
 218     return kFALSE ;
 219 }
 220 //____________________________________________________________________________
 221 Bool_t AliPHOSJet::IsInCone(const TParticle * p)const
 222 {
 223   //Say if  particle is inside the defined cone
 224   Double_t dEta ;
 225   Double_t dPhi ;
 226   if(!fEnergy){ //Final values not calculated yet, use intermediate
 227     if(fSumEnergy==0)
 228       return kTRUE ; //First particle
 229     if(p->Pt()/p->Energy() > 0.001)
 230       dEta=(p->Eta() - fSumEta/fSumEnergy) ;
 231     else
 232       dEta=(100.*p->Pz()/p->Energy() - fSumEta/fSumEnergy) ;
 233     dPhi=(p->Phi() - fSumPhi/fSumEnergy) ;
 234   }
 235   else{ //calculated with respect to final
 236     if(p->Pt()/p->Energy() > 0.001)
 237       dEta=(p->Eta() - fEta) ;
 238     else
 239       dEta=(100.*p->Pz()/p->Energy() - fEta) ;
 240     dPhi=(p->Phi() - fPhi) ;
 241   }
 242   if(TMath::Sqrt(dEta*dEta + dPhi*dPhi) < fConeRad)
 243     return kTRUE ;
 244   else
 245     return kFALSE ;
 246 }
 247 //____________________________________________________________________________
 248 Bool_t AliPHOSJet::IsInCone(Double_t eta, Double_t phi)const
 249 {
 250   //Says if particle is inside the defined cone
 251   Double_t dEta ;
 252   Double_t dPhi ;
 253   if(!fEnergy){ //Final values not calculated yet, use intermediate
 254     if(fSumEnergy==0)
 255       return kTRUE ; //First particle
 256     dEta=(eta - fSumEta/fSumEnergy) ;
 257     dPhi=(phi - fSumPhi/fSumEnergy) ;
 258   }
 259   else{ //calculated with respect to final
 260     dEta=(eta - fEta) ;
 261     dPhi=(phi - fPhi) ;
 262   }
 263   if(TMath::Sqrt(dEta*dEta + dPhi*dPhi) < fConeRad)
 264     return kTRUE ;
 265   else
 266     return kFALSE ;
 267 }
 268 //____________________________________________________________________________
 269 Double_t AliPHOSJet::DistanceToJet(const TParticle *p)const{
 270   //Calculate radius
 271   Double_t dEta ;
 272   Double_t dPhi ;
 273   if(!fEnergy){ //Final values not calculated yet, use intermediate
 274     if(fSumEnergy==0)
 275       return kTRUE ; //First particle
 276     if(p->Pt()/p->Energy() > 0.001)
 277       dEta=(p->Eta() - fSumEta/fSumEnergy) ;
 278     else
 279       dEta=(100.*p->Pz()/p->Energy() - fSumEta/fSumEnergy) ;
 280     dPhi=(p->Phi() - fSumPhi/fSumEnergy) ;
 281   }
 282   else{ //calculated with respect to final
 283     if(p->Pt()/p->Energy() > 0.001)
 284       dEta=(p->Eta() - fEta) ;
 285     else
 286       dEta=(100*p->Pz()/p->Energy() - fEta) ;
 287     dPhi=(p->Phi() - fPhi) ;
 288   }
 289   return TMath::Sqrt(dEta*dEta + dPhi*dPhi) ;
 290
 291 }
 292 //____________________________________________________________________________
 293 void AliPHOSJet::CalculateAll(void){
 294   //Calculate all jet parameters
 295   if(fSumEnergy==0)
 296     return  ; //Nothing to calculate
 297
 298   fEta = fSumEta/fSumEnergy ;
 299   fPhi = fSumPhi/fSumEnergy ;
 300   fEnergy = fSumEnergy ;
 301
 302   fSumEnergy = 0. ;
 303   fSumEta = 0. ;
 304   fSumPhi = 0. ;
 305 }
 306 //____________________________________________________________________________
 307 void AliPHOSJet::Print(const Option_t *) const {
 308   //Print jet parameters
 309   printf("-------------- AliPHOSJet ------------\n") ;
 310   printf(" Energy............. %f \n",fEnergy) ;
 311   printf(" Eta................ %f \n",fEta ) ;
 312   printf(" Phi................ %f \n",fPhi ) ;
 313   printf(" Leading Energy..... %f \n",fLEnergy) ;
 314   printf(" Leading Eta........ %f \n",fLEta) ;
 315   printf(" Leading Phi........ %f \n",fLPhi) ;
 316   printf(" N particles in jet  %d \n",fNpart) ;
 317   printf("----------------------------------\n") ;
 318 }
 319
 320
 321