TPC/TPCbase/AliTPCCalibGlobalMisalignment.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 ////////////////////////////////////////////////////////////////////////////
  18 //                                                                        //
  19 // AliTPCCalibGlobalMisalignment class                                        //
  20 // The class calculates the space point distortions due to simple         //
  21 // misalignments like shifts in caresian coordinates or a rotation        //
  22 // of the TPC read out planes (A and C side)                              //
  23 // Optionaly possible to use it for visualization of the alignemnt form the Alignment OCDB //
  24 // fUseGeoManager has to be set to kTRUE to enable this option
  25 //                                                                        //
  26 // date: 06/05/2010                                                       //
  27 // Authors: Stefan Rossegger, Jim Thomas, Magnus Mager                    //
  28 ////////////////////////////////////////////////////////////////////////////
  29
  30 #include "AliTPCCalibGlobalMisalignment.h"
  31 #include "TMath.h"
  32 #include "TGeoMatrix.h"
  33 #include "AliTPCROC.h"
  34 #include "AliTPCcalibDB.h"
  35 #include "AliTPCParam.h"
  36 #include <TGeoPhysicalNode.h>
  37 //
  38 #include "AliAlignObjParams.h"
  39 #include "AliGeomManager.h"
  40 #include "AliCDBManager.h"
  41 #include "AliCDBEntry.h"
  42
  43 AliTPCCalibGlobalMisalignment::AliTPCCalibGlobalMisalignment()
  44   : AliTPCCorrection("mialign","Misalignment"),
  45     fXShift(0.),fYShift(0.),fZShift(0.),
  46     fRotPhiA(0.),fRotPhiC(0.),
  47     fdRPhiOffsetA(0.),
  48     fdRPhiOffsetC(0.),
  49     fQuadrantQ0(0),   //OROC medium pads -delta ly+ - ly - shift
  50     fQuadrantRQ0(0),  //OROC medium pads -delta ly+ - ly - rotation
  51     fQuadrantQ1(0),   //OROC long   pads -delta ly+ - ly - shift
  52     fQuadrantQ2(0),   //OROC long   pads -shift
  53     fQuadrantRQ1(0),  //OROC long   pads -delta ly+ - ly - rotation
  54     fQuadrantRQ2(0),  //OROC long   pads -rotation
  55     fMatrixGlobal(0), // global Alignment common
  56     fMatrixGlobalDelta(0), // global Alignment Delta A side-c side
  57     fArraySector(0)   // fArraySector
  58 {
  59   //
  60   // default constructor
  61   //
  62 }
  63
  64 AliTPCCalibGlobalMisalignment::~AliTPCCalibGlobalMisalignment() {
  65   //
  66   // destructor
  67   //
  68   delete fQuadrantQ0;   //OROC medium pads -delta ly+ - ly - shift
  69   delete fQuadrantRQ0;  //OROC medium pads -delta ly+ - ly - rotation
  70   delete fQuadrantQ1;   //OROC long   pads -delta ly+ - ly - shift
  71   delete fQuadrantQ2;   //OROC long   pads -shift
  72   delete fQuadrantRQ1;  //OROC long   pads -delta ly+ - ly - rotation
  73   delete fQuadrantRQ2;  //OROC long   pads -rotation
  74   delete fMatrixGlobal; // global matrix
  75   delete fMatrixGlobal; // global matrix
  76   delete fArraySector;  // sector matrices
  77 }
  78
  79
  80
  81 Bool_t AliTPCCalibGlobalMisalignment::AddCorrectionCompact(AliTPCCorrection* corr, Double_t weight){
  82   //
  83   // Add correction  and make them compact
  84   // Assumptions:
  85   //  - origin of distortion/correction are additive
  86   //  - only correction ot the same type supported ()
  87   if (corr==NULL) {
  88     //AliError("Zerro pointer - correction");
  89     return kFALSE;
  90   }
  91   AliTPCCalibGlobalMisalignment* corrC = dynamic_cast<AliTPCCalibGlobalMisalignment *>(corr);
  92   if (corrC == NULL) {
  93     //AliError(TString::Format("Inconsistent class types: %s\%s",IsA()->GetName(),corr->IsA()->GetName()).Data());
  94     return kFALSE;
  95   }
  96   //
  97   AliTPCCalibGlobalMisalignment & add = *corrC;
  98   fXShift+=weight*add.fXShift;               // Shift in global X [cm]
  99   fYShift+=weight*add.fYShift;               // Shift in global Y [cm]
 100   fZShift+=weight*add.fZShift;               // Shift in global Z [cm]
 101
 102   fRotPhiA+=weight*add.fRotPhiA;      // simple rotation of A side read-out plane around the Z axis [rad]
 103   fRotPhiC+=weight*add.fRotPhiC;      // simple rotation of C side read-out plane around the Z axis [rad]
 104   fdRPhiOffsetA+=weight*add.fdRPhiOffsetA;  // add a constant offset of dRPhi (or local Y) in [cm]: purely for calibration purposes!
 105   fdRPhiOffsetC+=weight*add.fdRPhiOffsetC;  // add a constant offset of dRPhi (or local Y) in [cm]: purely for calibration purposes!
 106   //
 107   // Quadrant alignment
 108   //
 109   if (add.fQuadrantQ0) {
 110     if (fQuadrantQ0)  fQuadrantQ0->Add(weight*(*(add.fQuadrantQ0)));
 111     if (!fQuadrantQ0) {
 112       fQuadrantQ0 = (TVectorD*)(add.fQuadrantQ0->Clone());
 113       (*fQuadrantQ0)*=weight;
 114     }
 115   }
 116   if (add.fQuadrantRQ0) {
 117     if (fQuadrantRQ0) fQuadrantRQ0->Add(weight*(*(add.fQuadrantRQ0)));
 118     if (!fQuadrantRQ0) {
 119       fQuadrantRQ0 = (TVectorD*)(add.fQuadrantRQ0->Clone());
 120       (*fQuadrantRQ0)*=weight;
 121     }
 122   }
 123   //
 124   if (add.fQuadrantQ1) {
 125     if (fQuadrantQ1) fQuadrantQ1->Add(weight*(*(add.fQuadrantQ1)));
 126     if (!fQuadrantQ1) {
 127       fQuadrantQ1 = (TVectorD*)(add.fQuadrantQ1->Clone());
 128        (*fQuadrantQ1)*=weight;
 129     }
 130   }
 131   if (add.fQuadrantRQ1) {
 132     if (fQuadrantRQ1) fQuadrantRQ1->Add(weight*(*(add.fQuadrantRQ1)));
 133     if (!fQuadrantRQ1) {
 134       fQuadrantRQ1 = (TVectorD*)(add.fQuadrantRQ1->Clone());
 135        (*fQuadrantRQ1)*=weight;
 136     }
 137   }
 138   //
 139   if (add.fQuadrantQ2) {
 140     if (fQuadrantQ2) fQuadrantQ2->Add(weight*(*(add.fQuadrantQ2)));
 141     if (!fQuadrantQ2) {
 142       fQuadrantQ2 = (TVectorD*)(add.fQuadrantQ2->Clone());
 143       (*fQuadrantQ2)*=weight;
 144     }
 145   }
 146   if (add.fQuadrantRQ2) {
 147     if (fQuadrantRQ2) fQuadrantRQ2->Add(weight*(*(add.fQuadrantRQ2)));
 148     if (!fQuadrantRQ2) {
 149       fQuadrantRQ2 = (TVectorD*)(add.fQuadrantRQ2->Clone());
 150       (*fQuadrantQ2)*=weight;
 151     }
 152   }
 153   //
 154   // Global alignment - use native ROOT representation
 155   //
 156   Double_t delta[3]={0};
 157   if (add.fMatrixGlobal){
 158     TGeoHMatrix matrixW=*(add.fMatrixGlobal);
 159     TGeoHMatrix matrixScaled;
 160     const Double_t *rotMatrix   = matrixW.GetRotationMatrix();
 161     const Double_t *transMatrix = matrixW.GetTranslation();
 162     matrixScaled.RotateZ(-rotMatrix[1]*TMath::RadToDeg()*weight);
 163     matrixScaled.RotateY(rotMatrix[2]*TMath::RadToDeg()*weight);
 164     matrixScaled.RotateX(-rotMatrix[5]*TMath::RadToDeg()*weight);
 165     for (Int_t i=0; i<3; i++) delta[i]=weight*transMatrix[i];
 166     matrixScaled.SetTranslation(delta);
 167     // (matrixScaled*matrixW).Print(); in case weight -1 should be unit matrix
 168     //
 169     if (!fMatrixGlobal)  {
 170       fMatrixGlobal = new TGeoHMatrix(matrixScaled);
 171     }else{
 172       ((TGeoHMatrix*)fMatrixGlobal)->Multiply(&matrixScaled);
 173     }
 174   }
 175   if (add.fMatrixGlobalDelta){
 176     TGeoHMatrix matrixW=*(add.fMatrixGlobalDelta);
 177     TGeoHMatrix matrixScaled;
 178     const Double_t *rotMatrix   = matrixW.GetRotationMatrix();
 179     const Double_t *transMatrix = matrixW.GetTranslation();
 180     matrixScaled.RotateZ(-rotMatrix[1]*TMath::RadToDeg()*weight);
 181     matrixScaled.RotateY(rotMatrix[2]*TMath::RadToDeg()*weight);
 182     matrixScaled.RotateX(-rotMatrix[5]*TMath::RadToDeg()*weight);
 183     for (Int_t i=0; i<3; i++) delta[i]=weight*transMatrix[i];
 184     matrixScaled.SetTranslation(delta);
 185     // (matrixScaled*matrixW).Print(); in case weight -1 should be unit matrix
 186     //
 187     if (!fMatrixGlobalDelta)  {
 188       fMatrixGlobalDelta = new TGeoHMatrix(matrixScaled);
 189     }else{
 190       ((TGeoHMatrix*)fMatrixGlobalDelta)->Multiply(&matrixScaled);
 191     }
 192   }
 193   if (add.fArraySector){
 194     if (!fArraySector) {
 195       fArraySector = new TObjArray(72);
 196       for (Int_t isec=0; isec<72; isec++) fArraySector->AddAt(new TGeoHMatrix,isec);
 197     }
 198     for (Int_t isec=0; isec<72; isec++){
 199       TGeoHMatrix *mat0= (TGeoHMatrix*)fArraySector->At(isec);
 200       TGeoHMatrix *mat1= (TGeoHMatrix*)add.fArraySector->At(isec);
 201       if (mat0&&mat1){
 202         TGeoHMatrix matrixW=*(mat1);
 203         TGeoHMatrix matrixScaled;
 204         const Double_t *rotMatrix   = matrixW.GetRotationMatrix();
 205         const Double_t *transMatrix = matrixW.GetTranslation();
 206         matrixScaled.RotateZ(-rotMatrix[1]*TMath::RadToDeg()*weight);
 207         matrixScaled.RotateY(rotMatrix[2]*TMath::RadToDeg()*weight);
 208         matrixScaled.RotateX(-rotMatrix[5]*TMath::RadToDeg()*weight);
 209         for (Int_t i=0; i<3; i++) delta[i]=weight*transMatrix[i];
 210         matrixScaled.SetTranslation(delta);
 211           mat0->Multiply(&matrixScaled);
 212       }
 213     }
 214   }
 215   //
 216   return kTRUE;
 217 }
 218
 219
 220
 221
 222 void AliTPCCalibGlobalMisalignment::SetQuadranAlign(const TVectorD *quadrantQ0, const TVectorD *quadrantRQ0, const TVectorD *quadrantQ1,const TVectorD *quadrantRQ1,  const TVectorD *quadrantQ2,  const TVectorD *quadrantRQ2){
 223   //
 224   // Set quadrant alignment
 225   // 6 vectors for 36 (super) sectors
 226   //
 227   if (quadrantQ0) fQuadrantQ0   = new TVectorD(*quadrantQ0);
 228   if (quadrantRQ0) fQuadrantRQ0 = new TVectorD(*quadrantRQ0);
 229   //
 230   if (quadrantQ1) fQuadrantQ1   = new TVectorD(*quadrantQ1);
 231   if (quadrantQ1) fQuadrantRQ1  = new TVectorD(*quadrantRQ1);
 232   if (quadrantQ2) fQuadrantQ2   = new TVectorD(*quadrantQ2);
 233   if (quadrantQ2) fQuadrantRQ2  = new TVectorD(*quadrantRQ2);
 234 }
 235
 236
 237
 238 void AliTPCCalibGlobalMisalignment::SetAlignGlobal(const TGeoMatrix * matrixGlobal){
 239   //
 240   // Set global misalignment
 241   // Object is OWNER
 242   //
 243   if (fMatrixGlobal) delete fMatrixGlobal;
 244   fMatrixGlobal=0;
 245   if (matrixGlobal) fMatrixGlobal = new TGeoHMatrix(*matrixGlobal);
 246 }
 247
 248 void AliTPCCalibGlobalMisalignment::SetAlignGlobalDelta(const TGeoMatrix * matrixGlobalDelta){
 249   //
 250   // Set global misalignment
 251   // Object is OWNER
 252   //
 253   if (fMatrixGlobalDelta) delete fMatrixGlobalDelta;
 254   fMatrixGlobalDelta=0;
 255   if (matrixGlobalDelta) fMatrixGlobalDelta = new TGeoHMatrix(*matrixGlobalDelta);
 256 }
 257
 258 void AliTPCCalibGlobalMisalignment::SetAlignSectors(const TObjArray *arraySector){
 259   //
 260   // Set misalignment TObjArray of TGeoMatrices  - for each sector
 261   // Object is OWNER
 262   //
 263   if (fArraySector) delete fArraySector;
 264   fArraySector=0;
 265   if (arraySector) fArraySector = (TObjArray*)arraySector->Clone();
 266 }
 267
 268
 269 //void AliTPCCalibGlobalMisalignment::Init() {
 270 //  //
 271 // // Initialization funtion
 272 //  //
 273
 274 //  // nothing to be initialized, results of this calibration class will go to the global aligment structure
 275
 276 //}
 277
 278 //void AliTPCCalibGlobalMisalignment::Update(const TTimeStamp &/*timeStamp*/) {
 279 //  //
 280 //  // Update function
 281 //  //
 282 //
 283 //  // nothing to be updated, results of this calibration class will go to the global aligment structure
 284 //
 285 //}
 286
 287
 288
 289 void AliTPCCalibGlobalMisalignment::GetCorrection(const Float_t x[],const Short_t roc,Float_t dx[]) {
 290   //
 291   // Calculates the simple correction due to a shift (in x,y,z) or an rotation of the TPC (around z)
 292   //
 293   static AliTPCROC *tpcRoc =AliTPCROC::Instance();
 294   Double_t xref  = ( tpcRoc->GetPadRowRadii(0,0)+tpcRoc->GetPadRowRadii(36,tpcRoc->GetNRows(36)-1))*0.5;
 295   Double_t xquadrant  = tpcRoc->GetPadRowRadii(36,53); // row 53 from uli
 296   Double_t xIO  = ( tpcRoc->GetPadRowRadii(0,tpcRoc->GetNRows(0)-1)+tpcRoc->GetPadRowRadii(36,0))*0.5;
 297   Double_t r=0, phi=0;
 298   r   = TMath::Sqrt( x[0]*x[0] + x[1]*x[1] );
 299   phi = TMath::ATan2(x[1],x[0]);
 300   // Getsector number
 301   Double_t sec=TMath::Nint(-0.5+(phi*9./TMath::Pi()));
 302   if (sec<0) sec+=18;
 303   Int_t isec = TMath::Nint(sec);
 304   if  (roc%36>=18) isec+=18;
 305   //
 306   // Get the point on the local coordiante frame
 307   //
 308   Double_t alpha=(sec+0.5)*TMath::Pi()/9;
 309   Double_t pos[3]={0,0,x[2]};
 310   pos[0]=  TMath::Cos(alpha)*x[0]+TMath::Sin(alpha)*x[1];
 311   pos[1]= -TMath::Sin(alpha)*x[0]+TMath::Cos(alpha)*x[1];
 312   if (pos[0]>tpcRoc->GetPadRowRadiiUp(0))  isec+=36;
 313
 314   //
 315   // apply quadrant alignment if available - in local coordinate frame
 316   //
 317   //
 318   Double_t posQG[3]={x[0],x[1],x[2]};
 319   {
 320     Double_t dly=0;
 321     Bool_t isQ0 = (pos[0]<xquadrant)&&(pos[0]>xIO);
 322     Bool_t isQ1 = (pos[0]>xquadrant);
 323     Double_t  sign = (pos[1]>0)? 1.: -1.;
 324     if (isQ0){
 325       if (fQuadrantQ0)  dly+=sign*(*fQuadrantQ0)[isec%36];  // shift in cm
 326       if (fQuadrantRQ0) dly+=sign*(*fQuadrantRQ0)[isec%36]*(pos[0]-xref);
 327     }
 328     if (isQ1){
 329       if (fQuadrantQ1)  dly+=sign*(*fQuadrantQ1)[isec%36];  // shift in cm
 330       if (fQuadrantRQ1) dly+=sign*(*fQuadrantRQ1)[isec%36]*(pos[0]-xref);
 331       if (fQuadrantQ2)  dly+=(*fQuadrantQ2)[isec%36];  // shift in cm
 332       if (fQuadrantRQ2) dly+=(*fQuadrantRQ2)[isec%36]*(pos[0]-xref);
 333     }
 334     // Tranform the corrected point to the global frame
 335     posQG[0]=  TMath::Cos(alpha)*pos[0]-TMath::Sin(alpha)*(pos[1]+dly);
 336     posQG[1]=  TMath::Sin(alpha)*pos[0]+TMath::Cos(alpha)*(pos[1]+dly);
 337   }
 338   //
 339   // rotation of the read-out planes
 340   if  (roc%36<18) // A side
 341     phi += fRotPhiA;
 342   else         // C side
 343     phi += fRotPhiC;
 344
 345   // Simply adding a constant dRPHi residual. PURELY FOR CALIBRATION PURPOSES
 346   if  (roc%36<18) // A side
 347     phi += fdRPhiOffsetA/r;
 348   else         // C side
 349     phi += fdRPhiOffsetC/r;
 350
 351   dx[0] = r * TMath::Cos(phi) - x[0];
 352   dx[1] = r * TMath::Sin(phi) - x[1];
 353   dx[2] = 0.;
 354
 355   // Simple shifts
 356   dx[0] -= fXShift;
 357   dx[1] -= fYShift;
 358   dx[2] -= fZShift;
 359   // quadrant shifts
 360   dx[0] += (posQG[0]-x[0]);
 361   dx[1] += (posQG[1]-x[1]);
 362   //
 363   //
 364   if (fMatrixGlobal){
 365     // apply global alignment matrix
 366     Double_t ppos[3]={x[0],x[1],x[2]};
 367     Double_t pposC[3]={x[0],x[1],x[2]};
 368     fMatrixGlobal->LocalToMaster(ppos,pposC);
 369     dx[0]+=pposC[0]-ppos[0];
 370     dx[1]+=pposC[1]-ppos[1];
 371     dx[2]+=pposC[2]-ppos[2];
 372   }
 373   if (fMatrixGlobalDelta){
 374     // apply global alignment matrix A-C Side side
 375     Double_t ppos[3]={x[0],x[1],x[2]};
 376     Double_t pposC[3]={x[0],x[1],x[2]};
 377     fMatrixGlobalDelta->LocalToMaster(ppos,pposC);
 378     Double_t ssign=(roc%36<18) ? 1.:-1.;
 379     dx[0]+=ssign*(pposC[0]-ppos[0]);
 380     dx[1]+=ssign*(pposC[1]-ppos[1]);
 381     dx[2]+=ssign*(pposC[2]-ppos[2]);
 382   }
 383
 384   if (fArraySector){
 385     // apply global alignment matrix
 386     TGeoMatrix  *mat = (TGeoMatrix*)fArraySector->At(isec);
 387     if (mat){
 388       Double_t ppos[3]={x[0],x[1],x[2]};
 389       Double_t pposC[3]={x[0],x[1],x[2]};
 390       mat->LocalToMaster(ppos,pposC);
 391       dx[0]+=pposC[0]-ppos[0];
 392       dx[1]+=pposC[1]-ppos[1];
 393       dx[2]+=pposC[2]-ppos[2];
 394     }
 395   }
 396 }
 397
 398 void AliTPCCalibGlobalMisalignment::Print(Option_t* option )  const{
 399   //
 400   // Print function to check the settings
 401   //
 402   printf("%s",GetTitle());
 403   printf(" - Trivial Misalignments for calibration purposes: \n");
 404   printf(" - X-Shift: %1.3f cm, Y-Shift: %1.3f cm, Z-Shift: %1.3f cm \n",fXShift,fYShift,fZShift);
 405   printf(" - Phi-Rotations: A side: %1.5f rad, C side: %1.5f rad\n",fRotPhiA,fRotPhiC);
 406   printf(" - dRPhi offsets: A side: %1.5f cm, C side: %1.5f cm\n",fdRPhiOffsetA,fdRPhiOffsetC);
 407   TString opt = option; opt.ToLower();
 408   if (opt.Contains("a")){
 409     if (GetAlignGlobal()){
 410       printf("GetAlignGlobal()\n");
 411       GetAlignGlobal()->Print();
 412     }
 413     if (GetAlignGlobalDelta()){
 414       printf("GetAlignGlobalDelta()\n");
 415       GetAlignGlobalDelta()->Print();
 416     }
 417     if (GetAlignSectors()){
 418       printf("GetAlignSectors()\n");
 419       GetAlignSectors()->Print();
 420     }
 421   }
 422 }
 423
 424 void AliTPCCalibGlobalMisalignment::AddAlign(const  AliTPCCalibGlobalMisalignment & add){
 425   //
 426   // Add the alignmnet to current object
 427   //
 428   fXShift+=add.fXShift;               // Shift in global X [cm]
 429   fYShift+=add.fYShift;               // Shift in global Y [cm]
 430   fZShift+=add.fZShift;               // Shift in global Z [cm]
 431
 432   fRotPhiA+=add.fRotPhiA;      // simple rotation of A side read-out plane around the Z axis [rad]
 433   fRotPhiC+=add.fRotPhiC;      // simple rotation of C side read-out plane around the Z axis [rad]
 434   fdRPhiOffsetA+=add.fdRPhiOffsetA;  // add a constant offset of dRPhi (or local Y) in [cm]: purely for calibration purposes!
 435   fdRPhiOffsetC+=add.fdRPhiOffsetC;  // add a constant offset of dRPhi (or local Y) in [cm]: purely for calibration purposes!
 436   //
 437   // Quadrant alignment
 438   //
 439   if (add.fQuadrantQ0) {
 440     if (fQuadrantQ0) fQuadrantQ0->Add(*(add.fQuadrantQ0));
 441     if (!fQuadrantQ0) fQuadrantQ0 = (TVectorD*)(add.fQuadrantQ0->Clone());
 442   }
 443   if (add.fQuadrantRQ0) {
 444     if (fQuadrantRQ0) fQuadrantRQ0->Add(*(add.fQuadrantRQ0));
 445     if (!fQuadrantRQ0) fQuadrantRQ0 = (TVectorD*)(add.fQuadrantRQ0->Clone());
 446   }
 447   //
 448   if (add.fQuadrantQ1) {
 449     if (fQuadrantQ1) fQuadrantQ1->Add(*(add.fQuadrantQ1));
 450     if (!fQuadrantQ1) fQuadrantQ1 = (TVectorD*)(add.fQuadrantQ1->Clone());
 451   }
 452   if (add.fQuadrantRQ1) {
 453     if (fQuadrantRQ1) fQuadrantRQ1->Add(*(add.fQuadrantRQ1));
 454     if (!fQuadrantRQ1) fQuadrantRQ1 = (TVectorD*)(add.fQuadrantRQ1->Clone());
 455   }
 456   //
 457   if (add.fQuadrantQ2) {
 458     if (fQuadrantQ2) fQuadrantQ2->Add(*(add.fQuadrantQ2));
 459     if (!fQuadrantQ2) fQuadrantQ2 = (TVectorD*)(add.fQuadrantQ2->Clone());
 460   }
 461   if (add.fQuadrantRQ2) {
 462     if (fQuadrantRQ2) fQuadrantRQ2->Add(*(add.fQuadrantRQ2));
 463     if (!fQuadrantRQ2) fQuadrantRQ2 = (TVectorD*)(add.fQuadrantRQ2->Clone());
 464   }
 465   //
 466   // Global alignment - use native ROOT representation
 467   //
 468   if (add.fMatrixGlobal){
 469     if (!fMatrixGlobal)  fMatrixGlobal = new TGeoHMatrix(*(add.fMatrixGlobal));
 470     if (fMatrixGlobal)   ((TGeoHMatrix*)fMatrixGlobal)->Multiply(add.fMatrixGlobal);
 471   }
 472   if (add.fArraySector){
 473     if (!fArraySector) {SetAlignSectors(add.fArraySector);
 474     }else{
 475       for (Int_t isec=0; isec<72; isec++){
 476         TGeoHMatrix *mat0= (TGeoHMatrix*)fArraySector->At(isec);
 477         TGeoHMatrix *mat1= (TGeoHMatrix*)add.fArraySector->At(isec);
 478         if (mat0&&mat1) mat0->Multiply(mat1);
 479       }
 480     }
 481   }
 482 }
 483
 484
 485 AliTPCCalibGlobalMisalignment *  AliTPCCalibGlobalMisalignment::CreateOCDBAlign(){
 486   //
 487   // Create  AliTPCCalibGlobalMisalignment from OCDB Alignment entry
 488   // OCDB has to be initialized before in user code
 489   // All storages (defualt and specific)  and run number
 490   //
 491   AliCDBEntry * entry = AliCDBManager::Instance()->Get("TPC/Align/Data");
 492   if (!entry){
 493     printf("Missing alignmnet entry. OCDB not initialized?\n");
 494     return 0;
 495   }
 496   TClonesArray * array = (TClonesArray*)entry->GetObject();
 497   Int_t entries = array->GetEntries();
 498   TGeoHMatrix matrixGlobal;
 499   TObjArray *alignArrayOCDB= new TObjArray(73);  // sector misalignment + global misalignment
 500   //                                            // global is number 72
 501   //
 502   { for (Int_t i=0;i<entries; i++){
 503       //
 504       //
 505       TGeoHMatrix matrix;
 506       AliAlignObjParams *alignP = (AliAlignObjParams*)array->UncheckedAt(i);
 507       alignP->GetMatrix(matrix);
 508       Int_t imod;
 509       AliGeomManager::ELayerID ilayer;
 510       alignP->GetVolUID(ilayer, imod);
 511       if (ilayer==AliGeomManager::kInvalidLayer) {
 512         alignArrayOCDB->AddAt(matrix.Clone(),72);
 513         alignP->GetMatrix(matrixGlobal);
 514       }else{
 515         Int_t sector=imod;
 516         if (ilayer==AliGeomManager::kTPC2) sector+=36;
 517         alignArrayOCDB->AddAt(matrix.Clone(),sector);
 518       }
 519     }
 520   }
 521   AliTPCCalibGlobalMisalignment *align = new  AliTPCCalibGlobalMisalignment;
 522   align->SetAlignGlobal(&matrixGlobal);
 523   align->SetAlignSectors(alignArrayOCDB);
 524   return align;
 525 }
 526
 527
 528 AliTPCCalibGlobalMisalignment *  AliTPCCalibGlobalMisalignment::CreateMeanAlign(const AliTPCCalibGlobalMisalignment *alignIn){
 529   //
 530   // Create new object, disantangle common mean alignmnet and sector alignment
 531   //
 532   // 1. Try to get mean alignment
 533   // 2. Remove mean alignment from sector alignment
 534   // 3. Create new object
 535   //
 536   TObjArray * array = alignIn->GetAlignSectors();
 537   TObjArray * arrayNew = new TObjArray(72);
 538   //
 539   //Get mean transformation
 540   TGeoHMatrix matrix;
 541   {for (Int_t isec=0; isec<72; isec++){
 542       const TGeoMatrix* cmatrix=(TGeoMatrix*)array->At(isec);
 543       if (!cmatrix) continue;
 544       matrix.Multiply(cmatrix);
 545     }}
 546   TGeoHMatrix matrixMean(matrix);
 547   matrixMean.SetDx(matrix.GetTranslation()[0]/72.);
 548   matrixMean.SetDy(matrix.GetTranslation()[1]/72.);
 549   matrixMean.SetDz(matrix.GetTranslation()[2]/72.);
 550   Double_t rotation[12];
 551   {for (Int_t i=0; i<12; i++) {
 552       rotation[i]=1.0;
 553       if (TMath::Abs(matrix.GetRotationMatrix()[i]-1.)>0.1){
 554       rotation[i]=matrix.GetRotationMatrix()[i]/72.;
 555       }
 556     }}
 557   matrixMean.SetRotation(rotation);
 558   TGeoHMatrix matrixInv = matrixMean.Inverse();
 559   //
 560   {for (Int_t isec=0; isec<72; isec++){
 561       TGeoHMatrix* amatrix=(TGeoHMatrix*)(array->At(isec)->Clone());
 562       if (!amatrix) continue;
 563       amatrix->Multiply(&matrixInv);
 564       arrayNew->AddAt(amatrix,isec);
 565     }}
 566   if (alignIn->GetAlignGlobal()) matrixMean.Multiply((alignIn->GetAlignGlobal()));
 567   AliTPCCalibGlobalMisalignment *alignOut = new  AliTPCCalibGlobalMisalignment;
 568   alignOut->SetAlignGlobal(&matrixMean);
 569   alignOut->SetAlignSectors(arrayNew);
 570   /*
 571     Checks transformation:
 572     AliTPCCalibGlobalMisalignment *  alignIn =  AliTPCCalibGlobalMisalignment::CreateOCDBAlign()
 573     AliTPCCalibGlobalMisalignment * alignOut =  AliTPCCalibGlobalMisalignment::CreateMeanAlign(alignIn)
 574     alignOutM= (AliTPCCalibGlobalMisalignment*)alignOut->Clone();
 575     alignOutS= (AliTPCCalibGlobalMisalignment*)alignOut->Clone();
 576     alignOutS->SetAlignGlobal(0);
 577     alignOutM->SetAlignSectors(0);
 578     //
 579     AliTPCCorrection::AddVisualCorrection(alignOut,0);
 580     AliTPCCorrection::AddVisualCorrection(alignOutM,1);
 581     AliTPCCorrection::AddVisualCorrection(alignOutS,2);
 582     AliTPCCorrection::AddVisualCorrection(alignIn,3);
 583
 584     TF1 f0("f0","AliTPCCorrection::GetCorrSector(x,85,0.9,1,0)",0,18);
 585     TF1 f1("f1","AliTPCCorrection::GetCorrSector(x,85,0.9,1,1)",0,18);
 586     TF1 f2("f2","AliTPCCorrection::GetCorrSector(x,85,0.9,1,2)",0,18);
 587     TF1 f3("f3","AliTPCCorrection::GetCorrSector(x,85,0.9,1,3)",0,18);
 588     f0->SetLineColor(1);
 589     f1->SetLineColor(2);
 590     f2->SetLineColor(3);
 591     f3->SetLineColor(4);
 592     f0->Draw();
 593     f1->Draw("same");
 594     f2->Draw("same");
 595     f3->Draw("same");
 596
 597     TF2 f2D("f2D","AliTPCCorrection::GetCorrSector(x,y,0.9,1,0)-AliTPCCorrection::GetCorrSector(x,y,0.9,1,3)",0,18,85,245);
 598   */
 599   return alignOut;
 600 }
 601
 602
 603 void AliTPCCalibGlobalMisalignment::DumpAlignment( AliTPCCalibGlobalMisalignment* align, TTreeSRedirector *pcstream, const char *name){
 604   //
 605   // Dump alignment per sector into tree
 606   //
 607   TObjArray * array = align->GetAlignSectors();
 608   if (!array) return;
 609   //
 610   //Get mean transformation
 611   TGeoHMatrix matrix;
 612   {for (Int_t isec=0; isec<72; isec++){
 613       TGeoHMatrix* cmatrix=(TGeoHMatrix*)array->At(isec);
 614       TGeoHMatrix* cmatrixDown=(TGeoHMatrix*)array->At(isec%36);
 615       TGeoHMatrix* cmatrixUp=(TGeoHMatrix*)array->At(isec%36+36);
 616       TGeoHMatrix diff(*cmatrixDown);
 617       diff.Multiply(&(cmatrixUp->Inverse()));
 618       (*pcstream)<<name<<
 619         "isec="<<isec<<
 620         "m0.="<<cmatrix<<
 621         "diff.="<<&diff<<
 622         "\n";
 623     }
 624   }
 625
 626 }