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