STEER/ESD/AliTrackPointArray.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 //////////////////////////////////////////////////////////////////////////////
  17 //                          Class AliTrackPointArray                        //
  18 //   This class contains the ESD track space-points which are used during   //
  19 //   the alignment procedures. Each space-point consist of 3 coordinates    //
  20 //   (and their errors) and the index of the sub-detector which contains    //
  21 //   the space-point.                                                       //
  22 //   cvetan.cheshkov@cern.ch 3/11/2005                                      //
  23 //////////////////////////////////////////////////////////////////////////////
  24
  25 #include <TMath.h>
  26 #include <TMatrixD.h>
  27 #include <TMatrixDSym.h>
  28 #include <TGeoMatrix.h>
  29 #include <TMatrixDSymEigen.h>
  30 #include <TArrayI.h>
  31 #include "AliTrackPointArray.h"
  32 #include "AliPoolN.h"
  33
  34 ClassImp(AliTrackPointArray)
  35
  36 AliPoolN* AliTrackPointArray::fgPool = 0;
  37
  38 //______________________________________________________________________________
  39 AliTrackPointArray::AliTrackPointArray() :
  40   TObject(),
  41   fSorted(kFALSE),
  42   fNPoints(0),
  43   fX(0),
  44   fY(0),
  45   fZ(0),
  46   fCharge(0),
  47   fDriftTime(0),
  48   fChargeRatio(0),
  49   fClusterType(0),
  50   fIsExtra(0),
  51   fSize(0),
  52   fCov(0),
  53   fVolumeID(0)
  54 {
  55 }
  56
  57 //______________________________________________________________________________
  58 AliTrackPointArray::AliTrackPointArray(Int_t npoints):
  59   TObject(),
  60   fSorted(kFALSE),
  61   fNPoints(npoints),
  62   fX(           fgPool ? fgPool->BookF(npoints,&fX)           : new Float_t[npoints]),
  63   fY(           fgPool ? fgPool->BookF(npoints,&fY)           : new Float_t[npoints]),
  64   fZ(           fgPool ? fgPool->BookF(npoints,&fZ)           : new Float_t[npoints]),
  65   fCharge(      fgPool ? fgPool->BookF(npoints,&fCharge)      : new Float_t[npoints]),
  66   fDriftTime(   fgPool ? fgPool->BookF(npoints,&fDriftTime)   : new Float_t[npoints]),
  67   fChargeRatio( fgPool ? fgPool->BookF(npoints,&fChargeRatio) : new Float_t[npoints]),
  68   fClusterType( fgPool ? fgPool->BookI(npoints,&fClusterType) : new Int_t[npoints]),
  69   fIsExtra(     fgPool ? fgPool->BookB(npoints,&fIsExtra)     : new Bool_t[npoints]),
  70   fSize(6*npoints),
  71   fCov(         fgPool ? fgPool->BookF(fSize,&fCov)           : new Float_t[fSize]),
  72   fVolumeID(    fgPool ? fgPool->BookUS(npoints,&fVolumeID)   : new UShort_t[npoints])
  73 {
  74   // Constructor
  75   for (Int_t ip=npoints;ip--;){
  76     fX[ip]=0;
  77     fY[ip]=0;
  78     fZ[ip]=0;
  79     fCharge[ip]=0;
  80     fDriftTime[ip]=0;
  81     fChargeRatio[ip]=0;
  82     fClusterType[ip]=0;
  83     fIsExtra[ip]=kFALSE;
  84     fVolumeID[ip]=0;
  85     for (Int_t icov=6; icov--;) fCov[6*ip+icov]=0;
  86   }
  87 }
  88
  89 //______________________________________________________________________________
  90 AliTrackPointArray::AliTrackPointArray(const AliTrackPointArray &array):
  91   TObject(array),
  92   fSorted(array.fSorted),
  93   fNPoints(array.fNPoints),
  94   fX(           fgPool ? fgPool->BookF(fNPoints,&fX)           : new Float_t[fNPoints]),
  95   fY(           fgPool ? fgPool->BookF(fNPoints,&fY)           : new Float_t[fNPoints]),
  96   fZ(           fgPool ? fgPool->BookF(fNPoints,&fZ)           : new Float_t[fNPoints]),
  97   fCharge(      fgPool ? fgPool->BookF(fNPoints,&fCharge)      : new Float_t[fNPoints]),
  98   fDriftTime(   fgPool ? fgPool->BookF(fNPoints,&fDriftTime)   : new Float_t[fNPoints]),
  99   fChargeRatio( fgPool ? fgPool->BookF(fNPoints,&fChargeRatio) : new Float_t[fNPoints]),
 100   fClusterType( fgPool ? fgPool->BookI(fNPoints,&fClusterType) : new Int_t[fNPoints]),
 101   fIsExtra(     fgPool ? fgPool->BookB(fNPoints,&fIsExtra)     : new Bool_t[fNPoints]),
 102   fSize(6*fNPoints),
 103   fCov(         fgPool ? fgPool->BookF(fSize,&fCov)            : new Float_t[fSize]),
 104   fVolumeID(    fgPool ? fgPool->BookUS(fNPoints,&fVolumeID)   : new UShort_t[fNPoints])
 105 {
 106   // Copy constructor
 107   //
 108   memcpy(fX,array.fX,fNPoints*sizeof(Float_t));
 109   memcpy(fY,array.fY,fNPoints*sizeof(Float_t));
 110   memcpy(fZ,array.fZ,fNPoints*sizeof(Float_t));
 111   if (array.fCharge)      memcpy(fCharge,array.fCharge,fNPoints*sizeof(Float_t));
 112   else                    memset(fCharge, 0, fNPoints*sizeof(Float_t));
 113   //
 114   if (array.fDriftTime)   memcpy(fDriftTime,array.fDriftTime,fNPoints*sizeof(Float_t));
 115   else                    memset(fDriftTime, 0, fNPoints*sizeof(Float_t));
 116   //
 117   if (array.fChargeRatio) memcpy(fChargeRatio,array.fChargeRatio,fNPoints*sizeof(Float_t));
 118   else                    memset(fChargeRatio, 0, fNPoints*sizeof(Float_t));
 119   //
 120   if (array.fClusterType) memcpy(fClusterType,array.fClusterType,fNPoints*sizeof(Int_t));
 121   else                    memset(fClusterType, 0, fNPoints*sizeof(Int_t));
 122   //
 123   if (array.fIsExtra)     memcpy(fIsExtra,array.fIsExtra,fNPoints*sizeof(Bool_t));
 124   else                    memset(fIsExtra, 0, fNPoints*sizeof(Bool_t));
 125   //
 126   memcpy(fVolumeID,array.fVolumeID,fNPoints*sizeof(UShort_t));
 127   memcpy(fCov,array.fCov,fSize*sizeof(Float_t));
 128 }
 129
 130 //_____________________________________________________________________________
 131 AliTrackPointArray &AliTrackPointArray::operator =(const AliTrackPointArray& array)
 132 {
 133   // assignment operator
 134   //
 135   if(this==&array) return *this;
 136   ((TObject *)this)->operator=(array);
 137   fSorted = array.fSorted;
 138   fNPoints = array.fNPoints;
 139   fSize = array.fSize;
 140   //
 141   if (fgPool) Clear();
 142   if (!fgPool) {
 143     fX = new Float_t[fNPoints];
 144     fY = new Float_t[fNPoints];
 145     fZ = new Float_t[fNPoints];
 146     fCharge = new Float_t[fNPoints];
 147     fDriftTime = new Float_t[fNPoints];
 148     fChargeRatio = new Float_t[fNPoints];
 149     fClusterType = new Int_t[fNPoints];
 150     fIsExtra = new Bool_t[fNPoints];
 151     fVolumeID = new UShort_t[fNPoints];
 152     fCov = new Float_t[fSize];
 153   }
 154   else { // store in pools
 155     fX =           fgPool->BookF(fNPoints,&fX);
 156     fY =           fgPool->BookF(fNPoints,&fY);
 157     fZ =           fgPool->BookF(fNPoints,&fZ);
 158     fCharge =      fgPool->BookF(fNPoints,&fCharge);
 159     fDriftTime =   fgPool->BookF(fNPoints,&fDriftTime);
 160     fChargeRatio = fgPool->BookF(fNPoints,&fChargeRatio);
 161     fClusterType = fgPool->BookI(fNPoints,&fClusterType);
 162     fIsExtra =     fgPool->BookB(fNPoints,&fIsExtra);
 163     fCov =         fgPool->BookF(fSize,&fCov);
 164     fVolumeID =    fgPool->BookUS(fNPoints,&fVolumeID);
 165     //
 166   }
 167   //
 168   memcpy(fX,array.fX,fNPoints*sizeof(Float_t));
 169   memcpy(fY,array.fY,fNPoints*sizeof(Float_t));
 170   memcpy(fZ,array.fZ,fNPoints*sizeof(Float_t));
 171   memcpy(fCharge,array.fCharge,fNPoints*sizeof(Float_t));
 172   memcpy(fDriftTime,array.fDriftTime,fNPoints*sizeof(Float_t));
 173   memcpy(fChargeRatio,array.fChargeRatio,fNPoints*sizeof(Float_t));
 174   memcpy(fClusterType,array.fClusterType,fNPoints*sizeof(Int_t));
 175   memcpy(fIsExtra,array.fIsExtra,fNPoints*sizeof(Bool_t));
 176   memcpy(fVolumeID,array.fVolumeID,fNPoints*sizeof(UShort_t));
 177   memcpy(fCov,array.fCov,fSize*sizeof(Float_t));
 178
 179   return *this;
 180 }
 181
 182 //______________________________________________________________________________
 183 AliTrackPointArray::~AliTrackPointArray()
 184 {
 185   // Destructor
 186   //
 187   Clear();
 188 }
 189
 190 //______________________________________________________________________________
 191 void AliTrackPointArray::Clear(Option_t *)
 192 {
 193   // clean dynamical part
 194   TObject::Clear();
 195   if (!fX) return;
 196   if (!fgPool) {
 197     delete [] fX;
 198     delete [] fY;
 199     delete [] fZ;
 200     delete [] fCharge;
 201     delete [] fDriftTime;
 202     delete [] fChargeRatio;
 203     delete [] fClusterType;
 204     delete [] fIsExtra;
 205     delete [] fVolumeID;
 206     delete [] fCov;
 207   }
 208   else {
 209     if (!fgPool->IsReset()) { // UniqueID==0 means that the array was freed globally
 210       fgPool->FreeSlot(&fX);
 211       fgPool->FreeSlot(&fY);
 212       fgPool->FreeSlot(&fZ);
 213       fgPool->FreeSlot(&fCharge);
 214       fgPool->FreeSlot(&fDriftTime);
 215       fgPool->FreeSlot(&fChargeRatio);
 216       fgPool->FreeSlot(&fClusterType);
 217       fgPool->FreeSlot(&fIsExtra);
 218       fgPool->FreeSlot(&fCov);
 219       fgPool->FreeSlot(&fVolumeID);
 220     }
 221   }
 222 }
 223
 224 //______________________________________________________________________________
 225 Bool_t AliTrackPointArray::AddPoint(Int_t i, const AliTrackPoint *p)
 226 {
 227   // Add a point to the array at position i
 228   //
 229   if (i >= fNPoints) return kFALSE;
 230   fX[i] = p->GetX();
 231   fY[i] = p->GetY();
 232   fZ[i] = p->GetZ();
 233   fCharge[i] = p->GetCharge();
 234   fDriftTime[i] = p->GetDriftTime();
 235   fChargeRatio[i]=p->GetChargeRatio();
 236   fClusterType[i]=p->GetClusterType();
 237   fIsExtra[i] = p->IsExtra();
 238   fVolumeID[i] = p->GetVolumeID();
 239   memcpy(&fCov[6*i],p->GetCov(),6*sizeof(Float_t));
 240   return kTRUE;
 241 }
 242
 243
 244 //______________________________________________________________________________
 245 Bool_t AliTrackPointArray::GetPoint(AliTrackPoint &p, Int_t i) const
 246 {
 247   // Get the point at position i
 248   //
 249   if (i >= fNPoints) return kFALSE;
 250   p.SetXYZ(fX[i],fY[i],fZ[i],&fCov[6*i]);
 251   p.SetVolumeID(fVolumeID[i]);
 252   p.SetCharge(fCharge ? fCharge[i] : 0);
 253   p.SetDriftTime(fDriftTime ? fDriftTime[i] : 0);
 254   p.SetChargeRatio(fChargeRatio ? fChargeRatio[i] : 0);
 255   p.SetClusterType(fClusterType ? fClusterType[i] : 0);
 256   p.SetExtra(fIsExtra ? fIsExtra[i] : kFALSE);
 257   return kTRUE;
 258 }
 259
 260 //______________________________________________________________________________
 261 Bool_t AliTrackPointArray::HasVolumeID(UShort_t volid) const
 262 {
 263   // This method checks if the array
 264   // has at least one hit in the detector
 265   // volume defined by volid
 266   Bool_t check = kFALSE;
 267   for (Int_t ipoint = 0; ipoint < fNPoints; ipoint++)
 268     if (fVolumeID[ipoint] == volid) check = kTRUE;
 269
 270   return check;
 271 }
 272
 273 //______________________________________________________________________________
 274 void AliTrackPointArray::Sort(Bool_t down)
 275 {
 276   // Sort the array by the values of Y-coordinate of the track points.
 277   // The order is given by "down".
 278   // Optimized more for maintenance rather than for speed.
 279   if (fSorted) return;
 280   //
 281   static TArrayI indexAr;
 282   if (indexAr.GetSize()<fNPoints) indexAr.Set(fNPoints+10);
 283   //
 284   Int_t *index = indexAr.GetArray();
 285   TMath::Sort(fNPoints,fY,index,down);
 286
 287   AliTrackPoint p,p1;
 288   for (Int_t i = 0; i < fNPoints; i++) {
 289     GetPoint(p,index[i]); // swap points
 290     GetPoint(p1,i);       // and remember where point i was moved
 291     AddPoint(i,&p);
 292     AddPoint(index[i],&p1);
 293     for (int j=i;j<=fNPoints;j++) {if (index[j]==i) index[j]=index[i]; break;}
 294   }
 295   fSorted=kTRUE;
 296 }
 297
 298 ClassImp(AliTrackPoint)
 299
 300 //______________________________________________________________________________
 301 AliTrackPoint::AliTrackPoint() :
 302   TObject(),
 303   fX(0),
 304   fY(0),
 305   fZ(0),
 306   fCharge(0),
 307   fDriftTime(0),
 308   fChargeRatio(0),
 309   fClusterType(0),
 310   fIsExtra(kFALSE),
 311   fVolumeID(0)
 312 {
 313   // Default constructor
 314   //
 315   memset(fCov,0,6*sizeof(Float_t));
 316 }
 317
 318
 319 //______________________________________________________________________________
 320 AliTrackPoint::AliTrackPoint(Float_t x, Float_t y, Float_t z, const Float_t *cov, UShort_t volid, Float_t charge, Float_t drifttime,Float_t chargeratio, Int_t clutyp) :
 321   TObject(),
 322   fX(0),
 323   fY(0),
 324   fZ(0),
 325   fCharge(0),
 326   fDriftTime(0),
 327   fChargeRatio(0),
 328   fClusterType(0),
 329   fIsExtra(kFALSE),
 330   fVolumeID(0)
 331 {
 332   // Constructor
 333   //
 334   SetXYZ(x,y,z,cov);
 335   SetCharge(charge);
 336   SetDriftTime(drifttime);
 337   SetChargeRatio(chargeratio);
 338   SetClusterType(clutyp);
 339   SetVolumeID(volid);
 340 }
 341
 342 //______________________________________________________________________________
 343 AliTrackPoint::AliTrackPoint(const Float_t *xyz, const Float_t *cov, UShort_t volid, Float_t charge, Float_t drifttime,Float_t chargeratio, Int_t clutyp) :
 344   TObject(),
 345   fX(0),
 346   fY(0),
 347   fZ(0),
 348   fCharge(0),
 349   fDriftTime(0),
 350   fChargeRatio(0),
 351   fClusterType(0),
 352   fIsExtra(kFALSE),
 353   fVolumeID(0)
 354 {
 355   // Constructor
 356   //
 357   SetXYZ(xyz[0],xyz[1],xyz[2],cov);
 358   SetCharge(charge);
 359   SetDriftTime(drifttime);
 360   SetChargeRatio(chargeratio);
 361   SetVolumeID(volid);
 362   SetClusterType(clutyp);
 363 }
 364
 365 //______________________________________________________________________________
 366 AliTrackPoint::AliTrackPoint(const AliTrackPoint &p):
 367   TObject(p),
 368   fX(0),
 369   fY(0),
 370   fZ(0),
 371   fCharge(0),
 372   fDriftTime(0),
 373   fChargeRatio(0),
 374   fClusterType(0),
 375   fIsExtra(kFALSE),
 376   fVolumeID(0)
 377 {
 378   // Copy constructor
 379   //
 380   SetXYZ(p.fX,p.fY,p.fZ,&(p.fCov[0]));
 381   SetCharge(p.fCharge);
 382   SetDriftTime(p.fDriftTime);
 383   SetChargeRatio(p.fChargeRatio);
 384   SetClusterType(p.fClusterType);
 385   SetExtra(p.fIsExtra);
 386   SetVolumeID(p.fVolumeID);
 387 }
 388
 389 //_____________________________________________________________________________
 390 AliTrackPoint &AliTrackPoint::operator =(const AliTrackPoint& p)
 391 {
 392   // assignment operator
 393   //
 394   if(this==&p) return *this;
 395   ((TObject *)this)->operator=(p);
 396
 397   SetXYZ(p.fX,p.fY,p.fZ,&(p.fCov[0]));
 398   SetCharge(p.fCharge);
 399   SetDriftTime(p.fDriftTime);
 400   SetChargeRatio(p.fChargeRatio);
 401   SetClusterType(p.fClusterType);
 402   SetExtra(p.fIsExtra);
 403   SetVolumeID(p.fVolumeID);
 404
 405   return *this;
 406 }
 407
 408 //______________________________________________________________________________
 409 void AliTrackPoint::SetXYZ(Float_t x, Float_t y, Float_t z, const Float_t *cov)
 410 {
 411   // Set XYZ coordinates and their cov matrix
 412   //
 413   fX = x;
 414   fY = y;
 415   fZ = z;
 416   if (cov)
 417     memcpy(fCov,cov,6*sizeof(Float_t));
 418 }
 419
 420 //______________________________________________________________________________
 421 void AliTrackPoint::SetXYZ(const Float_t *xyz, const Float_t *cov)
 422 {
 423   // Set XYZ coordinates and their cov matrix
 424   //
 425   SetXYZ(xyz[0],xyz[1],xyz[2],cov);
 426 }
 427
 428 //______________________________________________________________________________
 429 void AliTrackPoint::SetCov(const Float_t *cov)
 430 {
 431   // Set XYZ cov matrix
 432   //
 433   if (cov)
 434     memcpy(fCov,cov,6*sizeof(Float_t));
 435 }
 436
 437 //______________________________________________________________________________
 438 void AliTrackPoint::GetXYZ(Float_t *xyz, Float_t *cov) const
 439 {
 440   xyz[0] = fX;
 441   xyz[1] = fY;
 442   xyz[2] = fZ;
 443   if (cov)
 444     memcpy(cov,fCov,6*sizeof(Float_t));
 445 }
 446
 447 //______________________________________________________________________________
 448 Float_t AliTrackPoint::GetResidual(const AliTrackPoint &p, Bool_t weighted) const
 449 {
 450   // This method calculates the track to space-point residuals. The track
 451   // interpolation is also stored as AliTrackPoint. Using the option
 452   // 'weighted' one can calculate the residual either with or without
 453   // taking into account the covariance matrix of the space-point and
 454   // track interpolation. The second case the residual becomes a pull.
 455
 456   Float_t res = 0;
 457
 458   if (!weighted) {
 459     Float_t xyz[3],xyzp[3];
 460     GetXYZ(xyz);
 461     p.GetXYZ(xyzp);
 462     res = (xyz[0]-xyzp[0])*(xyz[0]-xyzp[0])+
 463           (xyz[1]-xyzp[1])*(xyz[1]-xyzp[1])+
 464           (xyz[2]-xyzp[2])*(xyz[2]-xyzp[2]);
 465   }
 466   else {
 467     Float_t xyz[3],xyzp[3];
 468     Float_t cov[6],covp[6];
 469     GetXYZ(xyz,cov);
 470     TMatrixDSym mcov(3);
 471     mcov(0,0) = cov[0]; mcov(0,1) = cov[1]; mcov(0,2) = cov[2];
 472     mcov(1,0) = cov[1]; mcov(1,1) = cov[3]; mcov(1,2) = cov[4];
 473     mcov(2,0) = cov[2]; mcov(2,1) = cov[4]; mcov(2,2) = cov[5];
 474     p.GetXYZ(xyzp,covp);
 475     TMatrixDSym mcovp(3);
 476     mcovp(0,0) = covp[0]; mcovp(0,1) = covp[1]; mcovp(0,2) = covp[2];
 477     mcovp(1,0) = covp[1]; mcovp(1,1) = covp[3]; mcovp(1,2) = covp[4];
 478     mcovp(2,0) = covp[2]; mcovp(2,1) = covp[4]; mcovp(2,2) = covp[5];
 479     TMatrixDSym msum = mcov + mcovp;
 480     msum.Invert();
 481     //    mcov.Print(); mcovp.Print(); msum.Print();
 482     if (msum.IsValid()) {
 483       for (Int_t i = 0; i < 3; i++)
 484         for (Int_t j = 0; j < 3; j++)
 485           res += (xyz[i]-xyzp[i])*(xyz[j]-xyzp[j])*msum(i,j);
 486     }
 487   }
 488
 489   return res;
 490 }
 491
 492 //_____________________________________________________________________________
 493 Bool_t AliTrackPoint::GetPCA(const AliTrackPoint &p, AliTrackPoint &out) const
 494 {
 495   //
 496   // Get the intersection point between this point and
 497   // the point "p" belongs to.
 498   // The result is stored as a point 'out'
 499   // return kFALSE in case of failure.
 500   out.SetXYZ(0,0,0);
 501
 502   TMatrixD t(3,1);
 503   t(0,0)=GetX();
 504   t(1,0)=GetY();
 505   t(2,0)=GetZ();
 506
 507   TMatrixDSym tC(3);
 508   {
 509   const Float_t *cv=GetCov();
 510   tC(0,0)=cv[0]; tC(0,1)=cv[1]; tC(0,2)=cv[2];
 511   tC(1,0)=cv[1]; tC(1,1)=cv[3]; tC(1,2)=cv[4];
 512   tC(2,0)=cv[2]; tC(2,1)=cv[4]; tC(2,2)=cv[5];
 513   }
 514
 515   TMatrixD m(3,1);
 516   m(0,0)=p.GetX();
 517   m(1,0)=p.GetY();
 518   m(2,0)=p.GetZ();
 519
 520   TMatrixDSym mC(3);
 521   {
 522   const Float_t *cv=p.GetCov();
 523   mC(0,0)=cv[0]; mC(0,1)=cv[1]; mC(0,2)=cv[2];
 524   mC(1,0)=cv[1]; mC(1,1)=cv[3]; mC(1,2)=cv[4];
 525   mC(2,0)=cv[2]; mC(2,1)=cv[4]; mC(2,2)=cv[5];
 526   }
 527
 528   TMatrixDSym tmW(tC);
 529   tmW+=mC;
 530   tmW.Invert();
 531   if (!tmW.IsValid()) return kFALSE;
 532
 533   TMatrixD mW(tC,TMatrixD::kMult,tmW);
 534   TMatrixD tW(mC,TMatrixD::kMult,tmW);
 535
 536   TMatrixD mi(mW,TMatrixD::kMult,m);
 537   TMatrixD ti(tW,TMatrixD::kMult,t);
 538   ti+=mi;
 539
 540   TMatrixD iC(tC,TMatrixD::kMult,tmW);
 541   iC*=mC;
 542
 543   out.SetXYZ(ti(0,0),ti(1,0),ti(2,0));
 544   UShort_t id=p.GetVolumeID();
 545   out.SetVolumeID(id);
 546
 547   return kTRUE;
 548 }
 549
 550 //______________________________________________________________________________
 551 Float_t AliTrackPoint::GetAngle() const
 552 {
 553   // The method uses the covariance matrix of
 554   // the space-point in order to extract the
 555   // orientation of the detector plane.
 556   // The rotation in XY plane only is calculated.
 557
 558   Float_t phi= TMath::ATan2(TMath::Sqrt(fCov[0]),TMath::Sqrt(fCov[3]));
 559   if (fCov[1] > 0) {
 560     phi = TMath::Pi() - phi;
 561     if ((fY-fX) < 0) phi += TMath::Pi();
 562   }
 563   else {
 564     if ((fX+fY) < 0) phi += TMath::Pi();
 565   }
 566
 567   return phi;
 568
 569 }
 570
 571 //______________________________________________________________________________
 572 Bool_t AliTrackPoint::GetRotMatrix(TGeoRotation& rot) const
 573 {
 574   // Returns the orientation of the
 575   // sensitive layer (using cluster
 576   // covariance matrix).
 577   // Assumes that cluster has errors only in the layer's plane.
 578   // Return value is kTRUE in case of success.
 579
 580   TMatrixDSym mcov(3);
 581   {
 582     const Float_t *cov=GetCov();
 583     mcov(0,0)=cov[0]; mcov(0,1)=cov[1]; mcov(0,2)=cov[2];
 584     mcov(1,0)=cov[1]; mcov(1,1)=cov[3]; mcov(1,2)=cov[4];
 585     mcov(2,0)=cov[2]; mcov(2,1)=cov[4]; mcov(2,2)=cov[5];
 586   }
 587
 588   TMatrixDSymEigen eigen(mcov);
 589   TMatrixD eigenMatrix = eigen.GetEigenVectors();
 590
 591   rot.SetMatrix(eigenMatrix.GetMatrixArray());
 592
 593   return kTRUE;
 594 }
 595
 596
 597 //_____________________________________________________________________________
 598 AliTrackPoint& AliTrackPoint::Rotate(Float_t alpha) const
 599 {
 600   // Transform the space-point coordinates
 601   // and covariance matrix from global to
 602   // local (detector plane) coordinate system
 603   // XY plane rotation only
 604
 605   static AliTrackPoint p;
 606   p = *this;
 607
 608   Float_t xyz[3],cov[6];
 609   GetXYZ(xyz,cov);
 610
 611   Float_t sin = TMath::Sin(alpha), cos = TMath::Cos(alpha);
 612
 613   Float_t newxyz[3],newcov[6];
 614   newxyz[0] = cos*xyz[0] + sin*xyz[1];
 615   newxyz[1] = cos*xyz[1] - sin*xyz[0];
 616   newxyz[2] = xyz[2];
 617
 618   newcov[0] = cov[0]*cos*cos+
 619             2*cov[1]*sin*cos+
 620               cov[3]*sin*sin;
 621   newcov[1] = cov[1]*(cos*cos-sin*sin)+
 622              (cov[3]-cov[0])*sin*cos;
 623   newcov[2] = cov[2]*cos+
 624               cov[4]*sin;
 625   newcov[3] = cov[0]*sin*sin-
 626             2*cov[1]*sin*cos+
 627               cov[3]*cos*cos;
 628   newcov[4] = cov[4]*cos-
 629               cov[2]*sin;
 630   newcov[5] = cov[5];
 631
 632   p.SetXYZ(newxyz,newcov);
 633   p.SetVolumeID(GetVolumeID());
 634
 635   return p;
 636 }
 637
 638 //_____________________________________________________________________________
 639 AliTrackPoint& AliTrackPoint::MasterToLocal() const
 640 {
 641   // Transform the space-point coordinates
 642   // and the covariance matrix from the
 643   // (master) to the local (tracking)
 644   // coordinate system
 645
 646   Float_t alpha = GetAngle();
 647   return Rotate(alpha);
 648 }
 649
 650 //_____________________________________________________________________________
 651 void AliTrackPoint::Print(Option_t *) const
 652 {
 653   // Print the space-point coordinates and
 654   // covariance matrix
 655
 656   printf("VolumeID=%d\n", GetVolumeID());
 657   printf("X = %12.6f    Tx = %12.6f%12.6f%12.6f\n", fX, fCov[0], fCov[1], fCov[2]);
 658   printf("Y = %12.6f    Ty = %12.6f%12.6f%12.6f\n", fY, fCov[1], fCov[3], fCov[4]);
 659   printf("Z = %12.6f    Tz = %12.6f%12.6f%12.6f\n", fZ, fCov[2], fCov[4], fCov[5]);
 660   printf("Charge = %f\n", fCharge);
 661   printf("Drift Time = %f\n", fDriftTime);
 662   printf("Charge Ratio  = %f\n", fChargeRatio);
 663   printf("Cluster Type  = %d\n", fClusterType);
 664   if(fIsExtra) printf("This is an extra point\n");
 665
 666 }
 667
 668
 669 //________________________________
 670 void AliTrackPoint::SetAlignCovMatrix(const TMatrixDSym& alignparmtrx){
 671   // Add the uncertainty on the cluster position due to alignment
 672   // (using the 6x6 AliAlignObj Cov. Matrix alignparmtrx) to the already
 673   // present Cov. Matrix
 674
 675   TMatrixDSym cov(3);
 676   TMatrixD coval(3,3);
 677   TMatrixD jacob(3,6);
 678   Float_t newcov[6];
 679
 680   cov(0,0)=fCov[0];
 681   cov(1,0)=cov(0,1)=fCov[1];
 682   cov(2,0)=cov(0,2)=fCov[2];
 683   cov(1,1)=fCov[3];
 684   cov(2,1)=cov(1,2)=fCov[4];
 685   cov(2,2)=fCov[5];
 686
 687   jacob(0,0) = 1;      jacob(1,0) = 0;       jacob(2,0) = 0;
 688   jacob(0,1) = 0;      jacob(1,1) = 1;       jacob(2,1) = 0;
 689   jacob(0,2) = 0;      jacob(1,2) = 0;       jacob(2,2) = 1;
 690   jacob(0,3) = 0;      jacob(1,3) =-fZ;      jacob(2,3) = fY;
 691   jacob(0,4) = fZ;     jacob(1,4) = 0;       jacob(2,4) =-fX;
 692   jacob(0,5) = -fY;    jacob(1,5) = fX;      jacob(2,5) = 0;
 693
 694   TMatrixD jacobT=jacob.T();jacob.T();
 695
 696   coval=jacob*alignparmtrx*jacobT+cov;
 697
 698
 699   newcov[0]=coval(0,0);
 700   newcov[1]=coval(1,0);
 701   newcov[2]=coval(2,0);
 702   newcov[3]=coval(1,1);
 703   newcov[4]=coval(2,1);
 704   newcov[5]=coval(2,2);
 705
 706   SetXYZ(fX,fY,fZ,newcov);
 707
 708 }
 709
 710