ANALYSIS/AliAODParticle.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 /////////////////////////////////////////////////////////////
  20 //
  21 // class AliAODParticle
  22 //
  23 // Ali HBT Particle: simplified class TParticle
  24 // Simplified in order to minimize the size of object
  25 //  - we want to keep a lot of such a objects in memory
  26 // Additionaly adjusted for HBT Analysies purposes
  27 // + pointer to Track Points
  28 // + pointer to Cluster Map(s)
  29 //
  30 // Piotr.Skowronski@cern.ch
  31 //
  32 /////////////////////////////////////////////////////////////
  33 #include <TParticle.h>
  34
  35 #include "AliAODParticle.h"
  36 #include "AliClusterMap.h"
  37 #include "AliLog.h"
  38 #include "AliTrackPoints.h"
  39
  40 ClassImp(AliAODParticle)
  41
  42 //______________________________________________________________________________
  43 AliAODParticle::AliAODParticle():
  44  fPdgIdx(0), fIdxInEvent(0),fNPids(0),fPids(0x0),fPidProb(0x0),
  45  fCalcMass(0),fPx(0), fPy(0),fPz(0),fE(0), fVx(0), fVy(0),fVz(0),fVt(0),
  46  fTPCTrackPoints(0x0),fITSTrackPoints(0x0),fClusterMap(0x0)
  47 {//empty particle
  48 }
  49 //______________________________________________________________________________
  50
  51 AliAODParticle::AliAODParticle(Int_t pdg, Int_t idx,
  52                Double_t px, Double_t py, Double_t pz, Double_t etot,
  53                Double_t vx, Double_t vy, Double_t vz, Double_t time):
  54   fPdgIdx(0), fIdxInEvent(idx),fNPids(0),fPids(0x0),fPidProb(0x0),
  55   fCalcMass(0),
  56   fPx(px), fPy(py),fPz(pz),fE(etot),
  57   fVx(vx), fVy(vy),fVz(vz),fVt(time),
  58   fTPCTrackPoints(0x0),fITSTrackPoints(0x0),fClusterMap(0x0)
  59 {
  60 //mormal constructor
  61   SetPdgCode(pdg);
  62   if (GetPDG()) {
  63      fCalcMass    = GetPDG()->Mass();
  64   } else {
  65      Double_t a2 = fE*fE -fPx*fPx -fPy*fPy -fPz*fPz;
  66      if (a2 >= 0) fCalcMass =  TMath::Sqrt(a2);
  67      else         fCalcMass = -TMath::Sqrt(-a2);
  68   }
  69 }
  70 //______________________________________________________________________________
  71
  72 AliAODParticle::AliAODParticle(Int_t pdg, Float_t prob, Int_t idx,
  73                                Double_t px, Double_t py, Double_t pz, Double_t etot,
  74                                Double_t vx, Double_t vy, Double_t vz, Double_t time):
  75   fPdgIdx(0), fIdxInEvent(idx),fNPids(0),fPids(0x0),fPidProb(0x0),
  76   fCalcMass(0),
  77   fPx(px), fPy(py),fPz(pz),fE(etot),
  78   fVx(vx), fVy(vy),fVz(vz),fVt(time),
  79   fTPCTrackPoints(0x0),fITSTrackPoints(0x0),fClusterMap(0x0)
  80 {
  81 //mormal constructor
  82   SetPdgCode(pdg,prob);
  83   if (GetPDG()) {
  84      fCalcMass    = GetPDG()->Mass();
  85   } else {
  86      Double_t a2 = fE*fE -fPx*fPx -fPy*fPy -fPz*fPz;
  87      if (a2 >= 0) fCalcMass =  TMath::Sqrt(a2);
  88      else         fCalcMass = -TMath::Sqrt(-a2);
  89   }
  90 }
  91 //______________________________________________________________________________
  92
  93 AliAODParticle::AliAODParticle(const AliAODParticle& in):
  94    AliVAODParticle(in),
  95    fPdgIdx(in.fPdgIdx), fIdxInEvent(in.fIdxInEvent),
  96    fNPids(in.fNPids),fPids(new Int_t[fNPids]),fPidProb(new Float_t[fNPids]),
  97    fCalcMass(in.GetCalcMass()),
  98    fPx(in.Px()),fPy(in.Py()),fPz(in.Pz()),fE(in.E()),
  99    fVx(in.Vx()),fVy(in.Vy()),fVz(in.Vz()),fVt(in.T()),
 100    fTPCTrackPoints(0x0),fITSTrackPoints(0x0),fClusterMap(0x0)
 101 {
 102  //Copy constructor
 103 // Info("AliAODParticle(const AliAODParticle& in)","");
 104  for(Int_t i = 0; i<fNPids; i++)
 105   {
 106     fPids[i] =  in.fPids[i];
 107     fPidProb[i] = in.fPidProb[i];
 108   }
 109
 110  if (in.fTPCTrackPoints)
 111    fTPCTrackPoints = (AliTrackPoints*)in.fTPCTrackPoints->Clone();
 112  if (in.fITSTrackPoints)
 113    fITSTrackPoints = (AliTrackPoints*)in.fITSTrackPoints->Clone();
 114  if (in.fClusterMap)
 115    fClusterMap = (AliClusterMap*)in.fClusterMap->Clone();
 116 }
 117 //______________________________________________________________________________
 118
 119 AliAODParticle::AliAODParticle(const AliVAODParticle& in):
 120    AliVAODParticle(in),
 121    fPdgIdx(0), fIdxInEvent(in.GetUID()),
 122    fNPids(0),fPids(0x0),fPidProb(0x0),
 123    fCalcMass(-1.0),
 124    fPx(in.Px()),fPy(in.Py()),fPz(in.Pz()),fE(in.E()),
 125    fVx(in.Vx()),fVy(in.Vy()),fVz(in.Vz()),fVt(in.T()),
 126    fTPCTrackPoints(0x0),fITSTrackPoints(0x0),fClusterMap(0x0)
 127 {
 128  //Copy constructor
 129 // Info("AliAODParticle(const AliVAODParticle& in)","");
 130  for(Int_t i = 0; i<in.GetNumberOfPids(); i++)
 131   {
 132     SetPIDprobability(in.GetNthPid(i),in.GetNthPidProb(i));
 133   }
 134  SetPdgCode(in.GetPdgCode(),in.GetPidProb());
 135
 136  AliTrackPoints* tpts = in.GetTPCTrackPoints();
 137  if (tpts)  SetTPCTrackPoints((AliTrackPoints*)tpts->Clone());
 138
 139  tpts = in.GetITSTrackPoints();
 140  if (tpts) SetITSTrackPoints((AliTrackPoints*)tpts->Clone());
 141
 142  AliClusterMap* clmap = in.GetClusterMap();
 143  if (clmap) SetClusterMap((AliClusterMap*)clmap->Clone());
 144 }
 145 //______________________________________________________________________________
 146
 147 AliAODParticle::AliAODParticle(const TParticle &p,Int_t idx):
 148    fPdgIdx(0), fIdxInEvent(idx),
 149    fNPids(0),fPids(0x0),fPidProb(0x0),
 150    fCalcMass(p.GetCalcMass()),
 151    fPx(p.Px()),fPy(p.Py()),fPz(p.Pz()),fE(p.Energy()),
 152    fVx(p.Vx()),fVy(p.Vy()),fVz(p.Vz()),fVt(p.T()),
 153    fTPCTrackPoints(0x0),fITSTrackPoints(0x0),fClusterMap(0x0)
 154 {
 155  //all copied in the initialization
 156  SetPdgCode(p.GetPdgCode());
 157 }
 158 //______________________________________________________________________________
 159
 160 AliAODParticle::~AliAODParticle()
 161 {
 162 //dtor
 163   delete [] fPids;
 164   delete [] fPidProb;
 165   delete fTPCTrackPoints;
 166   delete fITSTrackPoints;
 167   delete fClusterMap;
 168 }
 169 //______________________________________________________________________________
 170
 171 void AliAODParticle::Clear(Option_t*)
 172 {
 173 //Must be implemented in order to store this object in Clones Array
 174   delete [] fPids;
 175   delete [] fPidProb;
 176   delete fTPCTrackPoints;
 177   delete fITSTrackPoints;
 178   delete fClusterMap;
 179
 180   fPids = 0x0;
 181   fPidProb = 0x0;
 182   fTPCTrackPoints = 0x0;
 183   fITSTrackPoints = 0x0;
 184   fClusterMap = 0x0;
 185 }
 186 //______________________________________________________________________________
 187
 188 AliAODParticle& AliAODParticle::operator=(const AliAODParticle& in)
 189 {
 190 //assigment operator
 191 //  Info("operator=(const AliAODParticle& in)","AliAODParticle");
 192   if (&in != this) {
 193     AliVAODParticle::operator=(in);
 194     fNPids = in.fNPids;
 195     delete [] fPids;
 196     delete [] fPidProb;
 197     fPids = new Int_t[fNPids];
 198     fPidProb = new Float_t[fNPids];
 199     for (Int_t i = 0; i < fNPids;i++)
 200       {
 201         fPids[i]    = in.fPids[i];
 202         fPidProb[i] = in.fPidProb[i];
 203       }
 204
 205     fPdgIdx = in.fPdgIdx;
 206     fIdxInEvent = in.fIdxInEvent;
 207     fCalcMass = in.GetCalcMass();
 208     fPx = in.Px();
 209     fPy = in.Py();
 210     fPz = in.Pz();
 211     fE = in.E();
 212     fVx = in.Vx();
 213     fVy = in.Vy();
 214     fVz = in.Vz();
 215     fVt = in.T();
 216
 217     delete fTPCTrackPoints;
 218     fTPCTrackPoints = (in.fTPCTrackPoints)?(AliTrackPoints*)in.fTPCTrackPoints->Clone():0x0;
 219
 220     delete fITSTrackPoints;
 221     fITSTrackPoints = (in.fITSTrackPoints)?(AliTrackPoints*)in.fITSTrackPoints->Clone():0x0;
 222
 223     delete fClusterMap;
 224     fClusterMap =  (in.fClusterMap)?(AliClusterMap*)in.fClusterMap->Clone():0x0;
 225   }
 226   return *this;
 227 }
 228 //______________________________________________________________________________
 229
 230 void AliAODParticle::SetPdgCode(Int_t pdg,Float_t prob)
 231 {
 232 //Set PDG Code
 233   SetPIDprobability(pdg,prob);
 234   fPdgIdx = GetPidSlot(pdg);
 235 }
 236
 237 //______________________________________________________________________________
 238 void AliAODParticle::SetPIDprobability(Int_t pdg, Float_t prob)
 239 {
 240 //Sets another pdg code and corresponding probabilty
 241 //Ids are set in decreasing order
 242 //Check if total probability is not overcoming unity is performed
 243 //in case, warning is printed
 244   AliDebug(9,Form("SetPIDprobability","Setting PID %d prob %f",pdg,prob));
 245
 246   Float_t totprob = 0.0;//sums up probabilities
 247   Int_t idx = GetPidSlot(pdg);
 248   Int_t i;
 249
 250   if (idx > -1)
 251    {
 252      fPidProb[idx] = prob;
 253      for (i = 0; i < fNPids;i++) totprob+=fPidProb[i];
 254      if (totprob > (1.0+0.000001))
 255        {
 256          Warning("SetPIDprobability","Total probability greater than unity (%f)",totprob);
 257        }
 258      AliDebug(9,Form("Current Total probability: %f",totprob));
 259      return;
 260    }
 261
 262   Int_t currentpid = GetPdgCode();
 263   fNPids++;
 264   Float_t* aPidProbNew = new Float_t[fNPids];
 265   Int_t* aPidsNew = new Int_t[fNPids];
 266
 267   for (i = 0; i < fNPids-1;i++)//find a slot
 268    {
 269      if ( fPidProb[i] > prob)
 270       {
 271         AliDebug(9,Form("Copying entry %d",i));
 272         aPidProbNew[i] = fPidProb[i];
 273         aPidsNew[i] = fPids[i];
 274         totprob+=fPidProb[i];
 275       }
 276      else break;
 277    }
 278
 279   AliDebug(9,Form("SetPID","Setting new PID on entry %d",i));
 280   aPidProbNew[i] = prob;
 281   aPidsNew[i] = pdg;
 282   totprob+=prob;
 283
 284
 285   for (Int_t j = fNPids-1; j > i ;j--)//copy rest of old arays
 286    {
 287      AliDebug(9,Form("SetPID","Copying from old entry %d to new entry %d",j-1,j));
 288      aPidProbNew[j] = fPidProb[j-1];
 289      aPidsNew[j] = fPids[j-1];
 290      totprob+=fPidProb[j-1];
 291    }
 292
 293   delete [] fPidProb;
 294   delete [] fPids;
 295
 296   fPidProb = aPidProbNew;
 297   fPids = aPidsNew;
 298
 299   fPdgIdx = GetPidSlot(currentpid);
 300   if (fPdgIdx == -1) fPdgIdx = 0;
 301
 302   if (totprob > (1.0+0.000001))//space for numerical error
 303    {
 304      Warning("SetId","Total probability is greater than unity (%f)!!!",totprob);
 305      Print();
 306    }
 307 }
 308 //______________________________________________________________________________
 309
 310 Float_t AliAODParticle::GetPIDprobability(Int_t pdg) const
 311 {
 312 //Returns probability that this particle is the type of pdg
 313   Int_t idx = GetPidSlot(pdg);
 314   if (idx < 0) return 0.0;//such pid was not specified for this particle
 315   return fPidProb[idx];
 316 }
 317 //______________________________________________________________________________
 318
 319 const Char_t* AliAODParticle::GetName() const
 320 {
 321   //returns name of this particle
 322    static char def[4] = "XXX";
 323    const TParticlePDG *ap = TDatabasePDG::Instance()->GetParticle(GetPdgCode());
 324    if (ap) return ap->GetName();
 325    else    return def;
 326 }
 327 //______________________________________________________________________________
 328
 329 Int_t AliAODParticle::GetPidSlot(Int_t pdg) const
 330 {
 331  //returns position of the given PID in fPids (and fPidProb) array.
 332  if (fPids == 0x0) return -1;
 333  for (Int_t i = 0; i< fNPids; i++)
 334   {
 335    if (fPids[i] == pdg) return i;
 336   }
 337  return -1;
 338 }
 339 //______________________________________________________________________________
 340
 341 Int_t AliAODParticle::GetNthPid(Int_t idx) const
 342 {
 343   //returns PID sitting on slot idx in fPids
 344   if ( (idx < 0) || (idx >= fNPids) )
 345    {
 346      Error("GetNthPid","Out Of Bounds");
 347      return 0;
 348    }
 349   return fPids[idx];
 350 }
 351 //______________________________________________________________________________
 352
 353 Float_t AliAODParticle::GetNthPidProb(Int_t idx) const
 354 {
 355   //returns PID sitting on slot idx in fPidProb
 356   if ( (idx < 0) || (idx >= fNPids) )
 357    {
 358      Error("GetNthPid","Out Of Bounds");
 359      return 0;
 360    }
 361   return fPidProb[idx];
 362 }
 363 //______________________________________________________________________________
 364
 365 void AliAODParticle::Print(const Option_t * /*opt*/) const
 366 {
 367 //prints information about particle
 368   printf("____________________________________________________\n");
 369   printf("Idx: %d  PID: %d  Name: ",fIdxInEvent,GetPdgCode());
 370   TParticlePDG *pdgp = TDatabasePDG::Instance()->GetParticle(GetPdgCode());
 371   if (pdgp)
 372    {
 373      printf("%s Mass: %f\n",pdgp->GetName(),pdgp->Mass());
 374    }
 375   else
 376    {
 377      printf("Not known\n");
 378    }
 379
 380   printf("Px: %+f Py: %+f Pz: %+f E: %+f Calculated Mass: %f\nVx: %+f Vy: %+f Vz: %+f T: %+f\n",
 381           Px(),Py(),Pz(),E(),GetCalcMass(),Vx(),Vy(),Vz(),T());
 382
 383   for (Int_t i = 0; i < fNPids; i++)
 384    {
 385      printf("# %d  PID: %d  Probability %f name ",i,fPids[i],fPidProb[i]);
 386      const TParticlePDG *ap = TDatabasePDG::Instance()->GetParticle(fPids[i]);
 387      if (ap)
 388       {
 389         printf("%s Mass %f\n",ap->GetName(),ap->Mass());
 390       }
 391      else
 392       {
 393         printf("Not known\n");
 394       }
 395    }
 396
 397   if (fITSTrackPoints) fITSTrackPoints->Print();
 398   if (fTPCTrackPoints) fTPCTrackPoints->Print();
 399
 400 }
 401
 402 //______________________________________________________________________________
 403
 404 //void AliAODParticle::Streamer(TBuffer &b)
 405 //{
 406 //     // Stream all objects in the array to or from the I/O buffer.
 407 //   UInt_t R__s, R__c;
 408 //   Int_t i;
 409 //   if (b.IsReading())
 410 //    {
 411 //      delete [] fPids;
 412 //      delete [] fPidProb;
 413 //
 414 //      Version_t v = b.ReadVersion(&R__s, &R__c);
 415 //      if (v == 1)
 416 //       {
 417 //         AliAODParticle::Class()->ReadBuffer(b, this);
 418 //      }
 419 //      else
 420 //       {
 421 //        TObject::Streamer(b);
 422 //       b >> fPdgIdx;
 423 //       b >> fIdxInEvent;
 424 //
 425 //       b >> fNPids;
 426 //       Int_t* fPids = new Int_t[fNPids];
 427 //        Float_t* fPidProb = new Float_t[fNPids];
 428 //        for (i = 0;i<fNPids;i++)
 429 //         {
 430 //           b >> fPids[i];
 431 //         }
 432 //        for (i = 0;i<fNPids;i++)
 433 //        {
 434 //          b >> fPidProb[i];
 435 //         }
 436 //        b >> fCalcMass;
 437 //
 438 //        b >> fPx;
 439 //       b >> fPy;
 440 //       b >> fPz;
 441 //        b >> fE;
 442 //
 443 //       b >> fVx;
 444 //        b >> fVy;
 445 //        b >> fVz;
 446 //       b >> fVt;
 447 //       Info("Streamer","Read data");
 448 //        Print();
 449 //       }
 450 //
 451 //      b.CheckByteCount(R__s, R__c,AliAODParticle::IsA());
 452 //    }
 453 //  else
 454 //   {
 455 //     R__c = b.WriteVersion(AliAODParticle::IsA(), kTRUE);
 456 //     TObject::Streamer(b);
 457 //     Info("Streamer","Read data");
 458 //     Print();
 459 //
 460 //     b << fPdgIdx;
 461 //     b << fIdxInEvent;
 462 //     b << fNPids;
 463 //     for (i = 0;i<fNPids;i++)
 464 //      {
 465 //        b << fPids[i];
 466 //      }
 467 //      {
 468 //      {
 469 //     for (i = 0;i<fNPids;i++)
 470 //     {
 471 //        b << fPidProb[i];
 472 //      }
 473 //     b << fCalcMass;
 474 //
 475 //     b << fPx;
 476 //     b << fPy;
 477 //     b << fPz;
 478 //     b << fE;
 479 //
 480 //     b << fVx;
 481 //     b << fVy;
 482 //     b << fVz;
 483 //     b << fVt;
 484 //
 485 //    b.SetByteCount(R__c, kTRUE);
 486 //   }
 487 //}