STEER/AliAODVertex.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2007, 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 //     AOD track base class
  20 //     Base class for Analysis Object Data
  21 //     Generic version
  22 //     Author: Markus Oldenburg, CERN
  23 //     Inheritance from AliVVertex: A. Dainese
  24 //-------------------------------------------------------------------------
  25
  26 #include "AliAODVertex.h"
  27 #include "AliAODTrack.h"
  28
  29 ClassImp(AliAODVertex)
  30
  31 //______________________________________________________________________________
  32 AliAODVertex::AliAODVertex() :
  33   AliVVertex(),
  34   fChi2perNDF(-999.),
  35   fID(-1),
  36   fType(kUndef),
  37   fNprong(0),
  38   fIprong(0),
  39   fCovMatrix(NULL),
  40   fParent(),
  41   fDaughters(),
  42   fProngs(NULL)
  43   {
  44   // default constructor
  45
  46   fPosition[0] = fPosition[1] = fPosition[2] = -999.;
  47 }
  48
  49 //______________________________________________________________________________
  50 AliAODVertex::AliAODVertex(const Double_t position[3],
  51                            const Double_t covMatrix[6],
  52                            Double_t  chi2perNDF,
  53                            TObject  *parent,
  54                            Short_t id,
  55                            Char_t vtype,
  56                            Int_t  nprong) :
  57   AliVVertex(),
  58   fChi2perNDF(chi2perNDF),
  59   fID(id),
  60   fType(vtype),
  61   fNprong(nprong),
  62   fIprong(0),
  63   fCovMatrix(NULL),
  64   fParent(parent),
  65   fDaughters(),
  66   fProngs(0)
  67 {
  68   // constructor
  69
  70   SetPosition(position);
  71   if (covMatrix) SetCovMatrix(covMatrix);
  72   MakeProngs();
  73 }
  74
  75 //______________________________________________________________________________
  76 AliAODVertex::AliAODVertex(const Float_t position[3],
  77                            const Float_t  covMatrix[6],
  78                            Double_t  chi2perNDF,
  79                            TObject  *parent,
  80                            Short_t id,
  81                            Char_t vtype,
  82                            Int_t nprong) :
  83
  84   AliVVertex(),
  85   fChi2perNDF(chi2perNDF),
  86   fID(id),
  87   fType(vtype),
  88   fNprong(nprong),
  89   fIprong(0),
  90   fCovMatrix(NULL),
  91   fParent(parent),
  92   fDaughters(),
  93   fProngs(0)
  94 {
  95   // constructor
  96
  97   SetPosition(position);
  98   if (covMatrix) SetCovMatrix(covMatrix);
  99   MakeProngs();
 100 }
 101
 102 //______________________________________________________________________________
 103 AliAODVertex::AliAODVertex(const Double_t position[3],
 104                            Double_t  chi2perNDF,
 105                            Char_t vtype,
 106                            Int_t nprong) :
 107   AliVVertex(),
 108   fChi2perNDF(chi2perNDF),
 109   fID(-1),
 110   fType(vtype),
 111   fNprong(nprong),
 112   fIprong(0),
 113   fCovMatrix(NULL),
 114   fParent(),
 115   fDaughters(),
 116   fProngs(0)
 117 {
 118   // constructor without covariance matrix
 119
 120   SetPosition(position);
 121   MakeProngs();
 122 }
 123
 124 //______________________________________________________________________________
 125 AliAODVertex::AliAODVertex(const Float_t position[3],
 126                            Double_t  chi2perNDF,
 127                            Char_t vtype, Int_t nprong) :
 128   AliVVertex(),
 129   fChi2perNDF(chi2perNDF),
 130   fID(-1),
 131   fType(vtype),
 132   fNprong(nprong),
 133   fIprong(0),
 134   fCovMatrix(NULL),
 135   fParent(),
 136   fDaughters(),
 137   fProngs(0)
 138 {
 139   // constructor without covariance matrix
 140
 141   SetPosition(position);
 142   MakeProngs();
 143 }
 144
 145 //______________________________________________________________________________
 146 AliAODVertex::~AliAODVertex()
 147 {
 148   // Destructor
 149
 150   delete fCovMatrix;
 151   if (fNprong > 0) delete[] fProngs;
 152 }
 153
 154 //______________________________________________________________________________
 155 AliAODVertex::AliAODVertex(const AliAODVertex& vtx) :
 156   AliVVertex(vtx),
 157   fChi2perNDF(vtx.fChi2perNDF),
 158   fID(vtx.fID),
 159   fType(vtx.fType),
 160   fNprong(vtx.fNprong),
 161   fIprong(vtx.fIprong),
 162   fCovMatrix(NULL),
 163   fParent(vtx.fParent),
 164   fDaughters(vtx.fDaughters),
 165   fProngs(0)
 166 {
 167   // Copy constructor.
 168
 169   for (int i = 0; i < 3; i++)
 170     fPosition[i] = vtx.fPosition[i];
 171
 172   if (vtx.fCovMatrix) fCovMatrix=new AliAODRedCov<3>(*vtx.fCovMatrix);
 173   MakeProngs();
 174   for (int i = 0; i < fNprong; i++) {
 175       fProngs[i] = vtx.fProngs[i];
 176   }
 177 }
 178
 179 //______________________________________________________________________________
 180 AliAODVertex& AliAODVertex::operator=(const AliAODVertex& vtx)
 181 {
 182   // Assignment operator
 183   if (this != &vtx) {
 184
 185     // name and type
 186     AliVVertex::operator=(vtx);
 187
 188     //momentum
 189     for (int i = 0; i < 3; i++)
 190       fPosition[i] = vtx.fPosition[i];
 191
 192     fChi2perNDF = vtx.fChi2perNDF;
 193     fID = vtx.fID;
 194     fType = vtx.fType;
 195
 196     //covariance matrix
 197     delete fCovMatrix;
 198     fCovMatrix = NULL;
 199     if (vtx.fCovMatrix) fCovMatrix=new AliAODRedCov<3>(*vtx.fCovMatrix);
 200
 201     //other stuff
 202     fParent = vtx.fParent;
 203     fDaughters = vtx.fDaughters;
 204     fNprong    = vtx.fNprong;
 205     fIprong    = vtx.fIprong;
 206
 207     MakeProngs();
 208     for (int i = 0; i < fNprong; i++) {
 209         fProngs[i] = vtx.fProngs[i];
 210     }
 211   }
 212
 213   return *this;
 214 }
 215
 216 //______________________________________________________________________________
 217 void AliAODVertex::AddDaughter(TObject *daughter)
 218 {
 219   // Add reference to daughter track
 220     if (!fProngs) {
 221         if (fDaughters.GetEntries()==0) {
 222             TRefArray* arr = &fDaughters;
 223             new(arr)TRefArray(TProcessID::GetProcessWithUID(daughter));
 224         }
 225         fDaughters.Add(daughter);
 226     } else {
 227         if (fIprong < fNprong) {
 228             fProngs[fIprong++] = daughter;
 229         } else {
 230             AliWarning("Number of daughters out of range !\n");
 231         }
 232     }
 233   return;
 234 }
 235
 236
 237 //______________________________________________________________________________
 238 template <class T> void AliAODVertex::GetSigmaXYZ(T sigma[3]) const
 239 {
 240   // Return errors on vertex position in thrust frame
 241
 242   if(fCovMatrix) {
 243     sigma[0]=fCovMatrix[3]; //GetCovXZ
 244     sigma[1]=fCovMatrix[4]; //GetCovYZ
 245     sigma[2]=fCovMatrix[5]; //GetCovZZ
 246   } else
 247     sigma[0]=sigma[1]=sigma[2]=-999.;
 248
 249   /*
 250   for (int i = 0, j = 6; i < 3; i++) {
 251     j -= i+1;
 252     sigma[2-i] = fCovMatrix ? TMath::Sqrt(fCovMatrix[j]) : -999.;
 253   }
 254   */
 255 }
 256
 257 //______________________________________________________________________________
 258 Int_t AliAODVertex::GetNContributors() const
 259 {
 260   // Returns the number of tracks used to fit this vertex.
 261   Int_t cont = 0;
 262
 263   if (!strcmp(GetTitle(), "vertexer: 3D")) {
 264     cont = fNContributors;
 265   } else {
 266     for (Int_t iDaug = 0; iDaug < GetNDaughters(); iDaug++) {
 267       if (((AliAODTrack*)fDaughters.At(iDaug))->GetUsedForVtxFit()) cont++;
 268     }
 269   }
 270   return cont;
 271 }
 272
 273 //______________________________________________________________________________
 274 Bool_t AliAODVertex::HasDaughter(TObject *daughter) const
 275 {
 276   // Checks if the given daughter (particle) is part of this vertex.
 277     if (!fProngs) {
 278         TRefArrayIter iter(&fDaughters);
 279         while (TObject *daugh = iter.Next()) {
 280             if (daugh == daughter) return kTRUE;
 281         }
 282         return kFALSE;
 283     } else {
 284         Bool_t has = kFALSE;
 285         for (int i = 0; i < fNprong; i++) {
 286             if (fProngs[i].GetObject() == daughter) has = kTRUE;
 287         }
 288         return has;
 289     }
 290 }
 291
 292 //______________________________________________________________________________
 293 Double_t AliAODVertex::RotatedCovMatrixXX(Double_t phi, Double_t theta) const
 294 {
 295   // XX term of covariance matrix after rotation by phi around z-axis
 296   // and, then, by theta around new y-axis
 297
 298   if (!fCovMatrix) {
 299     //AliFatal("Covariance matrix not set");
 300     return -999.;
 301   }
 302
 303   Double_t covMatrix[6];
 304
 305   GetCovMatrix(covMatrix);
 306
 307   Double_t cp = TMath::Cos(phi);
 308   Double_t sp = TMath::Sin(phi);
 309   Double_t ct = TMath::Cos(theta);
 310   Double_t st = TMath::Sin(theta);
 311   return
 312      covMatrix[0]*cp*cp*ct*ct  // GetCovXX
 313     +covMatrix[1]*2.*cp*sp*ct*ct  // GetCovXY
 314     +covMatrix[3]*2.*cp*ct*st  // GetCovXZ
 315     +covMatrix[2]*sp*sp*ct*ct  // GetCovYY
 316     +covMatrix[4]*2.*sp*ct*st  // GetCovYZ
 317     +covMatrix[5]*st*st;  // GetCovZZ
 318 }
 319
 320 //______________________________________________________________________________
 321 Double_t AliAODVertex::RotatedCovMatrixXY(Double_t phi, Double_t theta) const
 322 {
 323   // XY term of covariance matrix after rotation by phi around z-axis
 324   // and, then, by theta around new y-axis
 325
 326   if (!fCovMatrix) {
 327     //AliFatal("Covariance matrix not set");
 328     return -999.;
 329   }
 330
 331   Double_t covMatrix[6];
 332
 333   GetCovMatrix(covMatrix);
 334
 335   Double_t cp = TMath::Cos(phi);
 336   Double_t sp = TMath::Sin(phi);
 337   Double_t ct = TMath::Cos(theta);
 338   Double_t st = TMath::Sin(theta);
 339   return
 340     -covMatrix[0]*cp*sp*ct  // GetCovXX
 341     +covMatrix[1]*ct*(cp*cp-sp*sp)  // GetCovXY
 342     -covMatrix[3]*sp*st  // GetCovXZ
 343     +covMatrix[2]*cp*sp*ct  // GetCovYY
 344     +covMatrix[4]*cp*st;  // GetCovYZ
 345 }
 346
 347 //______________________________________________________________________________
 348 Double_t AliAODVertex::RotatedCovMatrixXZ(Double_t phi, Double_t theta) const
 349 {
 350   // XZ term of covariance matrix after rotation by phi around z-axis
 351   // and, then, by theta around new y-axis
 352
 353   if (!fCovMatrix) {
 354     //AliFatal("Covariance matrix not set");
 355     return -999.;
 356   }
 357
 358   Double_t covMatrix[6];
 359
 360   GetCovMatrix(covMatrix);
 361
 362   Double_t cp = TMath::Cos(phi);
 363   Double_t sp = TMath::Sin(phi);
 364   Double_t ct = TMath::Cos(theta);
 365   Double_t st = TMath::Sin(theta);
 366   return
 367     -covMatrix[0]*cp*cp*ct*st  // GetCovXX
 368     -covMatrix[1]*2.*cp*sp*ct*st  // GetCovXY
 369     +covMatrix[3]*cp*(ct*ct-st*st)  // GetCovXZ
 370     -covMatrix[2]*sp*sp*ct*st  // GetCovYY
 371     +covMatrix[4]*sp*(ct*ct-st*st)  // GetCovYZ
 372     +covMatrix[5]*ct*st;  // GetCovZZ
 373 }
 374
 375 //______________________________________________________________________________
 376 Double_t AliAODVertex::RotatedCovMatrixYY(Double_t phi) const
 377 {
 378   // YY term of covariance matrix after rotation by phi around z-axis
 379   // and, then, by theta around new y-axis
 380
 381   if (!fCovMatrix) {
 382     //AliFatal("Covariance matrix not set");
 383     return -999.;
 384   }
 385
 386   Double_t covMatrix[6];
 387
 388   GetCovMatrix(covMatrix);
 389
 390   Double_t cp = TMath::Cos(phi);
 391   Double_t sp = TMath::Sin(phi);
 392   return
 393      covMatrix[0]*sp*sp  // GetCovXX
 394     -covMatrix[1]*2.*cp*sp  // GetCovXY
 395     +covMatrix[2]*cp*cp;  // GetCovYY
 396 }
 397
 398 //______________________________________________________________________________
 399 Double_t AliAODVertex::RotatedCovMatrixYZ(Double_t phi, Double_t theta) const
 400 {
 401   // YZ term of covariance matrix after rotation by phi around z-axis
 402   // and, then, by theta around new y-axis
 403
 404   if (!fCovMatrix) {
 405     //AliFatal("Covariance matrix not set");
 406     return -999.;
 407   }
 408
 409   Double_t covMatrix[6];
 410
 411   GetCovMatrix(covMatrix);
 412
 413   Double_t cp = TMath::Cos(phi);
 414   Double_t sp = TMath::Sin(phi);
 415   Double_t ct = TMath::Cos(theta);
 416   Double_t st = TMath::Sin(theta);
 417   return
 418      covMatrix[0]*cp*sp*st  // GetCovXX
 419     +covMatrix[1]*st*(sp*sp-cp*cp)  // GetCovXY
 420     -covMatrix[3]*sp*ct  // GetCovXZ
 421     -covMatrix[2]*cp*sp*st  // GetCovYY
 422     +covMatrix[4]*cp*ct;  // GetCovYZ
 423 }
 424
 425 //______________________________________________________________________________
 426 Double_t AliAODVertex::RotatedCovMatrixZZ(Double_t phi, Double_t theta) const
 427 {
 428   // ZZ term of covariance matrix after rotation by phi around z-axis
 429   // and, then, by theta around new y-axis
 430
 431   if (!fCovMatrix) {
 432     //AliFatal("Covariance matrix not set");
 433     return -999.;
 434   }
 435
 436   Double_t covMatrix[6];
 437
 438   GetCovMatrix(covMatrix);
 439
 440   Double_t cp = TMath::Cos(phi);
 441   Double_t sp = TMath::Sin(phi);
 442   Double_t ct = TMath::Cos(theta);
 443   Double_t st = TMath::Sin(theta);
 444   return
 445      covMatrix[0]*cp*cp*st*st  // GetCovXX
 446     +covMatrix[1]*2.*cp*sp*st*st  // GetCovXY
 447     -covMatrix[3]*2.*cp*ct*st  // GetCovXZ
 448     +covMatrix[2]*sp*sp*st*st  // GetCovYY
 449     -covMatrix[4]*2.*sp*sp*ct*st  // GetCovYZ
 450     +covMatrix[5]*ct*ct;  // GetCovZZ
 451 }
 452
 453 //______________________________________________________________________________
 454 Double_t AliAODVertex::DistanceToVertex(AliAODVertex *vtx) const
 455 {
 456   // distance in 3D to another AliAODVertex
 457
 458   Double_t dx = GetX()-vtx->GetX();
 459   Double_t dy = GetY()-vtx->GetY();
 460   Double_t dz = GetZ()-vtx->GetZ();
 461
 462   return TMath::Sqrt(dx*dx+dy*dy+dz*dz);
 463 }
 464
 465 //______________________________________________________________________________
 466 Double_t AliAODVertex::DistanceXYToVertex(AliAODVertex *vtx) const
 467 {
 468   // distance in XY to another AliAODVertex
 469
 470   Double_t dx = GetX()-vtx->GetX();
 471   Double_t dy = GetY()-vtx->GetY();
 472
 473   return TMath::Sqrt(dx*dx+dy*dy);
 474 }
 475
 476 //______________________________________________________________________________
 477 Double_t AliAODVertex::ErrorDistanceToVertex(AliAODVertex *vtx) const
 478 {
 479   // error on the distance in 3D to another AliAODVertex
 480
 481   Double_t phi,theta;
 482   PhiAndThetaToVertex(vtx,phi,theta);
 483   // error2 due to this vertex
 484   Double_t error2 = RotatedCovMatrixXX(phi,theta);
 485   // error2 due to vtx vertex
 486   Double_t error2vtx = vtx->RotatedCovMatrixXX(phi,theta);
 487
 488   return TMath::Sqrt(error2+error2vtx);
 489 }
 490
 491 //______________________________________________________________________________
 492 Double_t AliAODVertex::ErrorDistanceXYToVertex(AliAODVertex *vtx) const
 493 {
 494   // error on the distance in XY to another AliAODVertex
 495
 496   Double_t phi,theta;
 497   PhiAndThetaToVertex(vtx,phi,theta);
 498   // error2 due to this vertex
 499   Double_t error2 = RotatedCovMatrixXX(phi);
 500   // error2 due to vtx vertex
 501   Double_t error2vtx = vtx->RotatedCovMatrixXX(phi);
 502
 503   return TMath::Sqrt(error2+error2vtx);
 504 }
 505
 506 //______________________________________________________________________________
 507 template <class T, class P>
 508 void AliAODVertex::PhiAndThetaToVertex(AliAODVertex *vtx, P &phi, T &theta) const
 509 {
 510   // rotation angles around z-axis (phi) and around new y-axis (theta)
 511   // with which vtx is seen (used by RotatedCovMatrix... methods)
 512
 513   phi = TMath::Pi()+TMath::ATan2(-vtx->GetY()+GetY(),-vtx->GetX()+GetX());
 514   Double_t vtxxphi = vtx->GetX()*TMath::Cos(phi)+vtx->GetY()*TMath::Sin(phi);
 515   Double_t xphi = GetX()*TMath::Cos(phi)+GetY()*TMath::Sin(phi);
 516   theta = TMath::ATan2(vtx->GetZ()-GetZ(),vtxxphi-xphi);
 517 }
 518
 519 //______________________________________________________________________________
 520 void AliAODVertex::PrintIndices() const
 521 {
 522   // Print indices of particles originating form this vertex
 523
 524   TRefArrayIter iter(&fDaughters);
 525   while (TObject *daugh = iter.Next()) {
 526     printf("Particle %p originates from this vertex.\n", static_cast<void*>(daugh));
 527   }
 528 }
 529
 530 //______________________________________________________________________________
 531 void AliAODVertex::Print(Option_t* /*option*/) const
 532 {
 533   // Print information of all data members
 534
 535   printf("Vertex position:\n");
 536   printf("     x = %f\n", fPosition[0]);
 537   printf("     y = %f\n", fPosition[1]);
 538   printf("     z = %f\n", fPosition[2]);
 539   printf(" parent particle: %p\n", static_cast<void*>(fParent.GetObject()));
 540   printf(" origin of %d particles\n", fDaughters.GetEntriesFast());
 541   printf(" vertex type %d\n", fType);
 542
 543   /*
 544   if (fCovMatrix) {
 545     printf("Covariance matrix:\n");
 546     printf(" %12.10f  %12.10f  %12.10f\n %12.10f  %12.10f  %12.10f\n %12.10f  %12.10f  %12.10f\n",
 547            fCovMatrix[0],
 548            fCovMatrix[1],
 549            fCovMatrix[3],
 550            fCovMatrix[1],
 551            fCovMatrix[2],
 552            fCovMatrix[4],
 553            fCovMatrix[3],
 554            fCovMatrix[4],
 555            fCovMatrix[5]);
 556            } */
 557   printf(" Chi^2/NDF = %f\n", fChi2perNDF);
 558 }
 559