]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TPC/AliTPCTransformation.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Updated list of files
[u/mrichter/AliRoot.git] / TPC / AliTPCTransformation.cxx
CommitLineData
94685752 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
22*/
23
24#include <string.h>
25#include "TRandom.h"
26#include "TMath.h"
27#include "TBits.h"
28#include "TFormula.h"
29#include "TF1.h"
30#include "TLinearFitter.h"
31#include "TFile.h"
32#include "TObjString.h"
33
34#include "TTreeStream.h"
35#include "AliTrackPointArray.h"
36#include "AliLog.h"
37#include "AliTPCTransformation.h"
38
39ClassImp(AliTPCTransformation)
40
41
42AliTPCTransformation::GenFuncG AliTPCTransformation::fgFormulas[10000];
43TObjArray* AliTPCTransformation::fgFormulasName = new TObjArray(10000);
44
45
46void AliTPCTransformation::RegisterFormula(const char * name, GenFuncG formula){
47 //
48 // Add Formula to the list of formulas
49 //
50 Int_t last= fgFormulasName->GetEntries();
51 fgFormulasName->AddLast(new TObjString(name));
52 fgFormulas[last]=formula;
53}
54
55Int_t AliTPCTransformation::BuildBasicFormulas(){
56 //
6d438146 57
58 //
59 //build list of basic TPC formulas - corrections
94685752 60 //
61 RegisterFormula("TPCscalingRPol",(GenFuncG)(AliTPCTransformation::TPCscalingRPol));
6d438146 62 RegisterFormula("TPCscalingRIFC",(GenFuncG)(AliTPCTransformation::TPCscalingRIFC));
63 RegisterFormula("TPCscalingROFC",(GenFuncG)(AliTPCTransformation::TPCscalingROFC));
64 //
94685752 65 RegisterFormula("TPCscalingZDr",(GenFuncG)(AliTPCTransformation::TPCscalingZDr));
66 RegisterFormula("TPCscalingPhiLocal",(GenFuncG)(AliTPCTransformation::TPCscalingPhiLocal));
0e9efcbe 67 //
68 // TPC Local X and Y misalignment
69 //
70 RegisterFormula("TPClocaldLxdGX",(GenFuncG)(AliTPCTransformation::TPClocaldLxdGX));
71 RegisterFormula("TPClocaldLxdGY",(GenFuncG)(AliTPCTransformation::TPClocaldLxdGY));
72 RegisterFormula("TPClocaldLydGX",(GenFuncG)(AliTPCTransformation::TPClocaldLydGX));
73 RegisterFormula("TPClocaldLydGY",(GenFuncG)(AliTPCTransformation::TPClocaldLydGY));
74 //
75 // Z offset
76 //
77 RegisterFormula("TPCDeltaZ",(GenFuncG)(AliTPCTransformation::TPCDeltaZ));
78 RegisterFormula("TPCDeltaZMediumLong",(GenFuncG)(AliTPCTransformation::TPCDeltaZMediumLong));
79 RegisterFormula("TPCTiltingZ",(GenFuncG)(AliTPCTransformation::TPCTiltingZ));
94685752 80 return 0;
81}
82
83AliTPCTransformation::GenFuncG AliTPCTransformation::FindFormula(const char * name){
84 //
85 // find formula - if registered
86 //
87 if (fgFormulasName->FindObject(name)==0) return 0;
88 Int_t entries = fgFormulasName->GetEntries();
89 for (Int_t i=0;i<entries;i++){
90 if (strcmp(fgFormulasName->At(i)->GetName(), name)==0){
91 return fgFormulas[i];
92 }
93 }
94 return 0;
95}
96
97Double_t AliTPCTransformation::Eval(const char * name, const Double_t*x,const Double_t*par){
98 //
99 // Only for test purposes - very slow
100 //
101 GenFuncG fun = FindFormula(name);
102 if (!fun) return 0;
103 return fun(x,par);
104}
105
106
107AliTPCTransformation::AliTPCTransformation():
108 TNamed(),
109 fNameX(0), // x formula name
110 fNameY(0), // y formula name
111 fNameZ(0), // z formula name
112 //
113 fBitMask(0), // bitmaps - transformation only for specified volID
114 fCoordSystem(0), // coord system of output deltas
115 fParam(0), // free parameter of transformation
116 fSigma(0), // uncertainty of the parameter
6d438146 117 fSigma2Time(0), // change of the error in time - (For kalman filter)
94685752 118 fFixedParam(0), // fixed parameters of tranformation
0e9efcbe 119 fIsActive(kTRUE), // swith - On/Off
94685752 120 //
121 fInit(kFALSE), // initialization flag - set to kTRUE if corresponding formulas found
122 fFormulaX(0), // x formula - pointer to the function
123 fFormulaY(0), // y formula - pointer to the function
124 fFormulaZ(0) // z formula - pointer to the function
125 //
126{
127 //
128 // default constructor
129 //
130}
131
132
133
6d438146 134AliTPCTransformation::AliTPCTransformation(const char *name, TBits *mask, const char *fx, const char *fy, const char *fz, Int_t coordSystem):
94685752 135 TNamed(name,name),
136 fNameX(0), // x formula name
137 fNameY(0), // y formula name
138 fNameZ(0), // z formula name
139 fBitMask(mask), // bitmaps - transformation only for specified volID
140 fCoordSystem(coordSystem), // coordinate system of output deltas
6d438146 141 fParam(0), // free parameter of transformation
142 fSigma(0),
143 fSigma2Time(0), // change of sigma in time
94685752 144 fFixedParam(0), // fixed parameters of tranformation
0e9efcbe 145 fIsActive(kTRUE), // swith - On/Off
94685752 146 //
147 fInit(kFALSE), // initialization flag - set to kTRUE if corresponding formulas found
148 fFormulaX(0), // x formula - pointer to the function
149 fFormulaY(0), // y formula - pointer to the function
150 fFormulaZ(0) // z formula - pointer to the function
151{
152 //
153 // non default constructor
154 //
155 if (fx) fNameX= new TString(fx);
156 if (fy) fNameY= new TString(fy);
157 if (fz) fNameZ= new TString(fz);
94685752 158 Init();
159}
160
6d438146 161AliTPCTransformation::AliTPCTransformation(const AliTPCTransformation&trafo):
162 TNamed(trafo),
163 fNameX(0), // x formula name
164 fNameY(0), // y formula name
165 fNameZ(0), // z formula name
166 //
167 fBitMask(0), // bitmaps - transformation only for specified volID
168 fCoordSystem(0), // coord system of output deltas
169 fParam(trafo.fParam), // free parameter of transformation
170 fSigma(trafo.fSigma), // uncertainty of the parameter
171 fSigma2Time(trafo.fSigma2Time), // change of the error in time - (For kalman filter)
e55e5512 172 fFixedParam(0), // fixed parameters of tranformation
173 fIsActive(trafo.fIsActive), // swith - On/Off
6d438146 174 //
175 fInit(kFALSE), // initialization flag - set to kTRUE if corresponding formulas found
176 fFormulaX(0), // x formula - pointer to the function
177 fFormulaY(0), // y formula - pointer to the function
178 fFormulaZ(0) // z formula - pointer to the function
179{
180 if (trafo.fNameX) fNameX = new TString(*(trafo.fNameX));
181 if (trafo.fNameY) fNameY = new TString(*(trafo.fNameY));
182 if (trafo.fNameZ) fNameZ = new TString(*(trafo.fNameZ));
183 if (trafo.fBitMask) fBitMask = new TBits(*(trafo.fBitMask));
184}
185
186AliTPCTransformation::~AliTPCTransformation(){
187 //
188 // destructor
189 //
190 delete fNameX;
191 delete fNameY;
192 delete fNameZ;
193 delete fBitMask;
194 delete fFixedParam;
195}
196void AliTPCTransformation::SetParams(Double_t param, Double_t sigma, Double_t sigma2Time, TVectorD* fixedParams){
197 //
198 // Set parameters of transformation
199 //
200 fParam = param;
201 fSigma = sigma;
202 fSigma2Time = sigma2Time;
203 if (fFixedParam) delete fFixedParam;
204 fFixedParam = new TVectorD(*fixedParams);
205 Init();
206}
94685752 207
208
209Bool_t AliTPCTransformation::Init(){
210 //
211 // associate formulas with pointer to the function
212 //
213 Bool_t isOK=kTRUE;
214 if (fNameX) {
215 fFormulaX=FindFormula(fNameX->Data());
216 if (fFormulaX==0) isOK=kFALSE;
217 }
218 if (fNameY) {
219 fFormulaY=FindFormula(fNameY->Data());
220 if (fFormulaY==0) isOK=kFALSE;
221 }
222 if (fNameZ) {
223 fFormulaZ=FindFormula(fNameZ->Data());
224 if (!fFormulaZ) isOK=kFALSE;
225 }
226 return isOK;
227}
228
0e9efcbe 229
230
94685752 231TBits * AliTPCTransformation::BitsSide(Bool_t aside){
232 //
233 TBits * bits = new TBits(72);
234 for (Int_t i=0; i<72;i++){
235 if (i%36<18 && aside) (*bits)[i]=kTRUE;
236 if (i%36<18 && (!aside)) (*bits)[i]=kFALSE;
237 if (i%36>=18 && aside) (*bits)[i]=kFALSE;
238 if (i%36>=18 && (!aside)) (*bits)[i]=kTRUE;
239 }
240 return bits;
241}
242
243TBits * AliTPCTransformation::BitsAll(){
244 //
245 //
246 //
247 TBits * bits = new TBits(72);
248 for (Int_t i=0; i<72;i++){
249 (*bits)[i]=kTRUE;
250 }
251 return bits;
252}
253
254Double_t AliTPCTransformation::GetDeltaXYZ(Int_t coord, Int_t volID, Double_t param, Double_t x, Double_t y, Double_t z){
255 //
256 //
257 //
0e9efcbe 258 if (!fIsActive) return 0;
94685752 259 if (fBitMask && (!(*fBitMask)[volID])) return 0;
0e9efcbe 260 Double_t xyz[5]={x,y,z, param,volID};
94685752 261 if (fCoordSystem==0){
262 // cartezian system
263 if (coord==0 && fFormulaX) return fFormulaX(xyz,fFixedParam->GetMatrixArray());
264 if (coord==1 && fFormulaY) return fFormulaY(xyz,fFixedParam->GetMatrixArray());
0e9efcbe 265 if (coord==2 && fFormulaZ) return fFormulaZ(xyz,fFixedParam->GetMatrixArray());
94685752 266 }
267 if (fCoordSystem==1){
268 // cylindrical system
269 if (coord==2) {
270 if (fFormulaZ==0) return 0;
271 return fFormulaZ(xyz,fFixedParam->GetMatrixArray());
272 }
273 Double_t rrphiz[3]={0,0,0};
274 if (fFormulaX) rrphiz[0] = fFormulaX(xyz,fFixedParam->GetMatrixArray());
275 if (fFormulaY) rrphiz[1] = fFormulaY(xyz,fFixedParam->GetMatrixArray());
276 Double_t alpha = TMath::ATan2(y,x);
277 Double_t ca = TMath::Cos(alpha);
278 Double_t sa = TMath::Sin(alpha);
279 if (coord==0) return ca*rrphiz[0]-sa*rrphiz[1];
280 if (coord==1) return sa*rrphiz[0]+ca*rrphiz[1];
281 }
282 return 0;
283}
284
285Double_t AliTPCTransformation::TPCscalingRPol(Double_t *xyz, Double_t * param){
286 //
287 // Scaling and shift of TPC radius
288 // xyz[0..2] - global xyz of point
289 // xyz[3] - scale parameter
6d438146 290 // param[0] - radial scaling power
291 // param[1] - drift scaling power
292 // radius from -1(at rInner) to 1 (rOuter)
293 // driftM from -1(at 0 drift) to 1 (250 cm drift)
294
295 Double_t rInner=78.8;
296 Double_t rOuter=258.0;
297 Double_t deltaR = rOuter-rInner;
298 Double_t radius = (TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1])-rInner)*2./deltaR;
299 Double_t driftM = (0.5 - TMath::Abs(xyz[2]/250.))*2.0;
300 Double_t delta = TMath::Power(radius,param[0])*TMath::Power(driftM,param[1]);
301 return delta*xyz[3];
94685752 302}
303
304
305Double_t AliTPCTransformation::TPCscalingZDr(Double_t *xyz, Double_t * param){
306 //
307 //
308 // Scaling and shift of TPC radius
309 // xyz[0..2] - global xyz of point
310 // xyz[3] - scale parameter
311 Double_t driftP = TMath::Power(1. - TMath::Abs(xyz[2]/250.), param[0]);
312 Double_t deltaZ = (xyz[2]>0) ? -driftP : driftP;
313 return deltaZ*xyz[3];
314}
315
316Double_t AliTPCTransformation::TPCscalingPhiLocal(Double_t *xyz, Double_t * param){
317 //
318 //
319 // Scaling if the local y -phi
320 // xyz[0..2] - global xyz of point
321 // xyz[3] - scale parameter
6d438146 322 // value = 1 for ful drift length and parameter 1
94685752 323 Double_t alpha = TMath::ATan2(xyz[1],xyz[0]);
324 Double_t sector = TMath::Nint(9*alpha/TMath::Pi()-0.5);
325 Double_t localAlpha = (alpha-(sector+0.5)*TMath::Pi()/9.);
326 Double_t radius = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1])/250.;
6d438146 327 //
328 Double_t deltaAlpha = radius*TMath::Power(2.*9.*localAlpha/TMath::Pi(),param[0]);
94685752 329 return deltaAlpha*xyz[3];
330}
331
332
6d438146 333Double_t AliTPCTransformation::TPCscalingRIFC(Double_t *xyz, Double_t * param){
334 //
335 // inner field cage r distorion - proportinal to 1 over distance to the IFC
336 // param[0] - drift polynom order
337 // distortion at first pad row - is normalized to
338 Double_t rInner=78.8;
339 Double_t rFirst=85.2;
340 Double_t deltaR = rFirst-rInner;
341 Double_t ndistR = (TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1])-rInner)/deltaR;
342 Double_t driftM = (0.5 - TMath::Abs(xyz[2]/250.))*2.;
343 Double_t value = TMath::Power(driftM,param[0])/ndistR;
344 return xyz[3]*value;
345}
346
347Double_t AliTPCTransformation::TPCscalingROFC(Double_t *xyz, Double_t * param){
348 //
349 // outer field cage r distorion - proportinal to 1 over distance to the OFC
350 // param[0] - drift polynom order
351 // driftM - from -1 to 1
352 //
353 Double_t rLast=245.8;
354 Double_t rOuter=258.0;
355 Double_t deltaR = rOuter-rLast;
356 Double_t ndistR = (rOuter-TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]))/deltaR;
357 Double_t driftM = (0.5 - TMath::Abs(xyz[2]/250.))*2.;
358 Double_t value = TMath::Power(driftM,param[0])/ndistR;
359 return xyz[3]*value;
360}
0e9efcbe 361
362
363//
364// TPC sector local misalignment
365//
366Double_t AliTPCTransformation:: TPClocaldLxdGX(Double_t *xyz, Double_t * /*param*/){
367 //
368 // xyz - [0..2] - position
369 // [3] - scale parameter
370 // [4] - volID
371 // return delta in global coordiante system
372 //
373 Int_t sector = TMath::Nint(xyz[4]);
374 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
375 Double_t ca = TMath::Cos(alpha);
376 // Double_t sa = TMath::Sin(alpha);
377 return ca*xyz[3];
378}
379
380Double_t AliTPCTransformation::TPClocaldLxdGY(Double_t *xyz, Double_t * /*param*/){
381 //
382 // xyz - [0..2] - position
383 // [3] - scale parameter
384 // [4] - volID
385 // return delta in global coordiante system
386 //
387 Int_t sector = TMath::Nint(xyz[4]);
388 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
389 //Double_t ca = TMath::Cos(alpha);
390 Double_t sa = TMath::Sin(alpha);
391 return sa*xyz[3];
392}
393
394Double_t AliTPCTransformation:: TPClocaldLydGX(Double_t *xyz, Double_t * /*param*/){
395 //
396 // xyz - [0..2] - position
397 // [3] - scale parameter
398 // [4] - volID
399 // return delta in global coordiante system
400 //
401 Int_t sector = TMath::Nint(xyz[4]);
402 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
403 //Double_t ca = TMath::Cos(alpha);
404 Double_t sa = TMath::Sin(alpha);
405 return -sa*xyz[3];
406}
407
408Double_t AliTPCTransformation::TPClocaldLydGY(Double_t *xyz, Double_t * /*param*/){
409 //
410 // xyz - [0..2] - position
411 // [3] - scale parameter
412 // [4] - volID
413 // return delta in global coordiante system
414 //
415 Int_t sector = TMath::Nint(xyz[4]);
416 Double_t alpha = TMath::Pi()*(sector+0.5)/9;
417 Double_t ca = TMath::Cos(alpha);
418 //Double_t sa = TMath::Sin(alpha);
419 return ca*xyz[3];
420}
421
422
423Double_t AliTPCTransformation::TPCDeltaZ(Double_t *xyz, Double_t * /*param*/){
424 //
425 // xyz - [0..2] - position
426 // [3] - scale parameter
427 // [4] - volID
428 // return delta in global coordiante system
429 //
430 Int_t sector = TMath::Nint(xyz[4]);
431 Double_t signZ = (sector%36<18) ? 1: -1; // drift direction
432 return signZ*xyz[3]; // IROC shift
433}
434
435Double_t AliTPCTransformation::TPCDeltaZMediumLong(Double_t *xyz, Double_t * /*param*/){
436 //
437 // xyz - [0..2] - position
438 // [3] - scale parameter
439 // [4] - volID
440 // return delta in global coordinate system
441 //
442 Int_t sector = TMath::Nint(xyz[4]);
443 Double_t signZ = (sector%36<18) ? 1: -1; // drift direction
444 if (sector<36) return 0;
445 //
446 Double_t radius = (TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]));
447 const Double_t radiusLong = 198.1;
448 Double_t sign = (radius<radiusLong) ? 1:-1;
449 return xyz[3]*sign*signZ;
450}
451
452Double_t AliTPCTransformation::TPCTiltingZ(Double_t *xyz, Double_t * param){
453 // xyz - [0..2] - position
454 // [3] - scale parameter
455 // [4] - volID
456 // param[0] - n for cos
457 // param[1] - n for sin
458 // return delta in global coordinate system
459 Double_t radius = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
460 const Double_t rFirst=85.2;
461 const Double_t rLast =245.8;
462 Double_t radiusC = (rFirst+rLast)*0.5;
463 Double_t deltaR = 2.0*(radius-radiusC)/(rLast-rFirst);
464 Double_t alpha = TMath::ATan2(xyz[1],xyz[0]);
465
466 if (param[0]>0) deltaR *= TMath::Cos(param[0]*alpha);
467 if (param[1]>0) deltaR *= TMath::Sin(param[1]*alpha);
468 return deltaR*xyz[3];
469}
470