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