STEER/AliTrackResidualsFast.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 //   Implementation of the derived class for track residuals
  18 //   based on linear chi2 minimization (in approximation of
  19 //   small alignment angles and translations)
  20 //
  21 //-----------------------------------------------------------------
  22
  23 #include <TMinuit.h>
  24 #include <TGeoMatrix.h>
  25
  26 #include "AliAlignObj.h"
  27 #include "AliTrackPointArray.h"
  28 #include "AliTrackResidualsFast.h"
  29
  30 ClassImp(AliTrackResidualsFast)
  31
  32 //______________________________________________________________________________
  33 AliTrackResidualsFast::AliTrackResidualsFast():AliTrackResiduals()
  34 {
  35   // Default constructor
  36   for (Int_t i = 0; i < 16; i++) fSum[i] = 0;
  37   fSumR = 0;
  38 }
  39
  40 //______________________________________________________________________________
  41 AliTrackResidualsFast::AliTrackResidualsFast(Int_t ntracks):
  42   AliTrackResiduals(ntracks)
  43 {
  44   // Constructor
  45   for (Int_t i = 0; i < 16; i++) fSum[i] = 0;
  46   fSumR = 0;
  47 }
  48
  49 //______________________________________________________________________________
  50 AliTrackResidualsFast::AliTrackResidualsFast(const AliTrackResidualsFast &res):
  51   AliTrackResiduals(res)
  52 {
  53   // Copy constructor
  54   for (Int_t i = 0; i < 16; i++) fSum[i] = res.fSum[i];
  55   fSumR = res.fSumR;
  56 }
  57
  58 //______________________________________________________________________________
  59 AliTrackResidualsFast &AliTrackResidualsFast::operator= (const AliTrackResidualsFast& res)
  60 {
  61   // Assignment operator
  62  ((AliTrackResiduals *)this)->operator=(res);
  63  for (Int_t i = 0; i < 16; i++) fSum[i] = res.fSum[i];
  64  fSumR = res.fSumR;
  65
  66  return *this;
  67 }
  68
  69 //______________________________________________________________________________
  70 Bool_t AliTrackResidualsFast::Minimize()
  71 {
  72   // Implementation of fast linear Chi2
  73   // based minimization of track residuals sum
  74   for (Int_t i = 0; i < 16; i++) fSum[i] = 0;
  75   fSumR = 0;
  76
  77   AliTrackPoint p1,p2;
  78
  79   for (Int_t itrack = 0; itrack < fLast; itrack++) {
  80     if (!fVolArray[itrack] || !fTrackArray[itrack]) continue;
  81     for (Int_t ipoint = 0; ipoint < fVolArray[itrack]->GetNPoints(); ipoint++) {
  82       fVolArray[itrack]->GetPoint(p1,ipoint);
  83       fTrackArray[itrack]->GetPoint(p2,ipoint);
  84       AddPoints(p1,p2);
  85     }
  86   }
  87
  88   return Update();
  89
  90   // debug info
  91 //   Float_t chi2 = 0;
  92 //   for (Int_t itrack = 0; itrack < fLast; itrack++) {
  93 //     if (!fVolArray[itrack] || !fTrackArray[itrack]) continue;
  94 //     for (Int_t ipoint = 0; ipoint < fVolArray[itrack]->GetNPoints(); ipoint++) {
  95 //       fVolArray[itrack]->GetPoint(p1,ipoint);
  96 //       fAlignObj->Transform(p1);
  97 //       fTrackArray[itrack]->GetPoint(p2,ipoint);
  98 //       Float_t residual = p2.GetResidual(p1,kFALSE);
  99 //       chi2 += residual;
 100 //     }
 101 //   }
 102 //   printf("Final chi2 = %f\n",chi2);
 103 }
 104
 105 //______________________________________________________________________________
 106 void AliTrackResidualsFast::AddPoints(AliTrackPoint &p, AliTrackPoint &pprime)
 107 {
 108   // Update the sums used for
 109   // the linear chi2 minimization
 110   Float_t xyz[3],xyzp[3];
 111   p.GetXYZ(xyz); pprime.GetXYZ(xyzp);
 112   Double_t weight = 1;
 113   weight *= weight;
 114   fSum[0] += weight;
 115   fSum[1] += xyz[0]*weight;
 116   fSum[2] += xyz[1]*weight;
 117   fSum[3] += xyz[2]*weight;
 118   fSum[4] += (xyz[0]*xyz[0]+xyz[1]*xyz[1])*weight;
 119   fSum[5] += (xyz[0]*xyz[0]+xyz[2]*xyz[2])*weight;
 120   fSum[6] += (xyz[1]*xyz[1]+xyz[2]*xyz[2])*weight;
 121   fSum[7] -= xyz[0]*xyz[1]*weight;
 122   fSum[8] -= xyz[0]*xyz[2]*weight;
 123   fSum[9] -= xyz[1]*xyz[2]*weight;
 124   fSum[10] += (xyzp[0]-xyz[0])*weight;
 125   fSum[11] += (xyzp[1]-xyz[1])*weight;
 126   fSum[12] += (xyzp[2]-xyz[2])*weight;
 127   fSum[13] += (xyzp[2]*xyz[1]-xyzp[1]*xyz[2])*weight;
 128   fSum[14] += (xyzp[0]*xyz[2]-xyzp[2]*xyz[0])*weight;
 129   fSum[15] += (xyzp[1]*xyz[0]-xyzp[0]*xyz[1])*weight;
 130
 131   fSumR += ((xyzp[0]-xyz[0])*(xyzp[0]-xyz[0])+
 132             (xyzp[1]-xyz[1])*(xyzp[1]-xyz[1])+
 133             (xyzp[2]-xyz[2])*(xyzp[2]-xyz[2]))*weight;
 134
 135   fNdf += 3;
 136 }
 137
 138 //______________________________________________________________________________
 139 Bool_t AliTrackResidualsFast::Update()
 140 {
 141   // Find the alignment parameters
 142   // by using the already accumulated
 143   // sums
 144   TMatrixDSym     smatrix(6);
 145   TMatrixD        sums(1,6);
 146
 147   smatrix(0,0) = fSum[0]; smatrix(1,1) = fSum[0]; smatrix(2,2)=fSum[0];
 148   smatrix(3,3) = fSum[6]; smatrix(4,4) = fSum[5]; smatrix(5,5)=fSum[4];
 149
 150   smatrix(0,1) = 0; smatrix(0,2) = 0;
 151   smatrix(0,3) = 0; smatrix(0,4) = fSum[3]; smatrix(0,5) = -fSum[2];
 152   smatrix(1,2) = 0;
 153   smatrix(1,3) = -fSum[3]; smatrix(1,4) = 0; smatrix(1,5) = fSum[1];
 154   smatrix(2,3) = fSum[2]; smatrix(2,4) = -fSum[1]; smatrix(2,5) = 0;
 155
 156   smatrix(3,4) = fSum[7]; smatrix(3,5) = fSum[8];
 157   smatrix(4,5) = fSum[9];
 158
 159   sums(0,0)    = fSum[10]; sums(0,1) = fSum[11]; sums(0,2) = fSum[12];
 160   sums(0,3)    = fSum[13]; sums(0,4) = fSum[14]; sums(0,5) = fSum[15];
 161
 162   smatrix.Invert();
 163   if (!smatrix.IsValid()) return kFALSE;
 164
 165   TMatrixD res = sums*smatrix;
 166   fAlignObj->SetPars(res(0,0),res(0,1),res(0,2),
 167                      TMath::RadToDeg()*res(0,3),
 168                      TMath::RadToDeg()*res(0,4),
 169                      TMath::RadToDeg()*res(0,5));
 170   TMatrixD  tmp = res*sums.T();
 171   fChi2 = fSumR - tmp(0,0);
 172   fNdf -= 6;
 173
 174   return kTRUE;
 175 }