TPC/AliTPCTransformation.cxx

   1 /*
   2
   3  Class AliTPCtransformation:
   4  Should represent general non linear transformation. Currently tune for TPConly.
   5  To be used:
   6  1. Simulation-Digitization
   7  2. Reconstruction - AliTPCTransform
   8  3. Calibration/Alignment (KalmanFilter, Milipedde)
   9  4. Set of transformation to be stored/retrieved as OCDB entry
  10
  11  Base functionality:
  12
  13  1. Double_t GetDeltaXYZ(Int_t coord, Int_t volID, Double_t param, Double_t x, Double_t y, Double_t z)
  14  Get correction - return the delta of coordinate coord dx or dy or dz for given volID for point at point (x,y,z)
  15  All coordinates are global
  16
  17  2. The transformation should work only for given volIDs and detector IDs
  18     Currently Bitmask is used for filtering
  19
  20
  21  Transformation - naming convention:
  22  //
  23  XXX(local)YYYZZZ
  24  TPClocaldLxdGX
  25  XXX   - detector if detector specific
  26  local - if local transforamtion
  27  YYY   - type of transformation
  28  ZZZ   - return type of transformation
  29
  30 */
  31
  32 #include <string.h>
  33 #include "TRandom.h"
  34 #include "TMath.h"
  35 #include "TBits.h"
  36 #include "TFormula.h"
  37 #include "TF1.h"
  38 #include "TLinearFitter.h"
  39 #include "TFile.h"
  40 #include "TObjString.h"
  41
  42 #include "TTreeStream.h"
  43 #include "AliTrackPointArray.h"
  44 #include "AliLog.h"
  45 #include "AliTPCTransformation.h"
  46
  47 ClassImp(AliTPCTransformation)
  48
  49
  50 AliTPCTransformation::GenFuncG    AliTPCTransformation::fgFormulas[10000];
  51 TObjArray*  AliTPCTransformation::fgFormulasName = new TObjArray(10000);
  52
  53
  54 void AliTPCTransformation::RegisterFormula(const char * name, GenFuncG formula){
  55   //
  56   // Add Formula to the list of formulas
  57   //
  58   Int_t last= fgFormulasName->GetEntries();
  59   fgFormulasName->AddLast(new TObjString(name));
  60   fgFormulas[last]=formula;
  61 }
  62
  63 Int_t  AliTPCTransformation::BuildBasicFormulas(){
  64   //
  65
  66   //
  67   //build list of basic TPC formulas - corrections
  68   //
  69   RegisterFormula("TPCscalingRPol",(GenFuncG)(AliTPCTransformation::TPCscalingRPol));
  70   RegisterFormula("TPCscalingRIFC",(GenFuncG)(AliTPCTransformation::TPCscalingRIFC));
  71   RegisterFormula("TPCscalingROFC",(GenFuncG)(AliTPCTransformation::TPCscalingROFC));
  72   //
  73   RegisterFormula("TPCdeltaFCROC",(GenFuncG)(AliTPCTransformation::TPCdeltaFCROC));
  74   RegisterFormula("TPCdeltaFCCE",(GenFuncG)(AliTPCTransformation::TPCdeltaFCCE));
  75   //
  76   RegisterFormula("TPCscalingZDr",(GenFuncG)(AliTPCTransformation::TPCscalingZDrift));
  77   RegisterFormula("TPCscalingZDrGy",(GenFuncG)(AliTPCTransformation::TPCscalingZDriftGy));
  78   RegisterFormula("TPCscalingZDriftT0",(GenFuncG)(AliTPCTransformation::TPCscalingZDriftT0));
  79   RegisterFormula("TPCscalingPhiLocal",(GenFuncG)(AliTPCTransformation::TPCscalingPhiLocal));
  80   RegisterFormula("TPClocalRPhiEdge",(GenFuncG)(AliTPCTransformation::TPClocalRPhiEdge));
  81   //
  82   // TPC Local X and Y misalignment + rotation
  83   //
  84   RegisterFormula("TPClocaldLxdGX",(GenFuncG)(AliTPCTransformation::TPClocaldLxdGX));
  85   RegisterFormula("TPClocaldLxdGY",(GenFuncG)(AliTPCTransformation::TPClocaldLxdGY));
  86   RegisterFormula("TPClocaldLydGX",(GenFuncG)(AliTPCTransformation::TPClocaldLydGX));
  87   RegisterFormula("TPClocaldLydGY",(GenFuncG)(AliTPCTransformation::TPClocaldLydGY));
  88   RegisterFormula("TPClocaldRzdGX",(GenFuncG)(AliTPCTransformation::TPClocaldRzdGX));
  89   RegisterFormula("TPClocaldRzdGY",(GenFuncG)(AliTPCTransformation::TPClocaldRzdGY));
  90
  91   //
  92   // Z offset
  93   //
  94   RegisterFormula("TPCDeltaZ",(GenFuncG)(AliTPCTransformation::TPCDeltaZ));
  95   RegisterFormula("TPCDeltaZMediumLong",(GenFuncG)(AliTPCTransformation::TPCDeltaZMediumLong));
  96   RegisterFormula("TPCTiltingZ",(GenFuncG)(AliTPCTransformation::TPCTiltingZ));
  97   return 0;
  98 }
  99
 100 AliTPCTransformation::GenFuncG  AliTPCTransformation::FindFormula(const char * name){
 101   //
 102   // find formula - if registered
 103   //
 104   if (fgFormulasName->FindObject(name)==0) return 0;
 105   Int_t entries = fgFormulasName->GetEntries();
 106   for (Int_t i=0;i<entries;i++){
 107     if (strcmp(fgFormulasName->At(i)->GetName(), name)==0){
 108       return fgFormulas[i];
 109     }
 110   }
 111   return 0;
 112 }
 113
 114 Double_t AliTPCTransformation::Eval(const char * name, const Double_t*x,const Double_t*par){
 115   //
 116   // Only for test purposes - very slow
 117   //
 118   GenFuncG fun = FindFormula(name);
 119   if (!fun) return 0;
 120   return fun(x,par);
 121 }
 122
 123
 124 AliTPCTransformation::AliTPCTransformation():
 125   TNamed(),
 126   fNameX(0),          // x formula name
 127   fNameY(0),          // y formula name
 128   fNameZ(0),          // z formula name
 129   //
 130   fBitMask(0),        // bitmaps - transformation only for specified volID
 131   fCoordSystem(0),    // coord system of  output deltas
 132   fParam(0),          // free parameter of transformation
 133   fSigma(0),          // uncertainty of the parameter
 134   fSigmaMax(0),       //maximal sigma (Not allowed to increase in propagate time by bigger factor)
 135   fSigma2Time(0),     // change of the error in time - (For kalman filter)
 136   fFixedParam(0),     // fixed parameters of tranformation
 137   fIsActive(kTRUE),   // swith - On/Off
 138   //
 139   fInit(kFALSE),      // initialization flag - set to kTRUE if corresponding formulas found
 140   fFormulaX(0),       // x formula - pointer to the function
 141   fFormulaY(0),       // y formula - pointer to the function
 142   fFormulaZ(0)       // z formula - pointer to the function
 143   //
 144 {
 145   //
 146   // default constructor
 147   //
 148 }
 149
 150
 151
 152 AliTPCTransformation::AliTPCTransformation(const char *name, TBits *mask, const char *fx, const char *fy, const char *fz, Int_t coordSystem):
 153   TNamed(name,name),
 154   fNameX(0),       // x formula name
 155   fNameY(0),       // y formula name
 156   fNameZ(0),       // z formula name
 157   fBitMask(mask),   // bitmaps - transformation only for specified volID
 158   fCoordSystem(coordSystem), // coordinate system of output deltas
 159   fParam(0),          // free parameter of transformation
 160   fSigma(0),
 161   fSigmaMax(0),       //maximal sigma (Not allowed to increase in propagate time by bigger factor)
 162   fSigma2Time(0),     // change of sigma in time
 163   fFixedParam(0),     // fixed parameters of tranformation
 164   fIsActive(kTRUE),   // swith - On/Off
 165   //
 166   fInit(kFALSE),      // initialization flag - set to kTRUE if corresponding formulas found
 167   fFormulaX(0),       // x formula - pointer to the function
 168   fFormulaY(0),       // y formula - pointer to the function
 169   fFormulaZ(0)       // z formula - pointer to the function
 170 {
 171   //
 172   // non default constructor
 173   //
 174   if (fx) fNameX= new TString(fx);
 175   if (fy) fNameY= new TString(fy);
 176   if (fz) fNameZ= new TString(fz);
 177   Init();
 178 }
 179
 180 AliTPCTransformation::AliTPCTransformation(const AliTPCTransformation&trafo):
 181   TNamed(trafo),
 182   fNameX(0),          // x formula name
 183   fNameY(0),          // y formula name
 184   fNameZ(0),          // z formula name
 185                      //
 186   fBitMask(0),        // bitmaps - transformation only for specified volID
 187   fCoordSystem(0),    // coord system of  output deltas
 188   fParam(trafo.fParam),          // free parameter of transformation
 189   fSigma(trafo.fSigma),          // uncertainty of the parameter
 190   fSigmaMax(trafo.fSigma),       //maximal sigma (Not allowed to increase in propagate time by bigger factor)
 191   fSigma2Time(trafo.fSigma2Time),     // change of the error in time - (For kalman filter)
 192   fFixedParam(0),     // fixed parameters of tranformation
 193   fIsActive(trafo.fIsActive),   // swith - On/Off
 194   //
 195   fInit(kFALSE),      // initialization flag - set to kTRUE if corresponding formulas found
 196   fFormulaX(0),       // x formula - pointer to the function
 197   fFormulaY(0),       // y formula - pointer to the function
 198   fFormulaZ(0)       // z formula - pointer to the function
 199 {
 200   if (trafo.fNameX) fNameX = new TString(*(trafo.fNameX));
 201   if (trafo.fNameY) fNameY = new TString(*(trafo.fNameY));
 202   if (trafo.fNameZ) fNameZ = new TString(*(trafo.fNameZ));
 203   if (trafo.fBitMask)  fBitMask = new TBits(*(trafo.fBitMask));
 204 }
 205
 206 AliTPCTransformation::~AliTPCTransformation(){
 207   //
 208   // destructor
 209   //
 210   delete fNameX;
 211   delete fNameY;
 212   delete fNameZ;
 213   delete fBitMask;
 214   delete fFixedParam;
 215 }
 216 void AliTPCTransformation::SetParams(Double_t param, Double_t sigma, Double_t sigma2Time, const TVectorD *const fixedParams){
 217   //
 218   // Set parameters of transformation
 219   //
 220   fParam = param;
 221   fSigma = sigma;
 222   fSigmaMax = sigma;
 223   fSigma2Time = sigma2Time;
 224   if (fFixedParam) delete fFixedParam;
 225   fFixedParam = new TVectorD(*fixedParams);
 226   Init();
 227 }
 228
 229
 230 Bool_t AliTPCTransformation::Init(){
 231   //
 232   // associate formulas with pointer to the function
 233   //
 234   Bool_t isOK=kTRUE;
 235   if (fNameX) {
 236     fFormulaX=FindFormula(fNameX->Data());
 237     if (fFormulaX==0) isOK=kFALSE;
 238   }
 239   if (fNameY) {
 240     fFormulaY=FindFormula(fNameY->Data());
 241     if (fFormulaY==0) isOK=kFALSE;
 242   }
 243   if (fNameZ) {
 244     fFormulaZ=FindFormula(fNameZ->Data());
 245     if (!fFormulaZ) isOK=kFALSE;
 246   }
 247   return isOK;
 248 }
 249
 250
 251
 252 TBits * AliTPCTransformation::BitsSide(Bool_t aside){
 253   //
 254   TBits * bits = new TBits(72);
 255   for (Int_t i=0; i<72;i++){
 256     if (i%36<18 && aside) (*bits)[i]=kTRUE;
 257     if (i%36<18 && (!aside)) (*bits)[i]=kFALSE;
 258     if (i%36>=18 && aside) (*bits)[i]=kFALSE;
 259     if (i%36>=18 && (!aside)) (*bits)[i]=kTRUE;
 260   }
 261   return bits;
 262 }
 263
 264 TBits * AliTPCTransformation::BitsAll(){
 265   //
 266   //
 267   //
 268   TBits * bits = new TBits(72);
 269   for (Int_t i=0; i<72;i++){
 270     (*bits)[i]=kTRUE;
 271   }
 272   return bits;
 273 }
 274
 275 // Double_t AliTPCTransformation::GetDeltaXYZ(Int_t coord, Int_t volID, Double_t param, Double_t x, Double_t y, Double_t z){
 276 //   //
 277 //   //
 278 //   // coord - type of coordinate
 279 //   //       - 0 -X
 280 //   //         1 -Y
 281 //   //         2 -Z
 282 //   //         3 -R
 283 //   //         4 -RPhi
 284 //   if (!fIsActive) return 0;
 285 //   if (fBitMask && (!(*fBitMask)[volID])) return 0;
 286 //   Double_t xyz[5]={x,y,z, param,volID};
 287 //   if (fCoordSystem==0){
 288 //     // cartezian system
 289 //     if (coord==0 && fFormulaX) return fFormulaX(xyz,fFixedParam->GetMatrixArray());
 290 //     if (coord==1 && fFormulaY) return fFormulaY(xyz,fFixedParam->GetMatrixArray());
 291 //     if (coord==2 && fFormulaZ) return fFormulaZ(xyz,fFixedParam->GetMatrixArray());
 292 //   }
 293 //   if (fCoordSystem==1){
 294 //     // cylindrical system
 295 //     if (coord==2) {
 296 //       if (fFormulaZ==0) return 0;
 297 //       return fFormulaZ(xyz,fFixedParam->GetMatrixArray());
 298 //     }
 299 //     Double_t rrphiz[3]={0,0,0};
 300 //     if (fFormulaX) rrphiz[0] = fFormulaX(xyz,fFixedParam->GetMatrixArray());
 301 //     if (fFormulaY) rrphiz[1] = fFormulaY(xyz,fFixedParam->GetMatrixArray());
 302 //     Double_t alpha = TMath::ATan2(y,x);
 303 //     Double_t ca    = TMath::Cos(alpha);
 304 //     Double_t sa    = TMath::Sin(alpha);
 305 //     if (coord==0) return ca*rrphiz[0]-sa*rrphiz[1];
 306 //     if (coord==1) return sa*rrphiz[0]+ca*rrphiz[1];
 307 //   }
 308 //   return 0;
 309 // }
 310
 311
 312 Double_t AliTPCTransformation::GetDeltaXYZ(Int_t coord, Int_t volID, Double_t param, Double_t x, Double_t y, Double_t z){
 313   //
 314   //
 315   // coord - type of coordinate
 316   //       - 0 -X
 317   //         1 -Y
 318   //         2 -Z
 319   //         3 -R
 320   //         4 -RPhi
 321   //         5 -Z
 322   if (!fIsActive) return 0;
 323   if (fBitMask && (!(*fBitMask)[volID])) return 0;
 324   Double_t xyz[5]={x,y,z, param,volID};
 325   Double_t alpha = TMath::ATan2(y,x);
 326   Double_t ca    = TMath::Cos(alpha);
 327   Double_t sa    = TMath::Sin(alpha);
 328
 329   if (fCoordSystem==0){
 330     // cartezian system
 331     Double_t dxdydz[3]={0,0,0};
 332     if(fFormulaX) dxdydz[0]=fFormulaX(xyz,fFixedParam->GetMatrixArray());
 333     if(fFormulaY) dxdydz[1]=fFormulaY(xyz,fFixedParam->GetMatrixArray());
 334     if(fFormulaZ) dxdydz[2]=fFormulaZ(xyz,fFixedParam->GetMatrixArray());
 335
 336     if (coord==0) return dxdydz[0];
 337     if (coord==1) return dxdydz[1];
 338     if (coord==2) return dxdydz[2];
 339     if (coord==3) return dxdydz[0]*ca+dxdydz[1]*sa;
 340     if (coord==4) return -dxdydz[0]*sa+dxdydz[1]*ca;
 341     if (coord==5) return dxdydz[2];
 342   }
 343   if (fCoordSystem==1){
 344     // cylindrical system
 345     if (coord==2||coord==5) {
 346       if (fFormulaZ==0) return 0;
 347       return fFormulaZ(xyz,fFixedParam->GetMatrixArray());
 348     }
 349     Double_t rrphiz[3]={0,0,0};
 350     if (fFormulaX) rrphiz[0] = fFormulaX(xyz,fFixedParam->GetMatrixArray());
 351     if (fFormulaY) rrphiz[1] = fFormulaY(xyz,fFixedParam->GetMatrixArray());
 352     alpha = TMath::ATan2(y,x);
 353     ca    = TMath::Cos(alpha);
 354     sa    = TMath::Sin(alpha);
 355     if (coord==0) return ca*rrphiz[0]-sa*rrphiz[1];
 356     if (coord==1) return sa*rrphiz[0]+ca*rrphiz[1];
 357     if (coord==3) return rrphiz[0];
 358     if (coord==4) return rrphiz[1];
 359   }
 360   return 0;
 361 }
 362
 363
 364
 365 Double_t  AliTPCTransformation::TPCscalingRPol(Double_t *xyz, Double_t * param){
 366   //
 367   // Scaling and shift of TPC radius
 368   // xyz[0..2] - global xyz of point
 369   // xyz[3]    - scale parameter
 370   // param[0]  - radial scaling power
 371   // param[1]  - drift  scaling power
 372   // radius  from -1(at rInner)   to 1 (rOuter)
 373   // driftM  from -1(at 0 drift)  to 1 (250 cm drift)
 374
 375   Double_t rInner=78.8;
 376   Double_t rOuter=258.0;
 377   Double_t deltaR  = rOuter-rInner;
 378   Double_t radius  = (TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1])-rInner)*2./deltaR;
 379   Double_t driftM  = (0.5 - TMath::Abs(xyz[2]/250.))*2.0;
 380   Double_t delta   = TMath::Power(radius,param[0])*TMath::Power(driftM,param[1]);
 381   return delta*xyz[3];
 382 }
 383
 384
 385 Double_t  AliTPCTransformation::TPCscalingZDrift(Double_t *xyz, Double_t * param){
 386   //
 387   //
 388   // Scaling and shift of TPC radius
 389   // xyz[0..2] - global xyz of point
 390   // xyz[3]    - scale parameter
 391   Double_t driftP  = TMath::Power(1. - TMath::Abs(xyz[2]/250.), param[0]);
 392   Int_t    sector = TMath::Nint(xyz[4]);
 393   Double_t deltaZ  = (sector%36<18) ? -driftP : driftP;
 394   return deltaZ*xyz[3];
 395 }
 396
 397 Double_t  AliTPCTransformation::TPCscalingZDriftT0(Double_t *xyz, Double_t * /*param*/){
 398   //
 399   //
 400   // Z shift because time 0 offset
 401   // opposite on A and C side
 402   //
 403   // xyz[0..2] - global xyz of point
 404   // xyz[3]    - scale parameter
 405   Int_t    sector = TMath::Nint(xyz[4]);
 406   Double_t sign  = (sector%36<18) ? -1 : 1;
 407   return sign*xyz[3];
 408 }
 409
 410
 411 Double_t  AliTPCTransformation::TPCscalingZDriftGy(Double_t *xyz, Double_t * param){
 412   //
 413   //
 414   // Scaling and shift of TPC radius
 415   // xyz[0..2] - global xyz of point
 416   // xyz[3]    - scale parameter
 417   Double_t driftP  = TMath::Power(1. - TMath::Abs(xyz[2]/250.), param[0]);
 418   Double_t gy      = xyz[1]/250.;
 419   Int_t    sector = TMath::Nint(xyz[4]);
 420   Double_t deltaZ  = (sector%36<18) ? -driftP : driftP;
 421   return deltaZ*xyz[3]*gy;
 422 }
 423
 424
 425
 426 Double_t  AliTPCTransformation::TPCscalingPhiLocal(Double_t *xyz, Double_t * param){
 427   //
 428   //
 429   // Scaling if the local y -phi
 430   // xyz[0..2] - global xyz of point
 431   // xyz[3]    - scale parameter
 432   // value = 1 for ful drift length and parameter 1
 433   Double_t alpha       = TMath::ATan2(xyz[1],xyz[0]);
 434   Double_t sector      = TMath::Nint(9*alpha/TMath::Pi()-0.5);
 435   Double_t localAlpha  = (alpha-(sector+0.5)*TMath::Pi()/9.);
 436   Double_t radius      = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1])/250.;
 437   //
 438   Double_t deltaAlpha  = radius*TMath::Power(2.*9.*localAlpha/TMath::Pi(),param[0]);
 439   return deltaAlpha*xyz[3];
 440 }
 441
 442 Double_t  AliTPCTransformation::TPClocalRPhiEdge(Double_t *xyz, const Double_t *const param){
 443   //
 444   //
 445   // Scaling if the local y -phi
 446   // xyz[0..2] - global xyz of point
 447   // xyz[3]    - scale parameter
 448   // param[0]  - dedge offset - should be around gap size/2.
 449   // param[1]  - dedge factor - should be around gap size/2.
 450   Double_t alpha       = TMath::ATan2(xyz[1],xyz[0]);
 451   Double_t sector      = TMath::Nint(9*alpha/TMath::Pi()-0.5);
 452   Double_t localAlpha  = (alpha-(sector+0.5)*TMath::Pi()/9.);
 453   Double_t radius      = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
 454   Double_t deltaAlpha  = TMath::Pi()/18.-TMath::Abs(localAlpha);
 455   Double_t distEdge    = (deltaAlpha*radius);
 456   Double_t factor      = 1./(1.+(distEdge-param[0])/param[1]);
 457   return factor*xyz[3]*((localAlpha>0)? -1.:1.);
 458 }
 459
 460
 461 Double_t       AliTPCTransformation::TPCscalingRIFC(Double_t *xyz, Double_t * param){
 462   //
 463   // inner field cage r distorion - proportinal to 1 over distance to the IFC
 464   // param[0] - drift polynom order
 465   // distortion at first pad row - is normalized to
 466   Double_t rInner=78.8;
 467   Double_t rFirst=85.2;
 468   Double_t deltaR  = rFirst-rInner;
 469   Double_t ndistR  = (TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1])-rInner)/deltaR;
 470   Double_t driftM  = (0.5 - TMath::Abs(xyz[2]/250.))*2.;
 471   Double_t value   = TMath::Power(driftM,param[0])/ndistR;
 472   return xyz[3]*value;
 473 }
 474
 475 Double_t       AliTPCTransformation::TPCscalingROFC(Double_t *xyz, Double_t * param){
 476   //
 477   // outer field cage r distorion - proportinal to 1 over distance to the OFC
 478   // param[0] - drift polynom order
 479   // driftM   - from -1 to 1
 480   //
 481   Double_t rLast=245.8;
 482   Double_t rOuter=258.0;
 483   Double_t deltaR  = rOuter-rLast;
 484   Double_t ndistR  = (rOuter-TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]))/deltaR;
 485   Double_t driftM  = (0.5 - TMath::Abs(xyz[2]/250.))*2.;
 486   Double_t value   = TMath::Power(driftM,param[0])/ndistR;
 487   return xyz[3]*value;
 488 }
 489
 490
 491 Double_t       AliTPCTransformation::TPCdeltaFCROC(Double_t *xyz, const Double_t *const param){
 492   //
 493   // delta R(Z) ROC induced
 494   // param[0] - switch  0 - use distance to IFC - 1 - distance to IFC
 495   // param[1] - kFC scaling factor  (multiplication factor  of (OFC-IFC))
 496   // param[2] - kROC scaling factor
 497   // parameters [1] and [2] should be obtained from the electric field
 498   //            simulation
 499   //
 500   Double_t rInner=78.8;
 501   Double_t rFirst=85.2;
 502   Double_t rLast=245.8;
 503   Double_t rOuter=258.0;
 504
 505   Double_t radius  = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
 506   //calculate distance to the FC - inner or outer
 507   Double_t deltaFC = (param[0]<0.5)? TMath::Abs(radius-rFirst) : TMath::Abs(radius-rLast);
 508   deltaFC/=(rOuter-rInner);
 509   Double_t scalingFC = 1./(1.+deltaFC/(param[1]));
 510   //
 511   Double_t drift = 1.-TMath::Abs(xyz[2]/250.);  // normalized drift length
 512   Double_t scalingROC = (1.-1./(1.+drift/param[2]));
 513   //
 514   return xyz[3]*scalingFC*scalingROC;
 515 }
 516
 517
 518 Double_t       AliTPCTransformation::TPCdeltaFCCE(Double_t *xyz, const Double_t *const param){
 519   //
 520   // delta R(Z) CE (central electrode) induced
 521   // param[0] - switch  0 - use distance to IFC - 1 - distance to IFC
 522   // param[1] - kFC scaling factor  (multiplication factor  of (OFC-IFC))
 523   // param[2] - kCE scaling factor
 524   // parameters [1] and [2] should be obtained from the electric field
 525   //            simulation
 526   Double_t rInner=78.8;
 527   Double_t rFirst=85.2;
 528   Double_t rLast =245.8;
 529   Double_t rOuter=258.0;
 530
 531   Double_t radius  = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
 532   //calculate distance to the FC - inner or outer
 533   Double_t deltaFC = (param[0]<0.5)? TMath::Abs(radius-rFirst) : TMath::Abs(radius-rLast);
 534   deltaFC/=(rOuter-rInner);
 535   Double_t scalingFC = 1./(1.+deltaFC/(param[1]));
 536   //
 537   Double_t drift     = 1.-TMath::Abs(xyz[2]/250.);  // normalized drift length
 538   Double_t scalingCE = 1/(1.+(1.-drift)/param[2]);  //
 539   //
 540   return xyz[3]*scalingFC*scalingCE;
 541 }
 542
 543
 544
 545
 546
 547
 548
 549 //
 550 // TPC sector local misalignment
 551 //
 552 //
 553 //
 554 Double_t AliTPCTransformation:: TPClocaldLxdGX(Double_t *xyz, const Double_t *const param){
 555   //
 556   // xyz - [0..2] - position
 557   //       [3]    - scale parameter
 558   //       [4]    - volID
 559   // param[0]= n  - cos(n *alpha)
 560   // param[1]= n  - sin(n *alpha)
 561   // param[2]     - indication - 0 - the same for IROC OROC 1 - opposite
 562   // return delta in global coordiante system
 563   //
 564   Int_t    sector = TMath::Nint(xyz[4]);
 565   Double_t alpha  = TMath::Pi()*(sector+0.5)/9;
 566   Double_t ca     = TMath::Cos(alpha);
 567   //  Double_t sa     = TMath::Sin(alpha);
 568   const Double_t xIROCOROC = 133.4;
 569   Double_t factor = xyz[3];
 570   if (param[0]>0)  factor*=TMath::Cos(alpha*param[0]);
 571   if (param[1]>0)  factor*=TMath::Sin(alpha*param[1]);
 572   if (param[2]>0.5 && TMath::Sqrt(xyz[1]*xyz[1]+xyz[0]*xyz[0])>xIROCOROC) factor*=-1;
 573   return ca*factor;
 574 }
 575
 576 Double_t AliTPCTransformation::TPClocaldLxdGY(Double_t *xyz, const Double_t *const param){
 577   //
 578   // xyz - [0..2] - position
 579   //       [3]    - scale parameter
 580   //       [4]    - volID
 581   // param[0]= n  - cos(n *alpha)
 582   // param[1]= n  - sin(n *alpha)
 583   // param[2]     - indication - 0 - the same for IROC OROC 1 - opposite
 584   // return delta in global coordiante system
 585   //
 586   Int_t    sector = TMath::Nint(xyz[4]);
 587   Double_t alpha  = TMath::Pi()*(sector+0.5)/9;
 588   //Double_t ca     = TMath::Cos(alpha);
 589   Double_t sa     = TMath::Sin(alpha);
 590   const Double_t xIROCOROC = 133.4;
 591   Double_t factor = xyz[3];
 592   if (param[0]>0)  factor*=TMath::Cos(alpha*param[0]);
 593   if (param[1]>0)  factor*=TMath::Sin(alpha*param[1]);
 594   if (param[2]>0.5 && TMath::Sqrt(xyz[1]*xyz[1]+xyz[0]*xyz[0])>xIROCOROC) factor*=-1;
 595   return   sa*factor;
 596 }
 597
 598 Double_t AliTPCTransformation:: TPClocaldLydGX(Double_t *xyz, const Double_t *const param){
 599   //
 600   // xyz - [0..2] - position
 601   //       [3]    - scale parameter
 602   //       [4]    - volID
 603   // param[0]= n  - cos(n *alpha)
 604   // param[1]= n  - sin(n *alpha)
 605   // param[2]     - indication - 0 - the same for IROC OROC 1 - opposite
 606   // return delta in global coordiante system
 607   //
 608   Int_t    sector = TMath::Nint(xyz[4]);
 609   Double_t alpha  = TMath::Pi()*(sector+0.5)/9;
 610   //Double_t ca     = TMath::Cos(alpha);
 611   Double_t sa     = TMath::Sin(alpha);
 612   const Double_t xIROCOROC = 133.4;
 613   Double_t factor = xyz[3];
 614   if (param[0]>0)  factor*=TMath::Cos(alpha*param[0]);
 615   if (param[1]>0)  factor*=TMath::Sin(alpha*param[1]);
 616   if (param[2]>0.5 && TMath::Sqrt(xyz[1]*xyz[1]+xyz[0]*xyz[0])>xIROCOROC) factor*=-1;
 617   return            -sa*factor;
 618 }
 619
 620 Double_t AliTPCTransformation::TPClocaldLydGY(Double_t *xyz, const Double_t *const param){
 621   //
 622   // xyz - [0..2] - position
 623   //       [3]    - scale parameter
 624   //       [4]    - volID
 625   // param[0]= n  - cos(n *alpha)
 626   // param[1]= n  - sin(n *alpha)
 627   // param[2]     - indication - 0 - the same for IROC OROC 1 - opposite
 628   // return delta in global coordiante system
 629   //
 630   Int_t    sector = TMath::Nint(xyz[4]);
 631   Double_t alpha  = TMath::Pi()*(sector+0.5)/9;
 632   Double_t ca     = TMath::Cos(alpha);
 633   //Double_t sa     = TMath::Sin(alpha);
 634   const Double_t xIROCOROC = 133.4;
 635   Double_t factor = xyz[3];
 636   if (param[0]>0)  factor*=TMath::Cos(alpha*param[0]);
 637   if (param[1]>0)  factor*=TMath::Sin(alpha*param[1]);
 638   if (param[2]>0.5 && TMath::Sqrt(xyz[1]*xyz[1]+xyz[0]*xyz[0])>xIROCOROC) factor*=-1;
 639   return   ca*factor;
 640 }
 641
 642
 643 Double_t AliTPCTransformation::TPClocaldRzdGX(Double_t *xyz, const Double_t *const param){
 644   //
 645   // xyz - [0..2] - position
 646   //       [3]    - scale parameter - rotation angle in mrad
 647   //       [4]    - volID
 648   // param[0]= n  - cos(n *alpha)
 649   // param[1]= n  - sin(n *alpha)
 650   // param[2]     - indication - 0 - the same for IROC OROC 1 - opposite
 651   // return delta in global coordiante system
 652   //
 653   Int_t    sector = TMath::Nint(xyz[4]);
 654   Double_t alpha  = TMath::Pi()*(sector+0.5)/9;
 655   Double_t ca     = TMath::Cos(alpha);
 656   Double_t sa     = TMath::Sin(alpha);
 657   Double_t lx     =  xyz[0]*ca + xyz[1]*sa;
 658   Double_t ly     = -xyz[0]*sa + xyz[1]*ca;
 659   //
 660   const Double_t xIROCOROC = 133.4;
 661   lx-=xIROCOROC;
 662   Double_t rot      =  xyz[3]*0.001;    // rotation in mrad
 663   if (param[0]>0)  rot*=TMath::Cos(alpha*param[0]);
 664   if (param[1]>0)  rot*=TMath::Sin(alpha*param[1]);
 665   if (param[2]>0.5 && lx>0) rot*=-1;
 666   //
 667   Double_t dlxR     =  - ly*rot;
 668   Double_t dlyR     =    lx*rot;
 669   Double_t dgxR     =  dlxR*ca - dlyR*sa;
 670   //Double_t dgyR     =  dlxR*sa + dlyR*ca;
 671   return  dgxR;
 672 }
 673
 674 Double_t AliTPCTransformation::TPClocaldRzdGY(Double_t *xyz, const Double_t *const param){
 675   //
 676   // xyz - [0..2] - position
 677   //       [3]    - scale parameter - rotation angle in mrad
 678   //       [4]    - volID
 679   // param[0]= n  - cos(n *alpha)
 680   // param[1]= n  - sin(n *alpha)
 681   // param[2]     - indication - 0 - the same for IROC OROC 1 - opposite
 682   // return delta in global coordiante system
 683   //
 684   Int_t    sector = TMath::Nint(xyz[4]);
 685   Double_t alpha  = TMath::Pi()*(sector+0.5)/9;
 686   Double_t ca     = TMath::Cos(alpha);
 687   Double_t sa     = TMath::Sin(alpha);
 688   Double_t lx     =  xyz[0]*ca + xyz[1]*sa;
 689   Double_t ly     = -xyz[0]*sa + xyz[1]*ca;
 690   //
 691   const Double_t xIROCOROC = 133.4;
 692   lx-=xIROCOROC;
 693   Double_t rot      =  xyz[3]*0.001;    // rotation in mrad
 694   if (param[0]>0)  rot*=TMath::Cos(alpha*param[0]);
 695   if (param[1]>0)  rot*=TMath::Sin(alpha*param[1]);
 696   if (param[2]>0.5 && lx>0) rot*=-1;
 697   Double_t dlxR     =  - ly*rot;
 698   Double_t dlyR     =    lx*rot;
 699   //Double_t dgxR     =  dlxR*ca - dlyR*sa;
 700   Double_t dgyR     =  dlxR*sa + dlyR*ca;
 701   return  dgyR;
 702 }
 703
 704
 705
 706 Double_t        AliTPCTransformation::TPCDeltaZMediumLong(Double_t *xyz, Double_t * /*param*/){
 707   //
 708   // xyz - [0..2] - position
 709   //        [3]    - scale parameter
 710   //        [4]    - volID
 711   // return delta in global coordinate system
 712   //
 713   Int_t    sector = TMath::Nint(xyz[4]);
 714   Double_t signZ  = (sector%36<18) ? 1: -1;  // drift direction
 715   if    (sector<36) return 0;
 716   //
 717   const Double_t radiusLong = 198.1;
 718   //
 719   Double_t alpha  = TMath::Pi()*(sector+0.5)/9;
 720   Double_t ca     = TMath::Cos(alpha);
 721   Double_t sa     = TMath::Sin(alpha);
 722   Double_t lx     =  xyz[0]*ca + xyz[1]*sa;
 723   Double_t sign   = (lx<radiusLong) ? 1:-1;
 724   return xyz[3]*sign*signZ;
 725 }
 726
 727 Double_t        AliTPCTransformation::TPCDeltaZ(Double_t *xyz, const Double_t *const param){
 728   //
 729   // xyz - [0..2] - position
 730   //        [3]    - scale parameter
 731   //        [4]    - volID
 732   // return delta in global coordiante system
 733   //
 734   Int_t    sector = TMath::Nint(xyz[4]);
 735   Double_t delta  = (sector%36<18) ? 1: -1;  // drift direction
 736   Double_t alpha  = TMath::Pi()*(sector+0.5)/9;
 737   Double_t ca     = TMath::Cos(alpha);
 738   Double_t sa     = TMath::Sin(alpha);
 739   Double_t lx     =  xyz[0]*ca + xyz[1]*sa;
 740   //
 741   const Double_t xIROCOROC = 133.4;
 742   if (param[0]>0) delta     *= TMath::Cos(param[0]*alpha);
 743   if (param[1]>0) delta     *= TMath::Sin(param[1]*alpha);
 744   if (param[2]>0.5 && lx >xIROCOROC) delta *=-1;
 745   return delta*xyz[3];     // IROC shift
 746 }
 747
 748
 749 Double_t       AliTPCTransformation::TPCTiltingZ(Double_t *xyz, const Double_t *const param){
 750   // xyz - [0..2] - position
 751   //        [3]    - scale parameter
 752   //        [4]    - volID
 753   // param[0]      - n for cos
 754   // param[1]      - n for sin
 755   // param[2]      - IROC-ORC relative (if >0.5 )
 756   // return delta in global coordinate system
 757   const Double_t rFirst=85.2;
 758   const Double_t rLast =245.8;
 759   const Double_t xIROCOROC = 133.4;
 760   //
 761   Int_t    sector = TMath::Nint(xyz[4]);
 762   Double_t alpha  = TMath::Pi()*(sector+0.5)/9;
 763   Double_t ca     = TMath::Cos(alpha);
 764   Double_t sa     = TMath::Sin(alpha);
 765   Double_t lx     =  xyz[0]*ca + xyz[1]*sa;
 766   Double_t deltaR = 2.0*(lx-xIROCOROC)/(rLast-rFirst);
 767   if (param[0]>0) deltaR *= TMath::Cos(param[0]*alpha);
 768   if (param[1]>0) deltaR *= TMath::Sin(param[1]*alpha);
 769   if (param[2]>0.5 && lx >xIROCOROC) deltaR *=-1;
 770   return deltaR*xyz[3];
 771 }
 772