3 Class AliTPCtransformation:
4 Should represent general non linear transformation. Currently tune for TPConly.
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
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
17 2. The transformation should work only for given volIDs and detector IDs
18 Currently Bitmask is used for filtering
21 Transformation - naming convention:
25 XXX - detector if detector specific
26 local - if local transforamtion
27 YYY - type of transformation
28 ZZZ - return type of transformation
38 #include "TLinearFitter.h"
40 #include "TObjString.h"
42 #include "TTreeStream.h"
43 #include "AliTrackPointArray.h"
45 #include "AliTPCTransformation.h"
47 ClassImp(AliTPCTransformation)
50 AliTPCTransformation::GenFuncG AliTPCTransformation::fgFormulas[10000];
51 TObjArray* AliTPCTransformation::fgFormulasName = new TObjArray(10000);
54 void AliTPCTransformation::RegisterFormula(const char * name, GenFuncG formula){
56 // Add Formula to the list of formulas
58 Int_t last= fgFormulasName->GetEntries();
59 fgFormulasName->AddLast(new TObjString(name));
60 fgFormulas[last]=formula;
63 Int_t AliTPCTransformation::BuildBasicFormulas(){
67 //build list of basic TPC formulas - corrections
69 RegisterFormula("TPCscalingRPol",(GenFuncG)(AliTPCTransformation::TPCscalingRPol));
70 RegisterFormula("TPCscalingRIFC",(GenFuncG)(AliTPCTransformation::TPCscalingRIFC));
71 RegisterFormula("TPCscalingROFC",(GenFuncG)(AliTPCTransformation::TPCscalingROFC));
73 RegisterFormula("TPCdeltaFCROC",(GenFuncG)(AliTPCTransformation::TPCdeltaFCROC));
74 RegisterFormula("TPCdeltaFCCE",(GenFuncG)(AliTPCTransformation::TPCdeltaFCCE));
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));
82 // TPC Local X and Y misalignment + rotation
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));
94 RegisterFormula("TPCDeltaZ",(GenFuncG)(AliTPCTransformation::TPCDeltaZ));
95 RegisterFormula("TPCDeltaZMediumLong",(GenFuncG)(AliTPCTransformation::TPCDeltaZMediumLong));
96 RegisterFormula("TPCTiltingZ",(GenFuncG)(AliTPCTransformation::TPCTiltingZ));
100 AliTPCTransformation::GenFuncG AliTPCTransformation::FindFormula(const char * name){
102 // find formula - if registered
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];
114 Double_t AliTPCTransformation::Eval(const char * name, const Double_t*x,const Double_t*par){
116 // Only for test purposes - very slow
118 GenFuncG fun = FindFormula(name);
124 AliTPCTransformation::AliTPCTransformation():
126 fNameX(0), // x formula name
127 fNameY(0), // y formula name
128 fNameZ(0), // z formula name
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
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
146 // default constructor
152 AliTPCTransformation::AliTPCTransformation(const char *name, TBits *mask, const char *fx, const char *fy, const char *fz, Int_t coordSystem):
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
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
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
172 // non default constructor
174 if (fx) fNameX= new TString(fx);
175 if (fy) fNameY= new TString(fy);
176 if (fz) fNameZ= new TString(fz);
180 AliTPCTransformation::AliTPCTransformation(const AliTPCTransformation&trafo):
182 fNameX(0), // x formula name
183 fNameY(0), // y formula name
184 fNameZ(0), // z formula name
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
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
203 if (trafo.fNameX) fNameX = new TString(*(trafo.fNameX));
204 if (trafo.fNameY) fNameY = new TString(*(trafo.fNameY));
205 if (trafo.fNameZ) fNameZ = new TString(*(trafo.fNameZ));
206 if (trafo.fBitMask) fBitMask = new TBits(*(trafo.fBitMask));
209 AliTPCTransformation::~AliTPCTransformation(){
219 void AliTPCTransformation::SetParams(Double_t param, Double_t sigma, Double_t sigma2Time, const TVectorD *const fixedParams){
221 // Set parameters of transformation
226 fSigma2Time = sigma2Time;
227 if (fFixedParam) delete fFixedParam;
228 fFixedParam = new TVectorD(*fixedParams);
233 Bool_t AliTPCTransformation::Init(){
235 // associate formulas with pointer to the function
239 fFormulaX=FindFormula(fNameX->Data());
240 if (fFormulaX==0) isOK=kFALSE;
243 fFormulaY=FindFormula(fNameY->Data());
244 if (fFormulaY==0) isOK=kFALSE;
247 fFormulaZ=FindFormula(fNameZ->Data());
248 if (!fFormulaZ) isOK=kFALSE;
255 TBits * AliTPCTransformation::BitsSide(Bool_t aside){
257 // Set bits for given side
259 TBits * bits = new TBits(72);
260 for (Int_t i=0; i<72;i++){
261 if (i%36<18 && aside) (*bits)[i]=kTRUE;
262 if (i%36<18 && (!aside)) (*bits)[i]=kFALSE;
263 if (i%36>=18 && aside) (*bits)[i]=kFALSE;
264 if (i%36>=18 && (!aside)) (*bits)[i]=kTRUE;
269 TBits * AliTPCTransformation::BitsAll(){
271 // Set all bits to kTRUE
273 TBits * bits = new TBits(72);
274 for (Int_t i=0; i<72;i++){
280 // Double_t AliTPCTransformation::GetDeltaXYZ(Int_t coord, Int_t volID, Double_t param, Double_t x, Double_t y, Double_t z){
283 // // coord - type of coordinate
289 // if (!fIsActive) return 0;
290 // if (fBitMask && (!(*fBitMask)[volID])) return 0;
291 // Double_t xyz[5]={x,y,z, param,volID};
292 // if (fCoordSystem==0){
293 // // cartezian system
294 // if (coord==0 && fFormulaX) return fFormulaX(xyz,fFixedParam->GetMatrixArray());
295 // if (coord==1 && fFormulaY) return fFormulaY(xyz,fFixedParam->GetMatrixArray());
296 // if (coord==2 && fFormulaZ) return fFormulaZ(xyz,fFixedParam->GetMatrixArray());
298 // if (fCoordSystem==1){
299 // // cylindrical system
301 // if (fFormulaZ==0) return 0;
302 // return fFormulaZ(xyz,fFixedParam->GetMatrixArray());
304 // Double_t rrphiz[3]={0,0,0};
305 // if (fFormulaX) rrphiz[0] = fFormulaX(xyz,fFixedParam->GetMatrixArray());
306 // if (fFormulaY) rrphiz[1] = fFormulaY(xyz,fFixedParam->GetMatrixArray());
307 // Double_t alpha = TMath::ATan2(y,x);
308 // Double_t ca = TMath::Cos(alpha);
309 // Double_t sa = TMath::Sin(alpha);
310 // if (coord==0) return ca*rrphiz[0]-sa*rrphiz[1];
311 // if (coord==1) return sa*rrphiz[0]+ca*rrphiz[1];
317 Double_t AliTPCTransformation::GetDeltaXYZ(Int_t coord, Int_t volID, Double_t param, Double_t x, Double_t y, Double_t z){
320 // coord - type of coordinate
327 if (!fIsActive) return 0;
328 if (fBitMask && (!(*fBitMask)[volID])) return 0;
329 Double_t xyz[5]={x,y,z, param,volID};
330 Double_t alpha = TMath::ATan2(y,x);
331 Double_t ca = TMath::Cos(alpha);
332 Double_t sa = TMath::Sin(alpha);
334 if (fCoordSystem==0){
336 Double_t dxdydz[3]={0,0,0};
337 if(fFormulaX) dxdydz[0]=fFormulaX(xyz,fFixedParam->GetMatrixArray());
338 if(fFormulaY) dxdydz[1]=fFormulaY(xyz,fFixedParam->GetMatrixArray());
339 if(fFormulaZ) dxdydz[2]=fFormulaZ(xyz,fFixedParam->GetMatrixArray());
341 if (coord==0) return dxdydz[0];
342 if (coord==1) return dxdydz[1];
343 if (coord==2) return dxdydz[2];
344 if (coord==3) return dxdydz[0]*ca+dxdydz[1]*sa;
345 if (coord==4) return -dxdydz[0]*sa+dxdydz[1]*ca;
346 if (coord==5) return dxdydz[2];
348 if (fCoordSystem==1){
349 // cylindrical system
350 if (coord==2||coord==5) {
351 if (fFormulaZ==0) return 0;
352 return fFormulaZ(xyz,fFixedParam->GetMatrixArray());
354 Double_t rrphiz[3]={0,0,0};
355 if (fFormulaX) rrphiz[0] = fFormulaX(xyz,fFixedParam->GetMatrixArray());
356 if (fFormulaY) rrphiz[1] = fFormulaY(xyz,fFixedParam->GetMatrixArray());
357 alpha = TMath::ATan2(y,x);
358 ca = TMath::Cos(alpha);
359 sa = TMath::Sin(alpha);
360 if (coord==0) return ca*rrphiz[0]-sa*rrphiz[1];
361 if (coord==1) return sa*rrphiz[0]+ca*rrphiz[1];
362 if (coord==3) return rrphiz[0];
363 if (coord==4) return rrphiz[1];
370 Double_t AliTPCTransformation::TPCscalingRPol(Double_t *xyz, const Double_t * const param){
372 // Scaling and shift of TPC radius
373 // xyz[0..2] - global xyz of point
374 // xyz[3] - scale parameter
375 // param[0] - radial scaling power
376 // param[1] - drift scaling power
377 // radius from -1(at rInner) to 1 (rOuter)
378 // driftM from -1(at 0 drift) to 1 (250 cm drift)
380 Double_t rInner=78.8;
381 Double_t rOuter=258.0;
382 Double_t deltaR = rOuter-rInner;
383 Double_t radius = (TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1])-rInner)*2./deltaR;
384 Double_t driftM = (0.5 - TMath::Abs(xyz[2]/250.))*2.0;
385 Double_t delta = TMath::Power(radius,param[0])*TMath::Power(driftM,param[1]);
390 Double_t AliTPCTransformation::TPCscalingZDrift(Double_t *xyz, const Double_t * const param){
393 // Scaling and shift of TPC radius
394 // xyz[0..2] - global xyz of point
395 // xyz[3] - scale parameter
396 Double_t driftP = TMath::Power(1. - TMath::Abs(xyz[2]/250.), param[0]);
397 Int_t sector = TMath::Nint(xyz[4]);
398 Double_t deltaZ = (sector%36<18) ? -driftP : driftP;
399 return deltaZ*xyz[3];
402 Double_t AliTPCTransformation::TPCscalingZDriftT0(Double_t *xyz, const Double_t * const /*param*/){
405 // Z shift because time 0 offset
406 // opposite on A and C side
408 // xyz[0..2] - global xyz of point
409 // xyz[3] - scale parameter
410 Int_t sector = TMath::Nint(xyz[4]);
411 Double_t sign = (sector%36<18) ? -1 : 1;
416 Double_t AliTPCTransformation::TPCscalingZDriftGy(Double_t *xyz, const Double_t * const param){
419 // Scaling and shift of TPC radius
420 // xyz[0..2] - global xyz of point
421 // xyz[3] - scale parameter
422 Double_t driftP = TMath::Power(1. - TMath::Abs(xyz[2]/250.), param[0]);
423 Double_t gy = xyz[1]/250.;
424 Int_t sector = TMath::Nint(xyz[4]);
425 Double_t deltaZ = (sector%36<18) ? -driftP : driftP;
426 return deltaZ*xyz[3]*gy;
431 Double_t AliTPCTransformation::TPCscalingPhiLocal(Double_t *xyz, const Double_t * const param){
434 // Scaling if the local y -phi
435 // xyz[0..2] - global xyz of point
436 // xyz[3] - scale parameter
437 // value = 1 for ful drift length and parameter 1
438 Double_t alpha = TMath::ATan2(xyz[1],xyz[0]);
439 Double_t sector = TMath::Nint(9*alpha/TMath::Pi()-0.5);
440 Double_t localAlpha = (alpha-(sector+0.5)*TMath::Pi()/9.);
441 Double_t radius = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1])/250.;
443 Double_t deltaAlpha = radius*TMath::Power(2.*9.*localAlpha/TMath::Pi(),param[0]);
444 return deltaAlpha*xyz[3];
447 Double_t AliTPCTransformation::TPClocalRPhiEdge(Double_t *xyz, const Double_t *const param){
450 // Scaling if the local y -phi
451 // xyz[0..2] - global xyz of point
452 // xyz[3] - scale parameter
453 // param[0] - dedge offset - should be around gap size/2.
454 // param[1] - dedge factor - should be around gap size/2.
455 Double_t alpha = TMath::ATan2(xyz[1],xyz[0]);
456 Double_t sector = TMath::Nint(9*alpha/TMath::Pi()-0.5);
457 Double_t localAlpha = (alpha-(sector+0.5)*TMath::Pi()/9.);
458 Double_t radius = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
459 Double_t deltaAlpha = TMath::Pi()/18.-TMath::Abs(localAlpha);
460 Double_t distEdge = (deltaAlpha*radius);
461 Double_t factor = 1./(1.+(distEdge-param[0])/param[1]);
462 return factor*xyz[3]*((localAlpha>0)? -1.:1.);
466 Double_t AliTPCTransformation::TPCscalingRIFC(Double_t *xyz, const Double_t * const param){
468 // inner field cage r distorion - proportinal to 1 over distance to the IFC
469 // param[0] - drift polynom order
470 // distortion at first pad row - is normalized to
471 Double_t rInner=78.8;
472 Double_t rFirst=85.2;
473 Double_t deltaR = rFirst-rInner;
474 Double_t ndistR = (TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1])-rInner)/deltaR;
475 Double_t driftM = (0.5 - TMath::Abs(xyz[2]/250.))*2.;
476 Double_t value = TMath::Power(driftM,param[0])/ndistR;
480 Double_t AliTPCTransformation::TPCscalingROFC(Double_t *xyz, const Double_t * const param){
482 // outer field cage r distorion - proportinal to 1 over distance to the OFC
483 // param[0] - drift polynom order
484 // driftM - from -1 to 1
486 Double_t rLast=245.8;
487 Double_t rOuter=258.0;
488 Double_t deltaR = rOuter-rLast;
489 Double_t ndistR = (rOuter-TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]))/deltaR;
490 Double_t driftM = (0.5 - TMath::Abs(xyz[2]/250.))*2.;
491 Double_t value = TMath::Power(driftM,param[0])/ndistR;
496 Double_t AliTPCTransformation::TPCdeltaFCROC(Double_t *xyz, const Double_t *const param){
498 // delta R(Z) ROC induced
499 // param[0] - switch 0 - use distance to IFC - 1 - distance to IFC
500 // param[1] - kFC scaling factor (multiplication factor of (OFC-IFC))
501 // param[2] - kROC scaling factor
502 // parameters [1] and [2] should be obtained from the electric field
505 Double_t rInner=78.8;
506 Double_t rFirst=85.2;
507 Double_t rLast=245.8;
508 Double_t rOuter=258.0;
510 Double_t radius = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
511 //calculate distance to the FC - inner or outer
512 Double_t deltaFC = (param[0]<0.5)? TMath::Abs(radius-rFirst) : TMath::Abs(radius-rLast);
513 deltaFC/=(rOuter-rInner);
514 Double_t scalingFC = 1./(1.+deltaFC/(param[1]));
516 Double_t drift = 1.-TMath::Abs(xyz[2]/250.); // normalized drift length
517 Double_t scalingROC = (1.-1./(1.+drift/param[2]));
519 return xyz[3]*scalingFC*scalingROC;
523 Double_t AliTPCTransformation::TPCdeltaFCCE(Double_t *xyz, const Double_t *const param){
525 // delta R(Z) CE (central electrode) induced
526 // param[0] - switch 0 - use distance to IFC - 1 - distance to IFC
527 // param[1] - kFC scaling factor (multiplication factor of (OFC-IFC))
528 // param[2] - kCE scaling factor
529 // parameters [1] and [2] should be obtained from the electric field
531 Double_t rInner=78.8;
532 Double_t rFirst=85.2;
533 Double_t rLast =245.8;
534 Double_t rOuter=258.0;
536 Double_t radius = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
537 //calculate distance to the FC - inner or outer
538 Double_t deltaFC = (param[0]<0.5)? TMath::Abs(radius-rFirst) : TMath::Abs(radius-rLast);
539 deltaFC/=(rOuter-rInner);
540 Double_t scalingFC = 1./(1.+deltaFC/(param[1]));
542 Double_t drift = 1.-TMath::Abs(xyz[2]/250.); // normalized drift length
543 Double_t scalingCE = 1/(1.+(1.-drift)/param[2]); //
545 return xyz[3]*scalingFC*scalingCE;
555 // TPC sector local misalignment
559 Double_t AliTPCTransformation:: TPClocaldLxdGX(Double_t *xyz, const Double_t *const param){
561 // xyz - [0..2] - position
562 // [3] - scale parameter
564 // param[0]= n - cos(n *alpha)
565 // param[1]= n - sin(n *alpha)
566 // param[2] - indication - 0 - the same for IROC OROC 1 - opposite
567 // return delta in global coordiante system
569 Int_t sector = TMath::Nint(xyz[4]);
570 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
571 Double_t ca = TMath::Cos(alpha);
572 // Double_t sa = TMath::Sin(alpha);
573 const Double_t xIROCOROC = 133.4;
574 Double_t factor = xyz[3];
575 if (param[0]>0) factor*=TMath::Cos(alpha*param[0]);
576 if (param[1]>0) factor*=TMath::Sin(alpha*param[1]);
577 if (param[2]>0.5 && TMath::Sqrt(xyz[1]*xyz[1]+xyz[0]*xyz[0])>xIROCOROC) factor*=-1;
581 Double_t AliTPCTransformation::TPClocaldLxdGY(Double_t *xyz, const Double_t *const param){
583 // xyz - [0..2] - position
584 // [3] - scale parameter
586 // param[0]= n - cos(n *alpha)
587 // param[1]= n - sin(n *alpha)
588 // param[2] - indication - 0 - the same for IROC OROC 1 - opposite
589 // return delta in global coordiante system
591 Int_t sector = TMath::Nint(xyz[4]);
592 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
593 //Double_t ca = TMath::Cos(alpha);
594 Double_t sa = TMath::Sin(alpha);
595 const Double_t xIROCOROC = 133.4;
596 Double_t factor = xyz[3];
597 if (param[0]>0) factor*=TMath::Cos(alpha*param[0]);
598 if (param[1]>0) factor*=TMath::Sin(alpha*param[1]);
599 if (param[2]>0.5 && TMath::Sqrt(xyz[1]*xyz[1]+xyz[0]*xyz[0])>xIROCOROC) factor*=-1;
603 Double_t AliTPCTransformation:: TPClocaldLydGX(Double_t *xyz, const Double_t *const param){
605 // xyz - [0..2] - position
606 // [3] - scale parameter
608 // param[0]= n - cos(n *alpha)
609 // param[1]= n - sin(n *alpha)
610 // param[2] - indication - 0 - the same for IROC OROC 1 - opposite
611 // return delta in global coordiante system
613 Int_t sector = TMath::Nint(xyz[4]);
614 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
615 //Double_t ca = TMath::Cos(alpha);
616 Double_t sa = TMath::Sin(alpha);
617 const Double_t xIROCOROC = 133.4;
618 Double_t factor = xyz[3];
619 if (param[0]>0) factor*=TMath::Cos(alpha*param[0]);
620 if (param[1]>0) factor*=TMath::Sin(alpha*param[1]);
621 if (param[2]>0.5 && TMath::Sqrt(xyz[1]*xyz[1]+xyz[0]*xyz[0])>xIROCOROC) factor*=-1;
625 Double_t AliTPCTransformation::TPClocaldLydGY(Double_t *xyz, const Double_t *const param){
627 // xyz - [0..2] - position
628 // [3] - scale parameter
630 // param[0]= n - cos(n *alpha)
631 // param[1]= n - sin(n *alpha)
632 // param[2] - indication - 0 - the same for IROC OROC 1 - opposite
633 // return delta in global coordiante system
635 Int_t sector = TMath::Nint(xyz[4]);
636 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
637 Double_t ca = TMath::Cos(alpha);
638 //Double_t sa = TMath::Sin(alpha);
639 const Double_t xIROCOROC = 133.4;
640 Double_t factor = xyz[3];
641 if (param[0]>0) factor*=TMath::Cos(alpha*param[0]);
642 if (param[1]>0) factor*=TMath::Sin(alpha*param[1]);
643 if (param[2]>0.5 && TMath::Sqrt(xyz[1]*xyz[1]+xyz[0]*xyz[0])>xIROCOROC) factor*=-1;
648 Double_t AliTPCTransformation::TPClocaldRzdGX(Double_t *xyz, const Double_t *const param){
650 // xyz - [0..2] - position
651 // [3] - scale parameter - rotation angle in mrad
653 // param[0]= n - cos(n *alpha)
654 // param[1]= n - sin(n *alpha)
655 // param[2] - indication - 0 - the same for IROC OROC 1 - opposite
656 // return delta in global coordiante system
658 Int_t sector = TMath::Nint(xyz[4]);
659 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
660 Double_t ca = TMath::Cos(alpha);
661 Double_t sa = TMath::Sin(alpha);
662 Double_t lx = xyz[0]*ca + xyz[1]*sa;
663 Double_t ly = -xyz[0]*sa + xyz[1]*ca;
665 const Double_t xIROCOROC = 133.4;
667 Double_t rot = xyz[3]*0.001; // rotation in mrad
668 if (param[0]>0) rot*=TMath::Cos(alpha*param[0]);
669 if (param[1]>0) rot*=TMath::Sin(alpha*param[1]);
670 if (param[2]>0.5 && lx>0) rot*=-1;
672 Double_t dlxR = - ly*rot;
673 Double_t dlyR = lx*rot;
674 Double_t dgxR = dlxR*ca - dlyR*sa;
675 //Double_t dgyR = dlxR*sa + dlyR*ca;
679 Double_t AliTPCTransformation::TPClocaldRzdGY(Double_t *xyz, const Double_t *const param){
681 // xyz - [0..2] - position
682 // [3] - scale parameter - rotation angle in mrad
684 // param[0]= n - cos(n *alpha)
685 // param[1]= n - sin(n *alpha)
686 // param[2] - indication - 0 - the same for IROC OROC 1 - opposite
687 // return delta in global coordiante system
689 Int_t sector = TMath::Nint(xyz[4]);
690 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
691 Double_t ca = TMath::Cos(alpha);
692 Double_t sa = TMath::Sin(alpha);
693 Double_t lx = xyz[0]*ca + xyz[1]*sa;
694 Double_t ly = -xyz[0]*sa + xyz[1]*ca;
696 const Double_t xIROCOROC = 133.4;
698 Double_t rot = xyz[3]*0.001; // rotation in mrad
699 if (param[0]>0) rot*=TMath::Cos(alpha*param[0]);
700 if (param[1]>0) rot*=TMath::Sin(alpha*param[1]);
701 if (param[2]>0.5 && lx>0) rot*=-1;
702 Double_t dlxR = - ly*rot;
703 Double_t dlyR = lx*rot;
704 //Double_t dgxR = dlxR*ca - dlyR*sa;
705 Double_t dgyR = dlxR*sa + dlyR*ca;
711 Double_t AliTPCTransformation::TPCDeltaZMediumLong(Double_t *xyz, Double_t * /*param*/){
713 // xyz - [0..2] - position
714 // [3] - scale parameter
716 // return delta in global coordinate system
718 Int_t sector = TMath::Nint(xyz[4]);
719 Double_t signZ = (sector%36<18) ? 1: -1; // drift direction
720 if (sector<36) return 0;
722 const Double_t radiusLong = 198.1;
724 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
725 Double_t ca = TMath::Cos(alpha);
726 Double_t sa = TMath::Sin(alpha);
727 Double_t lx = xyz[0]*ca + xyz[1]*sa;
728 Double_t sign = (lx<radiusLong) ? 1:-1;
729 return xyz[3]*sign*signZ;
732 Double_t AliTPCTransformation::TPCDeltaZ(Double_t *xyz, const Double_t *const param){
734 // xyz - [0..2] - position
735 // [3] - scale parameter
737 // return delta in global coordiante system
739 Int_t sector = TMath::Nint(xyz[4]);
740 Double_t delta = (sector%36<18) ? 1: -1; // drift direction
741 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
742 Double_t ca = TMath::Cos(alpha);
743 Double_t sa = TMath::Sin(alpha);
744 Double_t lx = xyz[0]*ca + xyz[1]*sa;
746 const Double_t xIROCOROC = 133.4;
747 if (param[0]>0) delta *= TMath::Cos(param[0]*alpha);
748 if (param[1]>0) delta *= TMath::Sin(param[1]*alpha);
749 if (param[2]>0.5 && lx >xIROCOROC) delta *=-1;
750 return delta*xyz[3]; // IROC shift
754 Double_t AliTPCTransformation::TPCTiltingZ(Double_t *xyz, const Double_t *const param){
755 // xyz - [0..2] - position
756 // [3] - scale parameter
758 // param[0] - n for cos
759 // param[1] - n for sin
760 // param[2] - IROC-ORC relative (if >0.5 )
761 // return delta in global coordinate system
762 const Double_t rFirst=85.2;
763 const Double_t rLast =245.8;
764 const Double_t xIROCOROC = 133.4;
766 Int_t sector = TMath::Nint(xyz[4]);
767 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
768 Double_t ca = TMath::Cos(alpha);
769 Double_t sa = TMath::Sin(alpha);
770 Double_t lx = xyz[0]*ca + xyz[1]*sa;
771 Double_t deltaR = 2.0*(lx-xIROCOROC)/(rLast-rFirst);
772 if (param[0]>0) deltaR *= TMath::Cos(param[0]*alpha);
773 if (param[1]>0) deltaR *= TMath::Sin(param[1]*alpha);
774 if (param[2]>0.5 && lx >xIROCOROC) deltaR *=-1;
775 return deltaR*xyz[3];