AliRoot Core  edcc906 (edcc906)
AliTPCPointCorrection.cxx
Go to the documentation of this file.
1 
3 /**************************************************************************
4  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
5  * *
6  * Author: The ALICE Off-line Project. *
7  * Contributors are mentioned in the code where appropriate. *
8  * *
9  * Permission to use, copy, modify and distribute this software and its *
10  * documentation strictly for non-commercial purposes is hereby granted *
11  * without fee, provided that the above copyright notice appears in all *
12  * copies and that both the copyright notice and this permission notice *
13  * appear in the supporting documentation. The authors make no claims *
14  * about the suitability of this software for any purpose. It is *
15  * provided "as is" without express or implied warranty. *
16  **************************************************************************/
17 
18 /*
19 
20  Unlinearities fitting:
21 
22  Model for Outer field cage:
23  Unlinearities at the edge aproximated using two exponential decays.
24 
25  dz = dz0(r,z) +dr(r,z)*tan(theta)
26  dy = +dr(r,z)*tan(phi)
27 
28 
29 
30 
31 */
32 
33 #include "AliTPCcalibDB.h"
34 #include "TLinearFitter.h"
35 #include "Riostream.h"
36 #include "TList.h"
37 #include "TMath.h"
38 #include "TCanvas.h"
39 #include "TFile.h"
40 #include "TF1.h"
41 #include "TVectorD.h"
42 #include "AliLog.h"
43 #include "AliTPCROC.h"
44 #include "AliTPCClusterParam.h"
45 #include "AliTPCPointCorrection.h"
46 
48 
50 ClassImp(AliTPCPointCorrection)
52 
54  TNamed(),
55  fParamsOutR(38),
56  fParamsOutZ(38),
57  fParamOutRVersion(0),
58  fErrorsOutR(38),
59  fErrorsOutZ(38),
60  fParamOutZVersion(0),
61  //
62  fXIO(0),
63  fXmiddle(0),
64  fXquadrant(0),
65  fArraySectorIntParam(36), // array of sector alignment parameters
66  fArraySectorIntCovar(36), // array of sector alignment covariances
67  //
68  // Kalman filter for global alignment
69  //
70  fSectorParam(0), // Kalman parameter
71  fSectorCovar(0) // Kalman covariance
72 
73 {
74  //
75  // Default constructor
76  //
78  fXquadrant = roc->GetPadRowRadii(36,53);
79  fXmiddle = ( roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5;
80  fXIO = ( roc->GetPadRowRadii(0,roc->GetNRows(0)-1)+roc->GetPadRowRadii(36,0))*0.5;
81 }
82 
83 AliTPCPointCorrection::AliTPCPointCorrection(const Text_t *name, const Text_t *title):
84  TNamed(name,title),
85  fParamsOutR(38),
86  fParamsOutZ(38),
87  fParamOutRVersion(0),
88  fErrorsOutR(38),
89  fErrorsOutZ(38),
90  fParamOutZVersion(0),
91  //
92  //
93  //
94  fXIO(0),
95  fXmiddle(0),
96  fXquadrant(0),
97  fArraySectorIntParam(36), // array of sector alignment parameters
98  fArraySectorIntCovar(36), // array of sector alignment covariances
99  //
100  // Kalman filter for global alignment
101  //
102  fSectorParam(0), // Kalman parameter for A side
103  fSectorCovar(0) // Kalman covariance for A side
104 
105 {
107 
108  AliTPCROC * roc = AliTPCROC::Instance();
109  fXquadrant = roc->GetPadRowRadii(36,53);
110  fXmiddle = ( roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5;
111  fXIO = ( roc->GetPadRowRadii(0,roc->GetNRows(0)-1)+roc->GetPadRowRadii(36,0))*0.5;
112 
113 }
114 
117 
118 }
119 
120 
122 {
125 
126  if (fgInstance == 0){
128  }
129  return fgInstance;
130 }
131 
132 
133 
134 Double_t AliTPCPointCorrection::GetDrOut(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector){
136 
137  if (fParamOutRVersion==0) return CorrectionOutR0(isGlobal, type,cx,cy,cz,sector);
138  return 0;
139 }
140 
141 Double_t AliTPCPointCorrection::SGetDrOut(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector){
143 
144  return fgInstance->GetDrOut(isGlobal, type,cx,cy,cz,sector);
145 }
146 
147 
148 
149 
150 Double_t AliTPCPointCorrection::GetDzOut(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector){
152 
153  if (fParamOutZVersion==0) return CorrectionOutZ0(isGlobal, type,cx,cy,cz,sector);
154  return 0;
155 }
156 
157 
158 Double_t AliTPCPointCorrection::SGetDzOut(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector){
160 
161  return fgInstance->GetDzOut(isGlobal, type,cx,cy,cz,sector);
162 }
163 
164 
165 
166 
167 Double_t AliTPCPointCorrection::CorrectionOutR0(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector){
174 
175  if (isGlobal){
176  // recalculate sector if in global frame
177  Double_t alpha = TMath::ATan2(cy,cx);
178  if (alpha<0) alpha+=TMath::Pi()*2;
179  sector = Int_t(18*alpha/(TMath::Pi()*2));
180  }
181 
182  if (type==kFALSE) sector=36+(sector%36>=18); // side level parameters
183  // dout
184  Double_t radius = TMath::Sqrt(cx*cx+cy*cy);
185  AliTPCROC * roc = AliTPCROC::Instance();
186  Double_t xout = roc->GetPadRowRadiiUp(roc->GetNRows(37)-1);
187  Double_t dout = xout-radius;
188  if (dout<0) return 0;
189  //drift
190  Double_t dr = 0.5 - TMath::Abs(cz/250.);
191  //
192  //
193  TVectorD * vec = GetParamOutR(sector);
194  if (!vec) return 0;
195  Double_t eout10 = TMath::Exp(-dout/10.);
196  Double_t eout5 = TMath::Exp(-dout/5.);
197  Double_t eoutp = (eout10+eout5)*0.5;
198  Double_t eoutm = (eout10-eout5)*0.5;
199  //
200  Double_t result=0;
201  result+=(*vec)[1]*eoutp;
202  result+=(*vec)[2]*eoutm;
203  result+=(*vec)[3]*eoutp*dr;
204  result+=(*vec)[4]*eoutm*dr;
205  result+=(*vec)[5]*eoutp*dr*dr;
206  result+=(*vec)[6]*eoutm*dr*dr;
207  return result;
208 }
209 
210 Double_t AliTPCPointCorrection::RPhiCOGCorrection(Int_t sector, Int_t padrow, Float_t pad, Float_t cy, Float_t y, Float_t z, Float_t ky,Float_t qMax, Float_t threshold){
227 
229  if (!clparam) {
230  AliFatal("TPC OCDB not initialized");
231  return 0;
232  }
233  Int_t padtype=0;
234  if (sector>=36) padtype = (padrow>64)?2:1;
235  Double_t padwidth=(padtype==0)? 0.4:0.6;
236 
237  Double_t sigma = clparam->GetRMS0(0,padtype,250-TMath::Abs(z),ky)/padwidth;
238  //
239  Int_t npads = AliTPCROC::Instance()->GetNPads(sector,padrow);
240  Float_t yshift = TMath::Abs(cy)-TMath::Abs((-npads*0.5+pad)*padwidth); // the clusters can be corrected before
241  // shift to undo correction
242 
243  y -= yshift*((y>0)?1.:-1.); // y in the sector coordinate
244  Double_t y0 = npads*0.5*padwidth;
245  Double_t dy = (TMath::Abs(y0)-TMath::Abs(y))/padwidth-0.5; // distance to the center of pad0
246  //
247  Double_t padcenter = TMath::Nint(dy);
248  Double_t sumw=0;
249  Double_t sumwi=0;
250  for (Double_t ip=padcenter-2.5; ip<=padcenter+2.5;ip++){
251  Double_t weightGaus = qMax*TMath::Exp(-(dy-ip)*(dy-ip)/(2*sigma*sigma));
252  Double_t ypad = (ip-npads*0.5)*padwidth;
253  Double_t weightGain = GetEdgeQ0(sector,padrow,ypad);
254  //
255  Double_t weight = TMath::Nint(weightGaus*weightGain);
256  if (ip<0 &&weight> threshold) weight = threshold; // this is done in cl finder
257  if (weight<0) continue;
258  sumw+=weight;
259  sumwi+=weight*(ip);
260  }
261  Double_t result =0;
262  if (sumw>0) result = sumwi/sumw;
263  result = (result-dy)*padwidth;
264  return result;
265 }
266 
267 
268 
269 
270 Double_t AliTPCPointCorrection::CorrectionOutZ0(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector){
277 
278  if (isGlobal){
279  // recalculate sector if in global frame
280  Double_t alpha = TMath::ATan2(cy,cx);
281  if (alpha<0) alpha+=TMath::Pi()*2;
282  sector = Int_t(18*alpha/(TMath::Pi()*2));
283  }
284 
285  if (type==kFALSE) sector=36+(sector%36>=18); // side level parameters
286  // dout
287  Double_t radius = TMath::Sqrt(cx*cx+cy*cy);
288  AliTPCROC * roc = AliTPCROC::Instance();
289  Double_t xout = roc->GetPadRowRadiiUp(roc->GetNRows(37)-1);
290  Double_t dout = xout-radius;
291  if (dout<0) return 0;
292  //drift
293  Double_t dr = 0.5 - TMath::Abs(cz/250.);
294  //
295  //
296  TVectorD * vec = GetParamOutR(sector);
297  if (!vec) return 0;
298  Double_t eout10 = TMath::Exp(-dout/10.);
299  Double_t eout5 = TMath::Exp(-dout/5.);
300  Double_t eoutp = (eout10+eout5)*0.5;
301  Double_t eoutm = (eout10-eout5)*0.5;
302  //
303  Double_t result=0;
304  result+=(*vec)[1]*eoutp;
305  result+=(*vec)[2]*eoutm;
306  result+=(*vec)[3]*eoutp*dr;
307  result+=(*vec)[4]*eoutm*dr;
308  result+=(*vec)[5]*eoutp*dr*dr;
309  result+=(*vec)[6]*eoutm*dr*dr;
310  return result;
311 
312 }
313 
314 Double_t AliTPCPointCorrection::GetEdgeQ0(Int_t sector, Int_t padrow, Float_t y){
316 
317  /* TF1 fexp("fexp","1-exp(-[0]*(x-[1]))",0,20)
318  | param [0] | param [1]
319  IROC | 4.71413 | 1.39558
320  OROC1| 2.11437 | 1.52643
321  OROC2| 2.15082 | 1.53537
322  */
323 
324  Double_t params[2]={100,0};
325  if (sector<36){
326  params[0]=4.71413; params[1]=1.39558;
327  }
328  if (sector>=36){
329  if (padrow<64) { params[0]=2.114; params[1]=1.526;}
330  if (padrow>=64){ params[0]=2.15; params[1]=1.535;}
331  }
332  Double_t result= 1;
333  Double_t xrow = AliTPCROC::Instance()->GetPadRowRadii(sector,padrow);
334  Double_t ymax = TMath::Tan(TMath::Pi()/18.)*xrow;
335  Double_t dedge = ymax-TMath::Abs(y);
336  if (dedge-params[1]<0) return 0;
337  if (dedge>10.*params[1]) return 1;
338  result = 1.-TMath::Exp(-params[0]*(dedge-params[1]));
339  return result;
340 }
341 
342 Double_t AliTPCPointCorrection::SRPhiCOGCorrection(Int_t sector, Int_t padrow, Float_t pad, Float_t cy, Float_t y, Float_t z, Float_t ky,Float_t qMax, Float_t threshold){
343  return fgInstance->RPhiCOGCorrection(sector, padrow, pad, cy, y,z, ky, qMax, threshold);
344 }
345 
346 Double_t AliTPCPointCorrection::SGetEdgeQ0(Int_t sector, Int_t padrow, Float_t y){
348 
349  return fgInstance->GetEdgeQ0(sector, padrow, y);
350 }
351 
352 void AliTPCPointCorrection::AddCorrectionSector(TObjArray & sideAPar, TObjArray &sideCPar, TObjArray & sideACov, TObjArray &sideCCov, Bool_t reset){
354 
355  for (Int_t isec=0;isec<36;isec++){
356  TMatrixD * mat1 = (TMatrixD*)(sideAPar.At(isec));
357  TMatrixD * mat1Covar = (TMatrixD*)(sideACov.At(isec));
358  if (!mat1) continue;
359  TMatrixD * mat0 = (TMatrixD*)(fArraySectorIntParam.At(isec));
360  TMatrixD * mat0Covar = (TMatrixD*)(fArraySectorIntCovar.At(isec));
361  if (!mat0) {
362  fArraySectorIntParam.AddAt(mat1->Clone(),isec);
363  fArraySectorIntCovar.AddAt(mat1Covar->Clone(),isec);
364  continue;
365  }
366  (*mat0Covar)=(*mat1Covar);
367  if (reset) (*mat0)=(*mat1);
368  if (!reset) (*mat0)+=(*mat1);
369  }
370 
371  for (Int_t isec=0;isec<36;isec++){
372  TMatrixD * mat1 = (TMatrixD*)(sideCPar.At(isec));
373  TMatrixD * mat1Covar = (TMatrixD*)(sideCCov.At(isec));
374  if (!mat1) continue;
375  TMatrixD * mat0 = (TMatrixD*)(fArraySectorIntParam.At(isec));
376  TMatrixD * mat0Covar = (TMatrixD*)(fArraySectorIntCovar.At(isec));
377  if (!mat0) {
378  fArraySectorIntParam.AddAt(mat1->Clone(),isec);
379  fArraySectorIntCovar.AddAt(mat1Covar->Clone(),isec);
380  continue;
381  }
382  (*mat0Covar)=(*mat1Covar);
383  if (reset) (*mat0)=(*mat1);
384  if (!reset) (*mat0)+=(*mat1);
385  }
386 }
387 
388 
389 Double_t AliTPCPointCorrection::GetCorrectionSector(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t /*lz*/, Int_t quadrant){
391 
392  TMatrixD * param = (TMatrixD*)fArraySectorIntParam.At(sector%36);
393  if (!param) return 0;
394  if (quadrant<0){ //recalc quadrant if not specified
395  if (lx<fXIO) quadrant=0;
396  if(lx>fXIO){
397  if (lx<fXquadrant) {
398  if (ly<0) quadrant=1;
399  if (ly>0) quadrant=2;
400  }
401  if (lx>=fXquadrant) {
402  if (ly<0) quadrant=3;
403  if (ly>0) quadrant=4;
404  }
405  }
406  }
407  Double_t a10 = (*param)(quadrant*6+0,0);
408  Double_t a20 = (*param)(quadrant*6+1,0);
409  Double_t a21 = (*param)(quadrant*6+2,0);
410  Double_t dx = (*param)(quadrant*6+3,0);
411  Double_t dy = (*param)(quadrant*6+4,0);
412  Double_t dz = (*param)(quadrant*6+5,0);
413  Double_t deltaX = lx-fXmiddle;
414  if (coord==0) return dx;
415  if (coord==1) return dy+deltaX*a10;
416  if (coord==2) return dz+deltaX*a20+ly*a21;
417  if (coord==3) return a10;
418  if (coord==4) return a20;
419  if (coord==5) return a21;
420  //
421  return 0;
422 }
423 
424 Double_t AliTPCPointCorrection::SGetCorrectionSector(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz, Int_t quadrant){
426 
427  if (!Instance()) return 0;
428  return Instance()->GetCorrectionSector(coord,sector,lx,ly,lz, quadrant);
429 }
430 
431 
432 
433 Double_t AliTPCPointCorrection::GetCorrection(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t /*lz*/){
435 
436  if (!fSectorParam) return 0;
437 
438  Double_t a10 = (*fSectorParam)(sector*6+0,0);
439  Double_t a20 = (*fSectorParam)(sector*6+1,0);
440  Double_t a21 = (*fSectorParam)(sector*6+2,0);
441  Double_t dx = (*fSectorParam)(sector*6+3,0);
442  Double_t dy = (*fSectorParam)(sector*6+4,0);
443  Double_t dz = (*fSectorParam)(sector*6+5,0);
444  Double_t deltaX = lx-fXIO;
445  //
446  if (coord==0) return dx;
447  if (coord==1) return dy+deltaX*a10;
448  if (coord==2) return dz+deltaX*a20+ly*a21;
449  if (coord==3) return a10;
450  if (coord==4) return a20;
451  if (coord==5) return a21;
452  return 0;
453 }
454 
455 Double_t AliTPCPointCorrection::SGetCorrection(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz){
457 
458  if (!Instance()) return 0;
459  return Instance()->GetCorrection(coord,sector,lx,ly,lz);
460 }
461 
462 
463 
464 
465 
466 
467 
static AliTPCcalibDB * Instance()
Double_t GetCorrection(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz)
Double_t CorrectionOutZ0(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector)
UInt_t GetNPads(UInt_t sector, UInt_t row) const
Definition: AliTPCROC.h:30
Double_t SRPhiCOGCorrection(Int_t sector, Int_t padrow, Float_t pad, Float_t cy, Float_t y, Float_t z, Float_t ky, Float_t qMax, Float_t threhsold)
void AddCorrectionSector(TObjArray &sideAPar, TObjArray &sideCPar, TObjArray &sideACov, TObjArray &sideCCov, Bool_t reset)
AliTPCClusterParam * GetClusterParam() const
TVectorD * GetParamOutR(Int_t sector)
TMatrixD * fSectorParam
Kalman parameter.
#define TObjArray
TVectorD vec
Definition: AnalyzeLaser.C:8
static AliTPCPointCorrection * fgInstance
UInt_t GetNRows(UInt_t sector) const
Definition: AliTPCROC.h:28
Double_t CorrectionOutR0(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector)
Double_t fXIO
OROC-IROC boundary.
static Double_t SGetDzOut(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector)
static Double_t SGetCorrection(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz)
Int_t fParamOutRVersion
version of the parameterization
Double_t GetDrOut(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector)
static Double_t SGetCorrectionSector(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz, Int_t quadrant=-1)
static Double_t SGetDrOut(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector)
TPC cluster error, shape and charge parameterization as function of drift length and inclination angl...
static AliTPCPointCorrection * Instance()
Float_t GetPadRowRadiiUp(UInt_t irow) const
Definition: AliTPCROC.h:55
Int_t fParamOutZVersion
version of the parameterization
Geometry class for a single ROC.
Definition: AliTPCROC.h:14
Double_t GetEdgeQ0(Int_t sector, Int_t padrow, Float_t y)
Float_t GetRMS0(Int_t dim, Int_t type, Float_t z, Float_t angle) const
#define AliFatal(message)
Definition: AliLog.h:640
TObjArray fArraySectorIntParam
array of sector alignment parameters
Double_t GetDzOut(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector)
Double_t RPhiCOGCorrection(Int_t sector, Int_t padrow, Float_t pad, Float_t cy, Float_t y, Float_t z, Float_t ky, Float_t qMax, Float_t threhsold)
static Double_t SGetEdgeQ0(Int_t sector, Int_t padrow, Float_t y)
static AliTPCROC * Instance()
Definition: AliTPCROC.cxx:34
class TVectorT< Double_t > TVectorD
Double_t fXquadrant
x quadrant
Float_t GetPadRowRadii(UInt_t isec, UInt_t irow) const
Definition: AliTPCROC.h:56
Double_t GetCorrectionSector(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz, Int_t quadrant=-1)
TObjArray fArraySectorIntCovar
array of sector alignment covariances
class TMatrixT< Double_t > TMatrixD
Double_t fXmiddle
center of the TPC sector local X