AliRoot Core  edcc906 (edcc906)
AliTPCROCVoltError3D.cxx
Go to the documentation of this file.
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 
52 
53 
54 #include "AliMagF.h"
55 #include "TGeoGlobalMagField.h"
56 #include "AliTPCcalibDB.h"
57 #include "AliTPCParam.h"
58 #include "AliLog.h"
59 #include "TMatrixD.h"
60 #include "TFile.h"
61 
62 #include "TMath.h"
63 #include "AliTPCROC.h"
64 #include "AliTPCROCVoltError3D.h"
65 
67 ClassImp(AliTPCROCVoltError3D)
69 
71  : AliTPCCorrection("ROCVoltErrors","ROC z alignment Errors"),
72  fC0(0.),fC1(0.),
73  fROCdisplacement(kTRUE),
74  fElectronArrivalCorrection(kTRUE),
75  fInitLookUp(kFALSE),
76  fROCDataFileName(""),
77  fdzDataLinFit(0)
78 {
79  //
80  // default constructor
81  //
82 
83  // Array which will contain the solution according to the setted boundary conditions
84  // main input: z alignment of the Read Out chambers
85  // see InitROCVoltError3D() function
86  for ( Int_t k = 0 ; k < kNPhi ; k++ ) {
87  fLookUpErOverEz[k] = new TMatrixF(kNR,kNZ);
88  fLookUpEphiOverEz[k] = new TMatrixF(kNR,kNZ);
89  fLookUpDeltaEz[k] = new TMatrixF(kNR,kNZ);
90  }
91  // fROCDataFileName="$ALICE_ROOT/TPC/Calib/maps/TPCROCdzSurvey.root";
92  // SetROCDataFileName(fROCDataFileName.Data()); // initialization of fdzDataLinFit is included
93 
94 }
95 
98 
99  for ( Int_t k = 0 ; k < kNPhi ; k++ ) {
100  delete fLookUpErOverEz[k];
101  delete fLookUpEphiOverEz[k];
102  delete fLookUpDeltaEz[k];
103  }
104 
105  delete fdzDataLinFit;
106 }
107 
108 
109 
115 
116  if (corr==NULL) {
117  AliError("Zerro pointer - correction");
118  return kFALSE;
119  }
120  AliTPCROCVoltError3D * corrC = dynamic_cast<AliTPCROCVoltError3D *>(corr);
121  if (corrC == NULL) {
122  AliError(TString::Format("Inconsistent class types: %s\%s",IsA()->GetName(),corr->IsA()->GetName()).Data());
123  return kFALSE;
124  }
125  //
126  TMatrixD matrixDz=*(corrC->fdzDataLinFit);
127  matrixDz*=weight;
128  if (fdzDataLinFit==NULL) {
129  fdzDataLinFit=new TMatrixD(matrixDz);
132  }
133  else{
134  (*fdzDataLinFit)+=matrixDz;
135  }
136  return kTRUE;
137 }
138 
139 
140 
141 
142 
146 
147  if (!fdzDataLinFit) fdzDataLinFit=new TMatrixD(*matrix);
148  else *fdzDataLinFit = *matrix;
149 }
150 
151 
154 
155  AliMagF* magF= (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
156  if (!magF) AliError("Magneticd field - not initialized");
157  Double_t bzField = magF->SolenoidField()/10.; //field in T
159  if (!param) AliError("Parameters - not initialized");
160  Double_t vdrift = param->GetDriftV()/1000000.; // [cm/us] // From dataBase: to be updated: per second (ideally)
161  Double_t ezField = 400; // [V/cm] // to be updated: never (hopefully)
162  Double_t wt = -10.0 * (bzField*10) * vdrift / ezField ;
163  // Correction Terms for effective omegaTau; obtained by a laser calibration run
164  SetOmegaTauT1T2(wt,fT1,fT2);
165 
167 }
168 
169 void AliTPCROCVoltError3D::Update(const TTimeStamp &/*timeStamp*/) {
171 
172  AliMagF* magF= (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
173  if (!magF) AliError("Magneticd field - not initialized");
174  Double_t bzField = magF->SolenoidField()/10.; //field in T
176  if (!param) AliError("Parameters - not initialized");
177  Double_t vdrift = param->GetDriftV()/1000000.; // [cm/us] // From dataBase: to be updated: per second (ideally)
178  Double_t ezField = 400; // [V/cm] // to be updated: never (hopefully)
179  Double_t wt = -10.0 * (bzField*10) * vdrift / ezField ;
180  // Correction Terms for effective omegaTau; obtained by a laser calibration run
181  SetOmegaTauT1T2(wt,fT1,fT2);
182 
183 }
184 
187 
189 
190  TFile f(fROCDataFileName.Data(),"READ");
191  TMatrixD *m = (TMatrixD*) f.Get("dzSurveyLinFitData");
192  TMatrixD &mf = *m;
193 
194  // prepare some space
195 
196  if (fdzDataLinFit) delete fdzDataLinFit;
197  fdzDataLinFit = new TMatrixD(72,3);
198  TMatrixD &dataIntern = *fdzDataLinFit;
199 
200  for (Int_t iroc=0;iroc<72;iroc++) {
201  dataIntern(iroc,0) = mf(iroc,0); // z0 offset
202  dataIntern(iroc,1) = mf(iroc,1); // slope in x
203  dataIntern(iroc,2) = mf(iroc,2); // slope in y
204  }
205 
206  f.Close();
207 
208  fInitLookUp = kFALSE;
209 
210 }
211 
212 void AliTPCROCVoltError3D::GetCorrection(const Float_t x[],const Short_t roc,Float_t dx[]) {
214 
215  const Double_t kEpsilon=Double_t(FLT_MIN);
216  if (!fInitLookUp) {
217  AliInfo("Lookup table was not initialized! Perform the inizialisation now ...");
219  }
220  static Bool_t forceInit=kTRUE; //temporary needed for back compatibility with old OCDB entries
221  if (forceInit&&fLookUpErOverEz[0]){
222  if (TMath::Abs(fLookUpErOverEz[0]->Sum())<kEpsilon){//temporary needed for back compatibility with old OCDB entries
224  }
225  forceInit=kFALSE;
226  }
227 
228 
229  Int_t order = 1 ; // FIXME: hardcoded? Linear interpolation = 1, Quadratic = 2
230 
231  Float_t intEr, intEphi, intDeltaEz;
232  Double_t r, phi, z ;
233  Int_t sign;
234 
235  r = TMath::Sqrt( x[0]*x[0] + x[1]*x[1] ) ;
236  phi = TMath::ATan2(x[1],x[0]) ;
237  if ( phi < 0 ) phi += TMath::TwoPi() ; // Table uses phi from 0 to 2*Pi
238  z = x[2] ; // Create temporary copy of x[2]
239 
240  if ( (roc%36) < 18 ) {
241  sign = 1; // (TPC A side)
242  } else {
243  sign = -1; // (TPC C side)
244  }
245 
246  if ( sign==1 && z < fgkZOffSet ) z = fgkZOffSet; // Protect against discontinuity at CE
247  if ( sign==-1 && z > -fgkZOffSet ) z = -fgkZOffSet; // Protect against discontinuity at CE
248 
249 
250  if ( (sign==1 && z<0) || (sign==-1 && z>0) ) // just a consistency check
251  AliError("ROC number does not correspond to z coordinate! Calculation of distortions is most likely wrong!");
252 
253  // Get the Er and Ephi field integrals plus the integral over DeltaEz
254  intEr = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi,
256  intEphi = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi,
258  intDeltaEz = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi,
260 
261  // printf("%lf %lf %lf\n",intEr,intEphi,intDeltaEz);
262 
263  // Calculate distorted position
264  if ( r > 0.0 ) {
265  phi = phi + ( fC0*intEphi - fC1*intEr ) / r;
266  r = r + ( fC0*intEr + fC1*intEphi );
267  }
268 
269  // Calculate correction in cartesian coordinates
270  dx[0] = r * TMath::Cos(phi) - x[0];
271  dx[1] = r * TMath::Sin(phi) - x[1];
272  dx[2] = intDeltaEz; // z distortion - (internally scaled with driftvelocity dependency
273  // on the Ez field plus the actual ROC misalignment (if set TRUE)
274 
275 
277 
278  // correction for the OROC (in average, a 0.014usec longer drift time
279  // due to different position of the anode wires) -> vd*dt -> 2.64*0.014 = 0.0369 cm
280  // FIXME: correction are token from Magboltz simulations
281  // should be loaded from a database
282 
283  AliTPCROC * rocInfo = AliTPCROC::Instance();
284  Double_t rCrossingROC = (rocInfo->GetPadRowRadii(0,62)+rocInfo->GetPadRowRadii(36,0))/2;
285 
286  if (r>rCrossingROC) {
287  if (sign==1)
288  dx[2] += 0.0369; // A side - negative correction
289  else
290  dx[2] -= 0.0369; // C side - positive correction
291  }
292 
293  }
294 
295 }
296 
302 
303  const Int_t order = 1 ; // Linear interpolation = 1, Quadratic = 2
304  const Float_t gridSizeR = (fgkOFCRadius-fgkIFCRadius) / (kRows-1) ;
305  const Float_t gridSizeZ = fgkTPCZ0 / (kColumns-1) ;
306  const Float_t gridSizePhi = TMath::TwoPi() / ( 18.0 * kPhiSlicesPerSector);
307 
308  // temporary arrays to create the boundary conditions
309  TMatrixD *arrayofArrayV[kPhiSlices], *arrayofCharge[kPhiSlices] ;
310  TMatrixD *arrayofEroverEz[kPhiSlices], *arrayofEphioverEz[kPhiSlices], *arrayofDeltaEz[kPhiSlices] ;
311 
312  for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) {
313  arrayofArrayV[k] = new TMatrixD(kRows,kColumns) ;
314  arrayofCharge[k] = new TMatrixD(kRows,kColumns) ;
315  arrayofEroverEz[k] = new TMatrixD(kRows,kColumns) ;
316  arrayofEphioverEz[k] = new TMatrixD(kRows,kColumns) ;
317  arrayofDeltaEz[k] = new TMatrixD(kRows,kColumns) ;
318  }
319 
320  // list of point as used in the poisson relation and the interpolation (during sum up)
321  Double_t rlist[kRows], zedlist[kColumns] , philist[kPhiSlices];
322  for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) {
323  philist[k] = gridSizePhi * k;
324  for ( Int_t i = 0 ; i < kRows ; i++ ) {
325  rlist[i] = fgkIFCRadius + i*gridSizeR ;
326  for ( Int_t j = 0 ; j < kColumns ; j++ ) { // Fill Vmatrix with Boundary Conditions
327  zedlist[j] = j * gridSizeZ ;
328  }
329  }
330  }
331 
332  // ==========================================================================
333  // Solve Poisson's equation in 3D cylindrical coordinates by relaxation technique
334  // Allow for different size grid spacing in R and Z directions
335 
336  const Int_t symmetry = 0;
337 
338  // Set bondaries and solve Poisson's equation --------------------------
339 
340  if ( !fInitLookUp ) {
341 
342  AliInfo(Form("Solving the poisson equation (~ %d sec)",2*10*(int)(kPhiSlices/10)));
343 
344  for ( Int_t side = 0 ; side < 2 ; side++ ) { // Solve Poisson3D twice; once for +Z and once for -Z
345 
346  for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) {
347  TMatrixD &arrayV = *arrayofArrayV[k] ;
348  TMatrixD &charge = *arrayofCharge[k] ;
349 
350  //Fill arrays with initial conditions. V on the boundary and Charge in the volume.
351  for ( Int_t i = 0 ; i < kRows ; i++ ) {
352  for ( Int_t j = 0 ; j < kColumns ; j++ ) { // Fill Vmatrix with Boundary Conditions
353  arrayV(i,j) = 0.0 ;
354  charge(i,j) = 0.0 ;
355 
356  Float_t radius0 = rlist[i] ;
357  Float_t phi0 = gridSizePhi * k ;
358 
359  // To avoid problems at sector boundaries, use an average of +- 1 degree from actual phi location
360  if ( j == (kColumns-1) ) {
361  arrayV(i,j) = 0.5* ( GetROCVoltOffset( side, radius0, phi0+0.02 ) + GetROCVoltOffset( side, radius0, phi0-0.02 ) ) ;
362 
363  if (side==1) // C side
364  arrayV(i,j) = -arrayV(i,j); // minus sign on the C side to allow a consistent usage of global z when setting the boundaries
365  }
366  }
367  }
368 
369  for ( Int_t i = 1 ; i < kRows-1 ; i++ ) {
370  for ( Int_t j = 1 ; j < kColumns-1 ; j++ ) {
371  charge(i,j) = 0.0 ;
372  }
373  }
374  }
375 
376  // Solve Poisson's equation in 3D cylindrical coordinates by relaxation technique
377  // Allow for different size grid spacing in R and Z directions
378 
379  PoissonRelaxation3D( arrayofArrayV, arrayofCharge,
380  arrayofEroverEz, arrayofEphioverEz, arrayofDeltaEz,
381  kRows, kColumns, kPhiSlices, gridSizePhi, kIterations,
382  symmetry, fROCdisplacement) ;
383 
384 
385  //Interpolate results onto a custom grid which is used just for these calculations.
386  Double_t r, phi, z ;
387  for ( Int_t k = 0 ; k < kNPhi ; k++ ) {
388  phi = fgkPhiList[k] ;
389 
390  TMatrixF &erOverEz = *fLookUpErOverEz[k] ;
391  TMatrixF &ephiOverEz = *fLookUpEphiOverEz[k];
392  TMatrixF &deltaEz = *fLookUpDeltaEz[k] ;
393 
394  for ( Int_t j = 0 ; j < kNZ ; j++ ) {
395 
396  z = TMath::Abs(fgkZList[j]) ; // Symmetric solution in Z that depends only on ABS(Z)
397 
398  if ( side == 0 && fgkZList[j] < 0 ) continue; // Skip rest of this loop if on the wrong side
399  if ( side == 1 && fgkZList[j] > 0 ) continue; // Skip rest of this loop if on the wrong side
400 
401  for ( Int_t i = 0 ; i < kNR ; i++ ) {
402  r = fgkRList[i] ;
403 
404  // Interpolate basicLookup tables; once for each rod, then sum the results
405  erOverEz(i,j) = Interpolate3DTable(order, r, z, phi, kRows, kColumns, kPhiSlices,
406  rlist, zedlist, philist, arrayofEroverEz );
407  ephiOverEz(i,j) = Interpolate3DTable(order, r, z, phi, kRows, kColumns, kPhiSlices,
408  rlist, zedlist, philist, arrayofEphioverEz);
409  deltaEz(i,j) = Interpolate3DTable(order, r, z, phi, kRows, kColumns, kPhiSlices,
410  rlist, zedlist, philist, arrayofDeltaEz );
411 
412  if (side == 1) deltaEz(i,j) = - deltaEz(i,j); // negative coordinate system on C side
413 
414  } // end r loop
415  }// end z loop
416  }// end phi loop
417 
418  if ( side == 0 ) AliInfo(" A side done");
419  if ( side == 1 ) AliInfo(" C side done");
420  } // end side loop
421  }
422 
423  // clear the temporary arrays lists
424  for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) {
425  delete arrayofArrayV[k];
426  delete arrayofCharge[k];
427  delete arrayofEroverEz[k];
428  delete arrayofEphioverEz[k];
429  delete arrayofDeltaEz[k];
430  }
431 
432 
433  fInitLookUp = kTRUE;
434 
435 }
436 
437 
438 Float_t AliTPCROCVoltError3D::GetROCVoltOffset(Int_t side, Float_t r0, Float_t phi0) {
441 
442  Float_t xp = r0*TMath::Cos(phi0);
443  Float_t yp = r0*TMath::Sin(phi0);
444 
445  // phi0 should be between 0 and 2pi
446  if (phi0<0) phi0+=TMath::TwoPi();
447  Int_t roc = (Int_t)TMath::Floor((TMath::RadToDeg()*phi0)/20);
448  if (side==1) roc+=18; // C side
449  if (r0>132) roc+=36; // OROC
450 
451  // linear-plane data: z = z0 + kx*lx + ky*ly (rotation in local coordinates)
452  TMatrixD &fitData = *fdzDataLinFit;
453 
454  // local coordinates
455  Double_t secAlpha = TMath::DegToRad()*(10.+20.*(((Int_t)roc)%18));
456  Float_t lx = xp*TMath::Cos(secAlpha)+yp*TMath::Sin(secAlpha);
457  Float_t ly = yp*TMath::Cos(secAlpha)-xp*TMath::Sin(secAlpha);
458 
459  // reference of rotation in lx is at the intersection between OROC and IROC
460  // necessary, since the Fitprozedure is otherwise useless
461 
462  AliTPCROC * rocInfo = AliTPCROC::Instance();
463  Double_t lxRef = (rocInfo->GetPadRowRadii(0,62)+rocInfo->GetPadRowRadii(36,0))/2;
464 
465  Float_t dz = fitData(roc,0)+fitData(roc,1)*(lx-lxRef) + fitData(roc,2)*ly; // value in cm
466 
467  // aproximated Voltage-offset-aquivalent to the z misalignment
468  // (linearly scaled with the z position)
469  Double_t ezField = (fgkCathodeV-fgkGG)/fgkTPCZ0; // = ALICE Electric Field (V/cm) Magnitude ~ -400 V/cm;
470  Float_t voltOff = dz*ezField; // values in "Volt equivalents"
471 
472  return voltOff;
473 }
474 
475 TH2F * AliTPCROCVoltError3D::CreateHistoOfZAlignment(Int_t side, Int_t nx, Int_t ny) {
478 
479  char hname[100];
480  if (side==0) snprintf(hname,100,"survey_dz_Aside");
481  if (side==1) snprintf(hname,100,"survey_dz_Cside");
482 
483  TH2F *h = new TH2F(hname,hname,nx,-250.,250.,ny,-250.,250.);
484 
485  for (Int_t iy=1;iy<=ny;++iy) {
486  Double_t yp = h->GetYaxis()->GetBinCenter(iy);
487  for (Int_t ix=1;ix<=nx;++ix) {
488  Double_t xp = h->GetXaxis()->GetBinCenter(ix);
489 
490  Float_t r0 = TMath::Sqrt(xp*xp+yp*yp);
491  Float_t phi0 = TMath::ATan2(yp,xp);
492 
493  Float_t dz = GetROCVoltOffset(side,r0,phi0); // in [volt]
494 
495  Double_t ezField = (fgkCathodeV-fgkGG)/fgkTPCZ0; // = ALICE Electric Field (V/cm) Magnitude ~ -400 V/cm;
496  dz = dz/ezField; // in [cm]
497 
498  if (85.<=r0 && r0<=245.) {
499  h->SetBinContent(ix,iy,dz);
500  } else {
501  h->SetBinContent(ix,iy,0.);
502  }
503  }
504  }
505 
506  h->GetXaxis()->SetTitle("x [cm]");
507  h->GetYaxis()->SetTitle("y [cm]");
508  h->GetZaxis()->SetTitle("dz [cm]");
509  h->SetStats(0);
510  // h->DrawCopy("colz");
511 
512  return h;
513 }
514 
515 void AliTPCROCVoltError3D::Print(const Option_t* option) const {
518 
519  TString opt = option; opt.ToLower();
520  printf("%s\n",GetTitle());
521  printf(" - z aligmnet of the TPC Read-Out chambers \n");
522  printf(" (linear interpolation within the chamber: dz = z0 + kx*(lx-133) + ky*ly [cm] ) \n");
523  printf(" Info: Check the following data-file for more details: %s \n",fROCDataFileName.Data());
524  printf(" fROCdisplacement: %s\n", fROCdisplacement ? "kTRUE" : "kFALSE");
525  printf(" fElectronArrivalCorrection: %s\n", fElectronArrivalCorrection ? "kTRUE" : "kFALSE");
526 
527  if (opt.Contains("a")) { // Print all details
528  TMatrixD &fitData = *fdzDataLinFit;
529  printf(" A side: IROC ROCX=(z0,kx,ky): \n");
530  for (Int_t roc = 0; roc<18; roc++)
531  printf("ROC%d:(%.2e,%.2e,%.2e) ",roc,fitData(roc,0),fitData(roc,1),fitData(roc,2));
532  printf("\n A side: OROC ROCX=(z0,kx,ky): \n");
533  for (Int_t roc = 36; roc<54; roc++)
534  printf("ROC%d:(%.2e,%.2e,%.2e) ",roc,fitData(roc,0),fitData(roc,1),fitData(roc,2));
535  printf("\n C side: IROC ROCX=(z0,kx,ky): \n");
536  for (Int_t roc = 18; roc<36; roc++)
537  printf("ROC%d:(%.2e,%.2e,%.2e) ",roc,fitData(roc,0),fitData(roc,1),fitData(roc,2));
538  printf("\n C side: OROC ROCX=(z0,kx,ky): \n");
539  for (Int_t roc = 54; roc<72; roc++)
540  printf("ROC%d:(%.2e,%.2e,%.2e) ",roc,fitData(roc,0),fitData(roc,1),fitData(roc,2));
541  printf("\n\n");
542  printf(" - T1: %1.4f, T2: %1.4f \n",fT1,fT2);
543  printf(" - C1: %1.4f, C0: %1.4f \n",fC1,fC0);
544  }
545 
546  if (!fInitLookUp) AliError("Lookup table was not initialized! You should do InitROCVoltError3D() ...");
547 
548 }
static AliTPCcalibDB * Instance()
printf("Chi2/npoints = %f\n", TMath::Sqrt(chi2/npoints))
const Double_t kEpsilon
TMatrixF * fLookUpErOverEz[kNPhi]
Array to store electric field integral (int Er/Ez)
Double_t fgkZList[kNZ]
points in the z direction (for the lookup table)
TString fROCDataFileName
filename of the survey data containing the lin Fit values
TMatrixF * fLookUpDeltaEz[kNPhi]
Array to store electric field integral (int Er/Ez)
Float_t GetROCVoltOffset(Int_t side, Float_t r0, Float_t phi0)
Manager and of geomety classes for set: TPC.
Definition: AliTPCParam.h:18
Float_t GetDriftV() const
Definition: AliTPCParam.h:342
Double_t Interpolate3DTable(Int_t order, Double_t x, Double_t y, Double_t z, Int_t nx, Int_t ny, Int_t nz, const Double_t xv[], const Double_t yv[], const Double_t zv[], TMatrixD **arrayofArrays)
static const Double_t fgkOFCRadius
Mean Radius of the Outer Field Cage (252.55 min, 256.45 max) (cm)
AliTPCParam * GetParameters() const
Double_t SolenoidField() const
Definition: AliMagF.h:67
The class calculates the space point distortions due to z offsets of the.
Double_t fT1
tensor term of wt - T1
static const Double_t fgkZOffSet
Offset from CE: calculate all distortions closer to CE as if at this point.
AliTPCCorrection class.
Bool_t fElectronArrivalCorrection
flag on wheter to consider the difference
virtual void GetCorrection(const Float_t x[], const Short_t roc, Float_t dx[])
virtual void Update(const TTimeStamp &timeStamp)
virtual void Print(const Option_t *option="") const
TMatrixF * fLookUpEphiOverEz[kNPhi]
Array to store electric field integral (int Er/Ez)
static const Double_t fgkTPCZ0
nominal gating grid position
Bool_t fInitLookUp
flag to check it the Look Up table was created (SUM)
Float_t fC1
coefficient C1 (compare Jim Thomas&#39;s notes for definitions)
#define AliInfo(message)
Definition: AliLog.h:484
Geometry class for a single ROC.
Definition: AliTPCROC.h:14
TMatrixD * fdzDataLinFit
Linear fits of dz survey points (each sector=72) (z0,slopeX,slopeY)
virtual void SetOmegaTauT1T2(Float_t omegaTau, Float_t t1, Float_t t2)
void SetROCData(TMatrixD *matrix)
Float_t fC0
coefficient C0 (compare Jim Thomas&#39;s notes for definitions)
TF1 * f
Definition: interpolTest.C:21
void SetROCDataFileName(const char *fname)
static const Double_t fgkCathodeV
Cathode Voltage (volts)
Double_t fgkPhiList[kNPhi]
points in the phi direction (for the lookup table)
static const Double_t fgkGG
Gating Grid voltage (volts)
Double_t fT2
tensor term of wt - T2
Double_t vdrift
Definition: DriftKalman.C:98
Double_t fgkRList[kNR]
points in the radial direction (for the lookup table)
static AliTPCROC * Instance()
Definition: AliTPCROC.cxx:34
#define AliError(message)
Definition: AliLog.h:591
Float_t GetPadRowRadii(UInt_t isec, UInt_t irow) const
Definition: AliTPCROC.h:56
char * fname
static const Double_t fgkIFCRadius
Mean Radius of the Inner Field Cage ( 82.43 min, 83.70 max) (cm)/hera/alice/wiechula/calib/guiTrees.
virtual Bool_t AddCorrectionCompact(AliTPCCorrection *corr, Double_t weight)
Bool_t fROCdisplacement
flag on wheter to consider the ROC displacement
class TMatrixT< Double_t > TMatrixD
TH2F * CreateHistoOfZAlignment(Int_t side, Int_t nx=250, Int_t ny=250)
void PoissonRelaxation3D(TMatrixD **arrayofArrayV, TMatrixD **arrayofChargeDensities, TMatrixD **arrayofEroverEz, TMatrixD **arrayofEPhioverEz, TMatrixD **arrayofEz, Int_t rows, Int_t columns, Int_t phislices, Float_t deltaphi, Int_t iterations, Int_t summetry, Bool_t rocDisplacement=kTRUE, IntegrationType integrationType=kIntegral)