HLT/TPCLib/AliHLTTPCClusterTransformation.cxx

   1 // $Id: AliHLTTPCClusterTransformation.cxx 41244 2010-05-14 08:13:35Z kkanaki $
   2
   3 //**************************************************************************
   4 //* This file is property of and copyright by the ALICE HLT Project        *
   5 //* ALICE Experiment at CERN, All rights reserved.                         *
   6 //*                                                                        *
   7 //* Primary Authors: Kalliopi Kanaki <Kalliopi.Kanaki@ift.uib.no>          *
   8 //*                  for The ALICE HLT Project.                            *
   9 //*                                                                        *
  10 //* Permission to use, copy, modify and distribute this software and its   *
  11 //* documentation strictly for non-commercial purposes is hereby granted   *
  12 //* without fee, provided that the above copyright notice appears in all   *
  13 //* copies and that both the copyright notice and this permission notice   *
  14 //* appear in the supporting documentation. The authors make no claims     *
  15 //* about the suitability of this software for any purpose. It is          *
  16 //* provided "as is" without express or implied warranty.                  *
  17 //**************************************************************************
  18
  19 /** @file   AliHLTTPCClusterTransformation.cxx
  20     @author Kalliopi Kanaki, Sergey Gorbubnov
  21     @date
  22     @brief
  23 */
  24
  25
  26 #include "AliHLTTPCClusterTransformation.h"
  27 #include "AliHLTTPCTransform.h"
  28 #include "AliHLTTPCFastTransform.h"
  29
  30 #include "AliTPCcalibDB.h"
  31 #include "AliTPCTransform.h"
  32 #include "AliTPCParam.h"
  33 #include "AliTPCRecoParam.h"
  34 #include "AliGeomManager.h"
  35 #include <iostream>
  36 #include <iomanip>
  37
  38 using namespace std;
  39
  40 ClassImp(AliHLTTPCClusterTransformation) //ROOT macro for the implementation of ROOT specific class methods
  41
  42 AliHLTTPCClusterTransformation::AliHLTTPCClusterTransformation()
  43 :
  44   fOfflineTPCParam( NULL ),
  45   fLastSector(-1)
  46 {
  47   // see header file for class documentation
  48   // or
  49   // refer to README to build package
  50   // or
  51   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
  52   fAliT[0] = 0.;
  53   fAliT[1] = 0.;
  54   fAliT[2] = 0.;
  55   SetRotationMatrix();
  56 }
  57 AliHLTTPCClusterTransformation::~AliHLTTPCClusterTransformation()
  58 {
  59   // see header file for class documentation
  60   AliHLTTPCFastTransform::Terminate();
  61 }
  62
  63
  64 int  AliHLTTPCClusterTransformation::Init( double FieldBz, UInt_t TimeStamp )
  65 {
  66   // Initialisation
  67
  68   fOfflineTPCParam = 0;
  69
  70   AliTPCcalibDB* pCalib=AliTPCcalibDB::Instance();
  71
  72   if(!pCalib ) return -1;
  73
  74   pCalib->SetExBField(FieldBz);
  75
  76   if(!AliGeomManager::GetGeometry()){
  77      AliGeomManager::LoadGeometry();
  78   }
  79
  80   if( !pCalib->GetTransform() ) return -2;
  81
  82   fOfflineTPCParam = pCalib->GetParameters();
  83   if( !fOfflineTPCParam ) return -3;
  84
  85   fOfflineTPCParam->Update();
  86   fOfflineTPCParam->ReadGeoMatrices();
  87
  88   fLastSector = -1;
  89
  90   fAliT[0] = 0.;
  91   fAliT[1] = 0.;
  92   fAliT[2] = 0.;
  93   SetRotationMatrix();
  94   SetCurrentTimeStamp( TimeStamp );
  95   return 0;
  96 }
  97
  98
  99 void AliHLTTPCClusterTransformation::SetCurrentTimeStamp( UInt_t TimeStamp )
 100 {
 101   // Set the current time stamp
 102   AliHLTTPCFastTransform::Instance()->SetCurrentTimeStamp( TimeStamp );
 103 }
 104
 105
 106 int  AliHLTTPCClusterTransformation::Transform( int Slice, int Row, float Pad, float Time, float XYZ[] )
 107 {
 108   // Convert row, pad, time to X Y Z
 109
 110   Int_t sector=-99, thisrow=-99;
 111   AliHLTTPCTransform::Slice2Sector( Slice, Row, sector, thisrow);
 112   Double_t x[3] = {0,0,0};
 113   AliHLTTPCFastTransform::Instance()->Transform(sector, thisrow, Pad, Time, x);
 114
 115   if( sector!= fLastSector ){
 116     if( fOfflineTPCParam && sector<fOfflineTPCParam->GetNSector() ){
 117       TGeoHMatrix  *alignment = fOfflineTPCParam->GetClusterMatrix( sector );
 118       if ( alignment ){
 119         const Double_t *tr = alignment->GetTranslation();
 120         const Double_t *rot = alignment->GetRotationMatrix();
 121         if( tr && rot ){
 122           for( int i=0; i<3; i++ ) fAliT[i] = tr[i];
 123           for( int i=0; i<9; i++ ) fAliR[i] = rot[i];
 124         }
 125       }
 126     } else {
 127       fAliT[0] = 0.;
 128       fAliT[1] = 0.;
 129       fAliT[2] = 0.;
 130       for( int i=0; i<9; i++ ) fAliR[i] = 0;
 131       fAliR[0] = 1.;
 132       fAliR[4] = 1.;
 133       fAliR[8] = 1.;
 134     }
 135     fLastSector = sector;
 136   }
 137   // alignment->LocalToMaster( x, y);
 138
 139   XYZ[0] = fAliT[0] + x[0]*fAliR[0] + x[1]*fAliR[1] + x[2]*fAliR[2];
 140   XYZ[1] = fAliT[1] + x[0]*fAliR[3] + x[1]*fAliR[4] + x[2]*fAliR[5];
 141   XYZ[2] = fAliT[2] + x[0]*fAliR[6] + x[1]*fAliR[7] + x[2]*fAliR[8];
 142
 143   return 0;
 144 }
 145
 146 int  AliHLTTPCClusterTransformation::ReverseAlignment( float XYZ[], int slice, int padrow)
 147 {
 148   // reverse the alignment correction
 149   Int_t sector=-99, thisrow=-99;
 150   AliHLTTPCTransform::Slice2Sector( slice, padrow, sector, thisrow);
 151   if( sector!= fLastSector ){
 152     if( fOfflineTPCParam && sector<fOfflineTPCParam->GetNSector() ){
 153       TGeoHMatrix  *alignment = fOfflineTPCParam->GetClusterMatrix( sector );
 154       if ( alignment ){
 155         const Double_t *tr = alignment->GetTranslation();
 156         const Double_t *rot = alignment->GetRotationMatrix();
 157         if(tr){
 158           for( int i=0; i<3; i++ ) fAliT[i] = tr[i];
 159         }
 160         SetRotationMatrix(rot, true);
 161       }
 162     } else {
 163       fAliT[0] = 0.;
 164       fAliT[1] = 0.;
 165       fAliT[2] = 0.;
 166       SetRotationMatrix(NULL, true);
 167     }
 168     fLastSector = sector;
 169   }
 170
 171   // correct for alignment: translation
 172   float xyz[3];
 173   xyz[0] = XYZ[0] - fAliT[0];
 174   xyz[1] = XYZ[1] - fAliT[1];
 175   xyz[2] = XYZ[2] - fAliT[2];
 176
 177   // correct for alignment: rotation
 178   XYZ[0]=xyz[0]*fAdjR[0] + xyz[1]*fAdjR[1] + xyz[2]*fAdjR[2];
 179   XYZ[1]=xyz[0]*fAdjR[3] + xyz[1]*fAdjR[4] + xyz[2]*fAdjR[5];
 180   XYZ[2]=xyz[0]*fAdjR[6] + xyz[1]*fAdjR[7] + xyz[2]*fAdjR[8];
 181
 182   return 0;
 183 }
 184
 185 void AliHLTTPCClusterTransformation::SetRotationMatrix(const Double_t *rot, bool bCalcAdjugate)
 186 {
 187   // set the rotation matrix and calculate the adjugate if requested
 188   if (rot) {
 189     for( int i=0; i<9; i++ ) fAliR[i] = rot[i];
 190     if (bCalcAdjugate) {
 191       CalcAdjugateRotation();
 192     }
 193     return;
 194   }
 195   for( int i=0; i<9; i++ ) {fAliR[i] = 0; fAdjR[i] = 0;}
 196   fAliR[0] = 1.;
 197   fAliR[4] = 1.;
 198   fAliR[8] = 1.;
 199   fAdjR[0] = 1.;
 200   fAdjR[4] = 1.;
 201   fAdjR[8] = 1.;
 202 }
 203
 204 bool AliHLTTPCClusterTransformation::CalcAdjugateRotation(bool bCheck)
 205 {
 206   // check rotation matrix and adjugate for consistency
 207   fAdjR[0]= fAliR[4]*fAliR[8]-fAliR[5]*fAliR[7];
 208   fAdjR[1]= fAliR[5]*fAliR[6]-fAliR[3]*fAliR[8];
 209   fAdjR[2]= fAliR[3]*fAliR[7]-fAliR[4]*fAliR[6];
 210
 211   fAdjR[3]= fAliR[2]*fAliR[7]-fAliR[1]*fAliR[8];
 212   fAdjR[4]= fAliR[0]*fAliR[8]-fAliR[2]*fAliR[6];
 213   fAdjR[5]= fAliR[2]*fAliR[6]-fAliR[0]*fAliR[7];
 214
 215   fAdjR[6]= fAliR[1]*fAliR[5]-fAliR[2]*fAliR[4];
 216   fAdjR[7]= fAliR[2]*fAliR[3]-fAliR[0]*fAliR[5];
 217   fAdjR[8]= fAliR[0]*fAliR[4]-fAliR[1]*fAliR[3];
 218
 219   if (bCheck) {
 220     for (int r=0; r<3; r++) {
 221       for (int c=0; c<3; c++) {
 222         float a=0.;
 223         float expected=0.;
 224         if (r==c) expected=1.;
 225         for (int i=0; i<3; i++) {
 226           a+=fAliR[3*r+i]*fAdjR[c+(3*i)];
 227         }
 228         if (TMath::Abs(a-expected)>0.00001) {
 229           std::cout << "inconsistent adjugate at " << r << c << ": " << a << " " << expected << std::endl;
 230           return false;
 231         }
 232       }
 233     }
 234   }
 235   return true;
 236 }
 237
 238 void AliHLTTPCClusterTransformation::Print(const char* /*option*/) const
 239 {
 240   // print info
 241   ios::fmtflags coutflags=std::cout.flags(); // backup cout status flags
 242   std::cout << "AliHLTTPCClusterTransformation for sector " << fLastSector << std::endl;
 243
 244   std::cout.setf(ios_base::showpos|ios_base::showpos|ios::right);
 245   std::cout << "  translation: " << std::endl;
 246   int r=0;
 247   for (r=0; r<3; r++) {
 248     std::cout << setw(7) << fixed << setprecision(2);
 249     cout << "  " << fAliT[r] << std::endl;
 250   }
 251   std::cout << "  rotation and adjugated rotation: " << std::endl;
 252   for (r=0; r<3; r++) {
 253     int c=0;
 254     std::cout << setw(7) << fixed << setprecision(2);
 255     for (c=0; c<3; c++) std::cout << "  " << fAliR[3*r+c];
 256     std::cout << "      ";
 257     for (c=0; c<3; c++) std::cout << "  " << fAdjR[3*r+c];
 258     std::cout << endl;
 259   }
 260   std::cout.flags(coutflags); // restore the original flags
 261 }