STEER/AliAlignObj.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 /* $Id$ */
  17
  18 //-----------------------------------------------------------------
  19 //  Implementation of the alignment object class, holding the alignment
  20 //  constants for a single volume, through the abstract class AliAlignObj.
  21 //  From it two derived concrete representation of alignment object class
  22 //  (AliAlignObjParams, AliAlignObjMatrix) are derived in separate files.
  23 //-----------------------------------------------------------------
  24
  25 #include <TGeoManager.h>
  26 #include <TGeoMatrix.h>
  27 #include <TGeoPhysicalNode.h>
  28 #include <TGeoOverlap.h>
  29 #include <TMath.h>
  30
  31 #include "AliAlignObj.h"
  32 #include "AliTrackPointArray.h"
  33 #include "AliLog.h"
  34
  35 ClassImp(AliAlignObj)
  36
  37 //_____________________________________________________________________________
  38 AliAlignObj::AliAlignObj():
  39   fVolPath(),
  40   fVolUID(0)
  41 {
  42   // default constructor
  43   for(Int_t i=0; i<6; i++) fDiag[i]=-999.;
  44   for(Int_t i=0; i<15; i++) fODia[i]=-999.;
  45 }
  46
  47 //_____________________________________________________________________________
  48 AliAlignObj::AliAlignObj(const char* symname, UShort_t voluid) :
  49   TObject(),
  50   fVolPath(symname),
  51   fVolUID(voluid)
  52 {
  53   // standard constructor
  54   //
  55   for(Int_t i=0; i<6; i++) fDiag[i]=-999.;
  56   for(Int_t i=0; i<15; i++) fODia[i]=-999.;
  57 }
  58
  59 //_____________________________________________________________________________
  60 AliAlignObj::AliAlignObj(const char* symname, UShort_t voluid, Double_t* cmat) :
  61   TObject(),
  62   fVolPath(symname),
  63   fVolUID(voluid)
  64 {
  65   // standard constructor
  66   //
  67   SetCorrMatrix(cmat);
  68 }
  69
  70 //_____________________________________________________________________________
  71 AliAlignObj::AliAlignObj(const AliAlignObj& theAlignObj) :
  72   TObject(theAlignObj),
  73   fVolPath(theAlignObj.GetSymName()),
  74   fVolUID(theAlignObj.GetVolUID())
  75 {
  76   //copy constructor
  77   for(Int_t i=0; i<6; i++) fDiag[i]=theAlignObj.fDiag[i];
  78   for(Int_t i=0; i<15; i++) fODia[i]=theAlignObj.fODia[i];
  79 }
  80
  81 //_____________________________________________________________________________
  82 AliAlignObj &AliAlignObj::operator =(const AliAlignObj& theAlignObj)
  83 {
  84   // assignment operator
  85   if(this==&theAlignObj) return *this;
  86   fVolPath = theAlignObj.GetSymName();
  87   fVolUID = theAlignObj.GetVolUID();
  88   for(Int_t i=0; i<6; i++) fDiag[i]=theAlignObj.fDiag[i];
  89   for(Int_t i=0; i<15; i++) fODia[i]=theAlignObj.fODia[i];
  90   return *this;
  91 }
  92
  93 //_____________________________________________________________________________
  94 AliAlignObj &AliAlignObj::operator*=(const AliAlignObj& theAlignObj)
  95 {
  96   // multiplication operator
  97   // The operator can be used to 'combine'
  98   // two alignment objects
  99   TGeoHMatrix m1;
 100   GetMatrix(m1);
 101   TGeoHMatrix m2;
 102   theAlignObj.GetMatrix(m2);
 103   m1.MultiplyLeft(&m2);
 104   SetMatrix(m1);
 105   // temporary solution: the covariance matrix of the resulting combined object
 106   // is set equal to the covariance matrix of the right operand
 107   // (not to be used for combining alignment objects for different levels)
 108   for(Int_t i=0; i<6; i++)  fDiag[i] = theAlignObj.fDiag[i];
 109   for(Int_t i=0; i<15; i++)  fODia[i] = theAlignObj.fODia[i];
 110   return *this;
 111 }
 112
 113 //_____________________________________________________________________________
 114 AliAlignObj::~AliAlignObj()
 115 {
 116   // dummy destructor
 117 }
 118
 119 //_____________________________________________________________________________
 120 void AliAlignObj::SetVolUID(AliGeomManager::ELayerID detId, Int_t modId)
 121 {
 122   // From detector name and module number (according to detector numbering)
 123   // build fVolUID, unique numerical identity of that volume inside ALICE
 124   // fVolUID is 16 bits, first 5 reserved for detID (32 possible values),
 125   // remaining 11 for module ID inside det (2048 possible values).
 126   //
 127   fVolUID = AliGeomManager::LayerToVolUID(detId,modId);
 128 }
 129
 130 //_____________________________________________________________________________
 131 void AliAlignObj::GetVolUID(AliGeomManager::ELayerID &layerId, Int_t &modId) const
 132 {
 133   // From the fVolUID, unique numerical identity of that volume inside ALICE,
 134   // (voluid is 16 bits, first 5 reserved for layerID (32 possible values),
 135   // remaining 11 for module ID inside det (2048 possible values)), sets
 136   // the argument layerId to the identity of the layer to which that volume
 137   // belongs and sets the argument modId to the identity of that volume
 138   // internally to the layer.
 139   //
 140   layerId = AliGeomManager::VolUIDToLayer(fVolUID,modId);
 141 }
 142
 143 //_____________________________________________________________________________
 144 Bool_t AliAlignObj::GetPars(Double_t tr[], Double_t angles[]) const
 145 {
 146   GetTranslation(tr);
 147   return GetAngles(angles);
 148 }
 149
 150 //_____________________________________________________________________________
 151 Int_t AliAlignObj::GetLevel() const
 152 {
 153   // Return the geometry level of the alignable volume to which
 154   // the alignment object is associated; this is the number of
 155   // slashes in the corresponding volume path
 156   //
 157   if(!gGeoManager){
 158     AliWarning("gGeoManager doesn't exist or it is still open: unable to return meaningful level value.");
 159     return (-1);
 160   }
 161   const char* symname = GetSymName();
 162   const char* path;
 163   TGeoPNEntry* pne = gGeoManager->GetAlignableEntry(symname);
 164   if(pne){
 165     path = pne->GetTitle();
 166   }else{
 167     path = symname;
 168   }
 169
 170   TString pathStr = path;
 171   if(pathStr[0]!='/') pathStr.Prepend('/');
 172   return pathStr.CountChar('/');
 173 }
 174
 175 //_____________________________________________________________________________
 176 Int_t AliAlignObj::Compare(const TObject *obj) const
 177 {
 178   // Compare the levels of two
 179   // alignment objects
 180   // Used in the sorting during
 181   // the application of alignment
 182   // objects to the geometry
 183   //
 184   Int_t level = GetLevel();
 185   Int_t level2 = ((AliAlignObj *)obj)->GetLevel();
 186   if (level == level2)
 187     return 0;
 188   else
 189     return ((level > level2) ? 1 : -1);
 190 }
 191
 192 //______________________________________________________________________________
 193 void AliAlignObj::GetCovMatrix(Double_t *cmat) const
 194 {
 195   // Fills the cmat argument with the coefficients of the external cov matrix (21 elements)
 196   // calculating them from the correlation matrix data member
 197   //
 198
 199   for(Int_t i=0; i<6; ++i) {
 200     // Off diagonal elements
 201     for(Int_t j=0; j<i; ++j) {
 202       cmat[i*(i+1)/2+j] = (fDiag[j] >= 0. && fDiag[i] >= 0.) ? fODia[(i-1)*i/2+j]*fDiag[j]*fDiag[i]: -999.;
 203     }
 204
 205     // Diagonal elements
 206     cmat[i*(i+1)/2+i] = (fDiag[i] >= 0.) ? fDiag[i]*fDiag[i] : -999.;
 207   }
 208
 209   return;
 210 }
 211
 212 //______________________________________________________________________________
 213 void AliAlignObj::GetCovMatrix(TMatrixDSym& mcov) const
 214 {
 215   // Fills the matrix m passed as argument as the covariance matrix calculated
 216   // from the coefficients of the reduced covariance matrix data members
 217   //
 218
 219   for(Int_t i=0; i<6; ++i) {
 220     // Off diagonal elements
 221     for(Int_t j=0; j<i; ++j) {
 222       mcov(j,i) = mcov(i,j) = (fDiag[j] >= 0. && fDiag[i] >= 0.) ? fODia[(i-1)*i/2+j]*fDiag[j]*fDiag[i]: -999.;
 223     }
 224
 225     // Diagonal elements
 226     mcov(i,i) = (fDiag[i] >= 0.) ? fDiag[i]*fDiag[i] : -999.;
 227   }
 228
 229 }
 230
 231 //______________________________________________________________________________
 232 Bool_t AliAlignObj::GetLocalCovMatrix(TMatrixDSym& lCov) const
 233 {
 234   // Calculates the covariance matrix (6x6) associated to the six parameters
 235   // defining the current alignment in the global coordinates system (and sets
 236   // in the internal data members) from the covariance matrix (6x6) for the six
 237   // parameters defining the alignment transformation in the local coordinates
 238   // system, passed as an argument.
 239   //
 240   TMatrixD mJ(6,6);// the jacobian of the transformation from local to global parameters
 241   if(!GetJacobian(mJ)) return kFALSE;
 242
 243   TMatrixDSym gCov(6);
 244   GetCovMatrix(gCov);
 245
 246   // Compute the local covariance matrix lcov = mJ^T gcov mJ
 247   TMatrixD gcovJ(gCov,TMatrixD::kMult,mJ);
 248   TMatrixD lCovM(mJ,TMatrixD::kTransposeMult,gcovJ);
 249   // To be done: somehow check that lCovM is close enough to be symmetric
 250   for(Int_t i=0; i<6; i++)
 251   {
 252     lCov(i,i) = lCovM(i,i);
 253     for(Int_t j=i+1; j<6; j++)
 254     {
 255       lCov(i,j)=lCovM(i,j);
 256       lCov(j,i)=lCovM(i,j);
 257     }
 258   }
 259
 260   return kTRUE;
 261
 262 }
 263
 264 //______________________________________________________________________________
 265 Bool_t AliAlignObj::GetLocalCovMatrix(Double_t *lCov) const
 266 {
 267   // Calculates the covariance matrix (6x6) associated to the six parameters
 268   // defining the current alignment in the global coordinates system (and sets
 269   // in the internal data members) from the covariance matrix (6x6) for the six
 270   // parameters defining the alignment transformation in the local coordinates
 271   // system, passed as an argument.
 272   //
 273   TMatrixDSym lCovMatrix(6);
 274   GetLocalCovMatrix(lCovMatrix);
 275
 276   Int_t k=0;
 277   for(Int_t i=0; i<6; i++)
 278     for(Int_t j=i; j<6; j++)
 279     {
 280        lCov[k++] = lCovMatrix(i,j);
 281     }
 282
 283   return kTRUE;
 284 }
 285
 286 //______________________________________________________________________________
 287 Bool_t AliAlignObj::GetJacobian(TMatrixD& mJ) const
 288 {
 289   // Compute the jacobian J of the transformation of the six local to the six global delta parameters
 290   //
 291   // R00 R01 R02 | (R01Rk2 - R02Rk1)Tk  (R02Rk0 - R00Rk2)Tk  (R00Rk1 - R01Rk0)Tk
 292   // R00 R01 R02 | (R11Rk2 - R12Rk1)Tk  (R12Rk0 - R10Rk2)Tk  (R10Rk1 - R11Rk0)Tk
 293   // R00 R01 R02 | (R21Rk2 - R22Rk1)Tk  (R22Rk0 - R20Rk2)Tk  (R20Rk1 - R21Rk0)Tk
 294   //  -  -   -   -   -   -   -  -  -  -  -  -  -  -  -  -  -  -  -  -  -  -  -
 295   //  0   0   0  |   R11R22 - R12R21      R12R20 - R10R22      R10R21 - R11R20
 296   //  0   0   0  |   R21R02 - R22R01      R22R00 - R20R02      R20R01 - R21R00
 297   //  0   0   0  |   R01R12 - R02R11      R02R10 - R00R12      R00R11 - R01R10
 298   //
 299   if (!gGeoManager || !gGeoManager->IsClosed()) {
 300     AliError("Can't compute the global covariance matrix from the local one without an open geometry!");
 301     return kFALSE;
 302   }
 303
 304   const char* symname = GetSymName();
 305   TGeoPhysicalNode* node;
 306   TGeoPNEntry* pne = gGeoManager->GetAlignableEntry(symname);
 307   if(pne){
 308     if(!pne->GetPhysicalNode()){
 309       node = gGeoManager->MakeAlignablePN(pne);
 310     }else{
 311       node = pne->GetPhysicalNode();
 312     }
 313   }else{
 314     AliWarning(Form("The symbolic volume name %s does not correspond to a physical entry. Using it as volume path!",symname));
 315     node = (TGeoPhysicalNode*) gGeoManager->MakePhysicalNode(symname);
 316   }
 317
 318   if (!node) {
 319     AliError(Form("Volume name or path %s not valid!",symname));
 320     return kFALSE;
 321   }
 322
 323   TGeoHMatrix gm; //global matrix
 324   gm = *node->GetMatrix();
 325   Double_t *tr  = gm.GetTranslation();
 326   Double_t *rot = gm.GetRotationMatrix();
 327
 328   TGeoHMatrix m; // global delta transformation matrix
 329   GetMatrix(m);
 330   // We should probably check that it's sufficinetly close to identity
 331   // if it's not return because the "small angles" approximation cannot hold
 332
 333   // 3x3 upper left part (global shifts derived w.r.t. local shifts)
 334   for(Int_t i=0; i<3; i++)
 335   {
 336     for(Int_t j=0; j<3; j++)
 337     {
 338       mJ(i,j) = rot[i+3*j];
 339     }
 340   }
 341
 342   // 3x3 lower left part (global angles derived w.r.t. local shifts)
 343   for(Int_t i=0; i<3; i++)
 344   {
 345     for(Int_t j=0; j<3; j++)
 346     {
 347       mJ(i+3,j) = 0.;
 348     }
 349   }
 350
 351   // 3x3 upper right part (global shifts derived w.r.t. local angles)
 352   for(Int_t i=0; i<3; i++)
 353   {
 354     for(Int_t j=0; j<3; j++)
 355     {
 356       Double_t mEl = 0.;
 357       Int_t b = (j+1)%3;
 358       Int_t d = (j+2)%3;
 359       for(Int_t k=0; k<3; k++)
 360       {
 361         mEl += (rot[3*i+b]*rot[3*k+d])*tr[k]-(rot[3*i+d]*rot[3*k+b])*tr[k];
 362       }
 363       mJ(i,j+3) = mEl;
 364     }
 365   }
 366
 367   // 3x3 lower right part (global angles derived w.r.t. local angles)
 368   for(Int_t i=0; i<3; i++)
 369     for(Int_t j=0; j<3; j++)
 370     {
 371       Int_t a = (i+1)%3;
 372       Int_t b = (j+1)%3;
 373       Int_t c = (i+2)%3;
 374       Int_t d = (j+2)%3;
 375       mJ(i+3,j+3) = rot[3*a+b]*rot[3*c+d]-rot[3*a+d]*rot[3*c+b];
 376     }
 377
 378   return kTRUE;
 379
 380 }
 381
 382 //______________________________________________________________________________
 383 Bool_t AliAlignObj::SetFromLocalCov(TMatrixDSym& lCov)
 384 {
 385   // Calculates the covariance matrix (6x6) associated to the six parameters
 386   // defining the current alignment in the global coordinates system (and sets
 387   // in the internal data members) from the covariance matrix (6x6) for the six
 388   // parameters defining the alignment transformation in the local coordinates
 389   // system, passed as an argument.
 390   //
 391   TMatrixD mJ(6,6);// the jacobian of the transformation from local to global parameters
 392   if(!GetJacobian(mJ)) return kFALSE;
 393
 394   // Compute the global covariance matrix gcov = mJ lcov mJ'
 395   TMatrixD trJ(TMatrixD::kTransposed, mJ);
 396   TMatrixD lcovTrJ(lCov,TMatrixD::kMult,trJ);
 397   TMatrixD gCovM(mJ,TMatrixD::kMult,lcovTrJ);
 398   // To be done: somehow check that gCovM is close enough to be symmetric
 399   TMatrixDSym gCov(6);
 400   for(Int_t i=0; i<6; i++)
 401   {
 402     gCov(i,i) = gCovM(i,i);
 403     for(Int_t j=i+1; j<6; j++)
 404     {
 405       gCov(i,j)=gCovM(i,j);
 406       gCov(j,i)=gCovM(i,j);
 407     }
 408   }
 409   SetCorrMatrix(gCov);
 410
 411   return kTRUE;
 412
 413 }
 414
 415 //______________________________________________________________________________
 416 Bool_t AliAlignObj::SetFromLocalCov(Double_t *lCov)
 417 {
 418   // Calculates the covariance matrix (6x6) associated to the six parameters
 419   // defining the current alignment in the global coordinates system, and sets
 420   // in the internal data members, from the 21 coefficients, passed as argument,
 421   // of the covariance matrix (6x6) for the six parameters defining the
 422   // alignment transformation in the local coordinates system.
 423   //
 424   TMatrixDSym lCovMatrix(6);
 425
 426   Int_t k=0;
 427   for(Int_t i=0; i<6; i++)
 428     for(Int_t j=i; j<6; j++)
 429     {
 430       lCovMatrix(i,j) = lCov[k++];
 431       if(j!=i) lCovMatrix(j,i) = lCovMatrix(i,j);
 432     }
 433
 434   return SetFromLocalCov(lCovMatrix);
 435
 436 }
 437
 438
 439 //______________________________________________________________________________
 440 void AliAlignObj::SetCorrMatrix(Double_t *cmat)
 441 {
 442   // Sets the correlation matrix data member from the coefficients of the external covariance
 443   // matrix (21 elements passed as argument).
 444   //
 445   if(cmat) {
 446
 447     // Diagonal elements first
 448     for(Int_t i=0; i<6; ++i) {
 449       fDiag[i] = (cmat[i*(i+1)/2+i] >= 0.) ? TMath::Sqrt(cmat[i*(i+1)/2+i]) : -999.;
 450     }
 451
 452     // ... then the ones off diagonal
 453     for(Int_t i=0; i<6; ++i)
 454       // Off diagonal elements
 455       for(Int_t j=0; j<i; ++j) {
 456         fODia[(i-1)*i/2+j] = (fDiag[i] > 0. && fDiag[j] > 0.) ? cmat[i*(i+1)/2+j]/(fDiag[j]*fDiag[i]) : 0.;       // check for division by zero (due to diagonal element of 0) and for fDiag != -999. (due to negative input diagonal element).
 457         if (fODia[(i-1)*i/2+j]>1.)  fODia[(i-1)*i/2+j] =  1.; // check upper boundary
 458         if (fODia[(i-1)*i/2+j]<-1.) fODia[(i-1)*i/2+j] = -1.; // check lower boundary
 459       }
 460   } else {
 461     for(Int_t i=0; i< 6; ++i) fDiag[i]=-999.;
 462     for(Int_t i=0; i< 6*(6-1)/2; ++i) fODia[i]=0.;
 463   }
 464
 465   return;
 466 }
 467
 468 //______________________________________________________________________________
 469 void AliAlignObj::SetCorrMatrix(TMatrixDSym& mcov)
 470 {
 471   // Sets the correlation matrix data member from the covariance matrix mcov passed
 472   // passed as argument.
 473   //
 474   if(mcov.IsValid()) {
 475
 476     // Diagonal elements first
 477     for(Int_t i=0; i<6; ++i) {
 478       fDiag[i] = (mcov(i,i) >= 0.) ? TMath::Sqrt(mcov(i,i)) : -999.;
 479     }
 480
 481     // ... then the ones off diagonal
 482     for(Int_t i=0; i<6; ++i)
 483       // Off diagonal elements
 484       for(Int_t j=0; j<i; ++j) {
 485         fODia[(i-1)*i/2+j] = (fDiag[i] > 0. && fDiag[j] > 0.) ? mcov(i,j)/(fDiag[j]*fDiag[i]) : 0.;       // check for division by zero (due to diagonal element of 0) and for fDiag != -999. (due to negative input diagonal element).
 486         if (fODia[(i-1)*i/2+j]>1.)  fODia[(i-1)*i/2+j] =  1.; // check upper boundary
 487         if (fODia[(i-1)*i/2+j]<-1.) fODia[(i-1)*i/2+j] = -1.; // check lower boundary
 488       }
 489   } else {
 490     for(Int_t i=0; i< 6; ++i) fDiag[i]=-999.;
 491     for(Int_t i=0; i< 6*(6-1)/2; ++i) fODia[i]=0.;
 492   }
 493
 494   return;
 495 }
 496
 497 //_____________________________________________________________________________
 498 void AliAlignObj::AnglesToMatrix(const Double_t *angles, Double_t *rot) const
 499 {
 500   // Calculates the rotation matrix using the
 501   // Euler angles in "x y z" notation
 502   //
 503   Double_t degrad = TMath::DegToRad();
 504   Double_t sinpsi = TMath::Sin(degrad*angles[0]);
 505   Double_t cospsi = TMath::Cos(degrad*angles[0]);
 506   Double_t sinthe = TMath::Sin(degrad*angles[1]);
 507   Double_t costhe = TMath::Cos(degrad*angles[1]);
 508   Double_t sinphi = TMath::Sin(degrad*angles[2]);
 509   Double_t cosphi = TMath::Cos(degrad*angles[2]);
 510
 511   rot[0] =  costhe*cosphi;
 512   rot[1] = -costhe*sinphi;
 513   rot[2] =  sinthe;
 514   rot[3] =  sinpsi*sinthe*cosphi + cospsi*sinphi;
 515   rot[4] = -sinpsi*sinthe*sinphi + cospsi*cosphi;
 516   rot[5] = -costhe*sinpsi;
 517   rot[6] = -cospsi*sinthe*cosphi + sinpsi*sinphi;
 518   rot[7] =  cospsi*sinthe*sinphi + sinpsi*cosphi;
 519   rot[8] =  costhe*cospsi;
 520 }
 521
 522 //_____________________________________________________________________________
 523 Bool_t AliAlignObj::MatrixToAngles(const Double_t *rot, Double_t *angles) const
 524 {
 525   // Calculates the Euler angles in "x y z" notation
 526   // using the rotation matrix
 527   // Returns false in case the rotation angles can not be
 528   // extracted from the matrix
 529   //
 530   if(TMath::Abs(rot[0])<1e-7 || TMath::Abs(rot[8])<1e-7) {
 531     AliError("Failed to extract roll-pitch-yall angles!");
 532     return kFALSE;
 533   }
 534   Double_t raddeg = TMath::RadToDeg();
 535   angles[0]=raddeg*TMath::ATan2(-rot[5],rot[8]);
 536   angles[1]=raddeg*TMath::ASin(rot[2]);
 537   angles[2]=raddeg*TMath::ATan2(-rot[1],rot[0]);
 538   return kTRUE;
 539 }
 540
 541 //______________________________________________________________________________
 542 void AliAlignObj::Transform(AliTrackPoint &p, Bool_t copycov) const
 543 {
 544   // The method transforms the space-point coordinates using the
 545   // transformation matrix provided by the AliAlignObj
 546   // In case the copycov flag is set to kTRUE, the covariance matrix
 547   // of the alignment object is copied into the space-point
 548   //
 549   if (fVolUID != p.GetVolumeID())
 550     AliWarning(Form("Alignment object ID is not equal to the space-point ID (%d != %d)",fVolUID,p.GetVolumeID()));
 551
 552   TGeoHMatrix m;
 553   GetMatrix(m);
 554   Double_t *rot = m.GetRotationMatrix();
 555   Double_t *tr  = m.GetTranslation();
 556
 557   Float_t xyzin[3],xyzout[3];
 558   p.GetXYZ(xyzin);
 559   for (Int_t i = 0; i < 3; i++)
 560     xyzout[i] = tr[i]+
 561                 xyzin[0]*rot[3*i]+
 562                 xyzin[1]*rot[3*i+1]+
 563                 xyzin[2]*rot[3*i+2];
 564   p.SetXYZ(xyzout);
 565
 566   if(copycov){
 567     TMatrixDSym covmat(6);
 568     GetCovMatrix(covmat);
 569     p.SetAlignCovMatrix(covmat);
 570   }
 571
 572 }
 573
 574 //_____________________________________________________________________________
 575 void AliAlignObj::Transform(AliTrackPointArray &array) const
 576 {
 577   // This method is used to transform all the track points
 578   // from the input AliTrackPointArray
 579   //
 580   AliTrackPoint p;
 581   for (Int_t i = 0; i < array.GetNPoints(); i++) {
 582     array.GetPoint(p,i);
 583     Transform(p);
 584     array.AddPoint(i,&p);
 585   }
 586 }
 587
 588 //_____________________________________________________________________________
 589 void AliAlignObj::Print(Option_t *) const
 590 {
 591   // Print the contents of the
 592   // alignment object in angles and
 593   // matrix representations
 594   //
 595   Double_t tr[3];
 596   GetTranslation(tr);
 597   Double_t angles[3];
 598   GetAngles(angles);
 599   TGeoHMatrix m;
 600   GetMatrix(m);
 601   const Double_t *rot = m.GetRotationMatrix();
 602
 603   printf("Volume=%s\n",GetSymName());
 604   if (GetVolUID() != 0) {
 605     AliGeomManager::ELayerID layerId;
 606     Int_t modId;
 607     GetVolUID(layerId,modId);
 608     printf("VolumeID=%d LayerID=%d ( %s ) ModuleID=%d\n", GetVolUID(),layerId,AliGeomManager::LayerName(layerId),modId);
 609   }
 610   printf("%12.8f%12.8f%12.8f    Tx = %12.8f    Psi   = %12.8f\n", rot[0], rot[1], rot[2], tr[0], angles[0]);
 611   printf("%12.8f%12.8f%12.8f    Ty = %12.8f    Theta = %12.8f\n", rot[3], rot[4], rot[5], tr[1], angles[1]);
 612   printf("%12.8f%12.8f%12.8f    Tz = %12.8f    Phi   = %12.8f\n", rot[6], rot[7], rot[8], tr[2], angles[2]);
 613
 614 }
 615
 616 //_____________________________________________________________________________
 617 void AliAlignObj::SetPars(Double_t x, Double_t y, Double_t z,
 618                           Double_t psi, Double_t theta, Double_t phi)
 619 {
 620   // Set the global delta transformation by passing 3 angles (expressed in
 621   // degrees) and 3 shifts (in centimeters)
 622   //
 623   SetTranslation(x,y,z);
 624   SetRotation(psi,theta,phi);
 625 }
 626
 627 //_____________________________________________________________________________
 628 Bool_t AliAlignObj::SetLocalPars(Double_t x, Double_t y, Double_t z,
 629                                  Double_t psi, Double_t theta, Double_t phi)
 630 {
 631   // Set the global delta transformation by passing the parameters
 632   // for the local delta transformation (3 shifts and 3 angles).
 633   // In case that the TGeo was not initialized or not closed,
 634   // returns false and the object parameters are not set.
 635   //
 636   TGeoHMatrix m;
 637   Double_t tr[3] = {x, y, z};
 638   m.SetTranslation(tr);
 639   Double_t angles[3] = {psi, theta, phi};
 640   Double_t rot[9];
 641   AnglesToMatrix(angles,rot);
 642   m.SetRotation(rot);
 643
 644   return SetLocalMatrix(m);
 645
 646 }
 647
 648 //_____________________________________________________________________________
 649 Bool_t AliAlignObj::SetLocalTranslation(Double_t x, Double_t y, Double_t z)
 650 {
 651   // Set the global delta transformation by passing the three shifts giving
 652   // the translation in the local reference system of the alignable
 653   // volume (known by TGeo geometry).
 654   // In case that the TGeo was not initialized or not closed,
 655   // returns false and the object parameters are not set.
 656   //
 657   TGeoHMatrix m;
 658   Double_t tr[3] = {x, y, z};
 659   m.SetTranslation(tr);
 660
 661   return SetLocalMatrix(m);
 662
 663 }
 664
 665 //_____________________________________________________________________________
 666 Bool_t AliAlignObj::SetLocalTranslation(const TGeoMatrix& m)
 667 {
 668   // Set the global delta transformation by passing the matrix of
 669   // the local delta transformation and taking its translational part
 670   // In case that the TGeo was not initialized or not closed,
 671   // returns false and the object parameters are not set.
 672   //
 673   const Double_t* tr = m.GetTranslation();
 674   TGeoHMatrix mtr;
 675   mtr.SetTranslation(tr);
 676
 677   return SetLocalMatrix(mtr);
 678
 679 }
 680
 681 //_____________________________________________________________________________
 682 Bool_t AliAlignObj::SetLocalRotation(Double_t psi, Double_t theta, Double_t phi)
 683 {
 684   // Set the global delta transformation by passing the three angles giving
 685   // the rotation in the local reference system of the alignable
 686   // volume (known by TGeo geometry).
 687   // In case that the TGeo was not initialized or not closed,
 688   // returns false and the object parameters are not set.
 689   //
 690   TGeoHMatrix m;
 691   Double_t angles[3] = {psi, theta, phi};
 692   Double_t rot[9];
 693   AnglesToMatrix(angles,rot);
 694   m.SetRotation(rot);
 695
 696   return SetLocalMatrix(m);
 697
 698 }
 699
 700 //_____________________________________________________________________________
 701 Bool_t AliAlignObj::SetLocalRotation(const TGeoMatrix& m)
 702 {
 703   // Set the global delta transformation by passing the matrix of
 704   // the local delta transformation and taking its rotational part
 705   // In case that the TGeo was not initialized or not closed,
 706   // returns false and the object parameters are not set.
 707   //
 708   TGeoHMatrix rotm;
 709   const Double_t* rot = m.GetRotationMatrix();
 710   rotm.SetRotation(rot);
 711
 712   return SetLocalMatrix(rotm);
 713
 714 }
 715
 716 //_____________________________________________________________________________
 717 Bool_t AliAlignObj::SetLocalMatrix(const TGeoMatrix& m)
 718 {
 719   // Set the global delta transformation by passing the TGeo matrix
 720   // for the local delta transformation.
 721   // In case that the TGeo was not initialized or not closed,
 722   // returns false and the object parameters are not set.
 723   //
 724   if (!gGeoManager || !gGeoManager->IsClosed()) {
 725     AliError("Can't set the local alignment object parameters! gGeoManager doesn't exist or it is still open!");
 726     return kFALSE;
 727   }
 728
 729   const char* symname = GetSymName();
 730   TGeoHMatrix gprime,gprimeinv;
 731   TGeoPhysicalNode* pn = 0;
 732   TGeoPNEntry* pne = gGeoManager->GetAlignableEntry(symname);
 733   if(pne)
 734   {
 735     pn = pne->GetPhysicalNode();
 736     if(pn){
 737       if (pn->IsAligned())
 738         AliWarning(Form("Volume %s has been misaligned already!",symname));
 739       gprime = *pn->GetMatrix();
 740     }else{
 741       gprime = pne->GetGlobalOrig();
 742     }
 743   }else{
 744     AliWarning(Form("The symbolic volume name %s does not correspond to a physical entry. Using it as volume path!",symname));
 745     if(!gGeoManager->cd(symname)) {
 746       AliError(Form("Volume name or path %s not valid!",symname));
 747       return kFALSE;
 748     }
 749     gprime = *gGeoManager->GetCurrentMatrix();
 750   }
 751
 752   TGeoHMatrix m1; // the TGeoHMatrix copy of the local delta "m"
 753   const Double_t *tr = m.GetTranslation();
 754   m1.SetTranslation(tr);
 755   const Double_t* rot = m.GetRotationMatrix();
 756   m1.SetRotation(rot);
 757
 758   gprimeinv = gprime.Inverse();
 759   m1.Multiply(&gprimeinv);
 760   m1.MultiplyLeft(&gprime);
 761
 762   return SetMatrix(m1);
 763 }
 764
 765 //_____________________________________________________________________________
 766 Bool_t AliAlignObj::SetMatrix(const TGeoMatrix& m)
 767 {
 768   // Set the global delta transformation by passing the TGeoMatrix
 769   // for it
 770   //
 771   SetTranslation(m);
 772   return SetRotation(m);
 773 }
 774
 775 //_____________________________________________________________________________
 776 Bool_t AliAlignObj::GetLocalPars(Double_t transl[], Double_t angles[]) const
 777 {
 778   // Get the translations and angles (in degrees) expressing the
 779   // local delta transformation.
 780   // In case that the TGeo was not initialized or not closed,
 781   // returns false and the object parameters are not set.
 782   //
 783   if(!GetLocalTranslation(transl)) return kFALSE;
 784   return GetLocalAngles(angles);
 785 }
 786
 787 //_____________________________________________________________________________
 788 Bool_t AliAlignObj::GetLocalTranslation(Double_t* tr) const
 789 {
 790   // Get the 3 shifts giving the translational part of the local
 791   // delta transformation.
 792   // In case that the TGeo was not initialized or not closed,
 793   // returns false and the object parameters are not set.
 794   //
 795   TGeoHMatrix ml;
 796   if(!GetLocalMatrix(ml)) return kFALSE;
 797   const Double_t* transl;
 798   transl = ml.GetTranslation();
 799   tr[0]=transl[0];
 800   tr[1]=transl[1];
 801   tr[2]=transl[2];
 802   return kTRUE;
 803 }
 804
 805 //_____________________________________________________________________________
 806 Bool_t AliAlignObj::GetLocalAngles(Double_t* angles) const
 807 {
 808   // Get the 3 angles giving the rotational part of the local
 809   // delta transformation.
 810   // In case that the TGeo was not initialized or not closed,
 811   // returns false and the object parameters are not set.
 812   //
 813   TGeoHMatrix ml;
 814   if(!GetLocalMatrix(ml)) return kFALSE;
 815   const Double_t *rot = ml.GetRotationMatrix();
 816   return MatrixToAngles(rot,angles);
 817 }
 818
 819 //_____________________________________________________________________________
 820 Bool_t AliAlignObj::GetLocalMatrix(TGeoHMatrix& m) const
 821 {
 822   // Get the matrix for the local delta transformation.
 823   // In case that the TGeo was not initialized or not closed,
 824   // returns false and the object parameters are not set.
 825   //
 826   if (!gGeoManager || !gGeoManager->IsClosed()) {
 827     AliError("Can't get the local alignment object parameters! gGeoManager doesn't exist or it is still open!");
 828     return kFALSE;
 829   }
 830
 831   const char* symname = GetSymName();
 832   TGeoPhysicalNode* node;
 833   TGeoPNEntry* pne = gGeoManager->GetAlignableEntry(symname);
 834   if(pne){
 835     if(!pne->GetPhysicalNode()){
 836       node = gGeoManager->MakeAlignablePN(pne);
 837     }else{
 838       node = pne->GetPhysicalNode();
 839     }
 840   }else{
 841     AliWarning(Form("The symbolic volume name %s does not correspond to a physical entry. Using it as volume path!",symname));
 842     node = (TGeoPhysicalNode*) gGeoManager->MakePhysicalNode(symname);
 843   }
 844
 845   if (!node) {
 846     AliError(Form("Volume name or path %s not valid!",symname));
 847     return kFALSE;
 848   }
 849 //  if (node->IsAligned())
 850 //    AliWarning(Form("Volume %s has been misaligned already!",symname));
 851
 852   GetMatrix(m);
 853   TGeoHMatrix gprime,gprimeinv;
 854   gprime = *node->GetMatrix();
 855   gprimeinv = gprime.Inverse();
 856   m.Multiply(&gprime);
 857   m.MultiplyLeft(&gprimeinv);
 858
 859   return kTRUE;
 860 }
 861
 862 //_____________________________________________________________________________
 863 Bool_t AliAlignObj::ApplyToGeometry(Bool_t ovlpcheck)
 864 {
 865   // Apply the current alignment object to the TGeo geometry
 866   // This method returns FALSE if the symname of the object was not
 867   // valid neither to get a TGeoPEntry nor as a volume path
 868   //
 869   if (!gGeoManager || !gGeoManager->IsClosed()) {
 870     AliError("Can't apply the alignment object! gGeoManager doesn't exist or it is still open!");
 871     return kFALSE;
 872   }
 873
 874   if (gGeoManager->IsLocked()){
 875     AliError("Can't apply the alignment object! Geometry is locked!");
 876     return kFALSE;
 877   }
 878
 879   const char* symname = GetSymName();
 880   const char* path;
 881   TGeoPhysicalNode* node;
 882   TGeoPNEntry* pne = gGeoManager->GetAlignableEntry(symname);
 883   if(pne){
 884     path = pne->GetTitle();
 885     node = gGeoManager->MakeAlignablePN(pne);
 886   }else{
 887     AliDebug(1,Form("The symbolic volume name %s does not correspond to a physical entry. Using it as a volume path!",symname));
 888     path=symname;
 889     if (!gGeoManager->CheckPath(path)) {
 890       AliDebug(1,Form("Volume path %s not valid!",path));
 891       return kFALSE;
 892     }
 893     if (gGeoManager->GetListOfPhysicalNodes()->FindObject(path)) {
 894       AliError(Form("Volume %s has been misaligned already!",path));
 895       return kFALSE;
 896     }
 897     node = (TGeoPhysicalNode*) gGeoManager->MakePhysicalNode(path);
 898   }
 899
 900   if (!node) {
 901     AliError(Form("Volume path %s not valid!",path));
 902     return kFALSE;
 903   }
 904
 905   //  Double_t threshold = 0.001;
 906
 907   TGeoHMatrix align,gprime;
 908   gprime = *node->GetMatrix();
 909   GetMatrix(align);
 910   gprime.MultiplyLeft(&align);
 911   TGeoHMatrix *ginv = new TGeoHMatrix;
 912   TGeoHMatrix *g = node->GetMatrix(node->GetLevel()-1);
 913   *ginv = g->Inverse();
 914   *ginv *= gprime;
 915   AliGeomManager::ELayerID layerId; // unique identity for layer in the alobj
 916   Int_t modId; // unique identity for volume inside layer in the alobj
 917   GetVolUID(layerId, modId);
 918   AliDebug(2,Form("Aligning volume %s of detector layer %d with local ID %d",symname,layerId,modId));
 919   if(ovlpcheck){
 920     node->Align(ginv,0,kTRUE); //(trunk of root takes threshold as additional argument)
 921   }else{
 922     node->Align(ginv,0,kFALSE);
 923   }
 924   if(ovlpcheck)
 925   {
 926     TObjArray* ovlpArray =  gGeoManager->GetListOfOverlaps();
 927     Int_t nOvlp = ovlpArray->GetEntriesFast();
 928     if(nOvlp)
 929     {
 930       AliInfo(Form("Misalignment of node %s generated the following overlaps/extrusions:",node->GetName()));
 931       for(Int_t i=0; i<nOvlp; i++)
 932         ((TGeoOverlap*)ovlpArray->UncheckedAt(i))->PrintInfo();
 933     }
 934   }
 935
 936   return kTRUE;
 937 }
 938
 939