]>
Commit | Line | Data |
---|---|---|
c9cbd2f2 | 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 | ||
16 | ////////////////////////////////////////////////////////////////////////////// | |
17 | // // | |
18 | // AliTPCROCVoltError3D class // | |
19 | // The class calculates the space point distortions due to residual voltage // | |
20 | // errors on Read Out Chambers of the TPC in 3D. // | |
21 | // // | |
22 | // The class allows "effective Omega Tau" corrections. // | |
23 | // // | |
24 | // NOTE: This class is capable of calculating z distortions due to // | |
25 | // misalignment and the vd dependency on the residual drift field // | |
26 | // // | |
27 | // date: 08/08/2010 // | |
28 | // Authors: Jim Thomas, Stefan Rossegger // | |
29 | // // | |
30 | // Example usage : // | |
31 | // AliTPCROCVoltError3D ROCerror; // | |
32 | ////////////////////////////////////////////////////////////////////////////// | |
33 | ||
34 | #include "AliMagF.h" | |
35 | #include "TGeoGlobalMagField.h" | |
36 | #include "AliTPCcalibDB.h" | |
37 | #include "AliTPCParam.h" | |
38 | #include "AliLog.h" | |
39 | #include "TMatrixD.h" | |
40 | #include "TFile.h" | |
41 | ||
42 | #include "TMath.h" | |
43 | #include "AliTPCROC.h" | |
44 | #include "AliTPCROCVoltError3D.h" | |
45 | ||
46 | ClassImp(AliTPCROCVoltError3D) | |
47 | ||
48 | AliTPCROCVoltError3D::AliTPCROCVoltError3D() | |
49 | : AliTPCCorrection("ROCVoltErrors","ROC z alignment Errors"), | |
50 | fC0(0.),fC1(0.), | |
51 | fROCdisplacement(kTRUE), | |
52 | fInitLookUp(kFALSE), | |
2b68ab9c | 53 | fROCDataFileName((char*)"$(ALICE_ROOT)/TPC/Calib/maps/TPCROCdzSurvey.root"), |
c9cbd2f2 | 54 | fdzDataLinFit(0) |
55 | { | |
56 | // | |
57 | // default constructor | |
58 | // | |
59 | ||
60 | // Array which will contain the solution according to the setted boundary conditions | |
61 | // main input: z alignment of the Read Out chambers | |
62 | // see InitROCVoltError3D() function | |
63 | for ( Int_t k = 0 ; k < kNPhi ; k++ ) { | |
64 | fLookUpErOverEz[k] = new TMatrixD(kNR,kNZ); | |
65 | fLookUpEphiOverEz[k] = new TMatrixD(kNR,kNZ); | |
66 | fLookUpDeltaEz[k] = new TMatrixD(kNR,kNZ); | |
67 | } | |
68 | ||
69 | SetROCDataFileName(fROCDataFileName); // initialization of fdzDataLinFit is included | |
70 | ||
71 | } | |
72 | ||
73 | AliTPCROCVoltError3D::~AliTPCROCVoltError3D() { | |
74 | // | |
75 | // destructor | |
76 | // | |
77 | ||
78 | for ( Int_t k = 0 ; k < kNPhi ; k++ ) { | |
79 | delete fLookUpErOverEz[k]; | |
80 | delete fLookUpEphiOverEz[k]; | |
81 | delete fLookUpDeltaEz[k]; | |
82 | } | |
83 | ||
84 | delete fdzDataLinFit; | |
85 | } | |
86 | ||
2bbac918 | 87 | void AliTPCROCVoltError3D::SetDZMap(TMatrixD * matrix){ |
88 | // | |
89 | // | |
90 | // | |
91 | if (!fdzDataLinFit) fdzDataLinFit=new TMatrixD(*matrix); | |
92 | else *fdzDataLinFit = *matrix; | |
93 | } | |
94 | ||
95 | ||
c9cbd2f2 | 96 | void AliTPCROCVoltError3D::Init() { |
97 | // | |
98 | // Initialization funtion | |
99 | // | |
100 | ||
101 | AliMagF* magF= (AliMagF*)TGeoGlobalMagField::Instance()->GetField(); | |
102 | if (!magF) AliError("Magneticd field - not initialized"); | |
103 | Double_t bzField = magF->SolenoidField()/10.; //field in T | |
104 | AliTPCParam *param= AliTPCcalibDB::Instance()->GetParameters(); | |
105 | if (!param) AliError("Parameters - not initialized"); | |
106 | Double_t vdrift = param->GetDriftV()/1000000.; // [cm/us] // From dataBase: to be updated: per second (ideally) | |
107 | Double_t ezField = 400; // [V/cm] // to be updated: never (hopefully) | |
108 | Double_t wt = -10.0 * (bzField*10) * vdrift / ezField ; | |
109 | // Correction Terms for effective omegaTau; obtained by a laser calibration run | |
110 | SetOmegaTauT1T2(wt,fT1,fT2); | |
111 | ||
35ae345f | 112 | if (!fInitLookUp) InitROCVoltError3D(); |
c9cbd2f2 | 113 | } |
114 | ||
115 | void AliTPCROCVoltError3D::Update(const TTimeStamp &/*timeStamp*/) { | |
116 | // | |
117 | // Update function | |
118 | // | |
119 | AliMagF* magF= (AliMagF*)TGeoGlobalMagField::Instance()->GetField(); | |
120 | if (!magF) AliError("Magneticd field - not initialized"); | |
121 | Double_t bzField = magF->SolenoidField()/10.; //field in T | |
122 | AliTPCParam *param= AliTPCcalibDB::Instance()->GetParameters(); | |
123 | if (!param) AliError("Parameters - not initialized"); | |
124 | Double_t vdrift = param->GetDriftV()/1000000.; // [cm/us] // From dataBase: to be updated: per second (ideally) | |
125 | Double_t ezField = 400; // [V/cm] // to be updated: never (hopefully) | |
126 | Double_t wt = -10.0 * (bzField*10) * vdrift / ezField ; | |
127 | // Correction Terms for effective omegaTau; obtained by a laser calibration run | |
128 | SetOmegaTauT1T2(wt,fT1,fT2); | |
129 | ||
130 | } | |
131 | ||
132 | void AliTPCROCVoltError3D::SetROCDataFileName(char *const fname) { | |
133 | // | |
134 | // Set / load the ROC data (linear fit of ROC misalignments) | |
135 | // | |
136 | ||
137 | fROCDataFileName = fname; | |
138 | ||
139 | TFile f(fROCDataFileName,"READ"); | |
140 | TMatrixD *m = (TMatrixD*) f.Get("dzSurveyLinFitData"); | |
141 | TMatrixD &mf = *m; | |
142 | ||
143 | // prepare some space | |
144 | ||
145 | if (fdzDataLinFit) delete fdzDataLinFit; | |
146 | fdzDataLinFit = new TMatrixD(72,3); | |
147 | TMatrixD &dataIntern = *fdzDataLinFit; | |
148 | ||
149 | for (Int_t iroc=0;iroc<72;iroc++) { | |
150 | dataIntern(iroc,0) = mf(iroc,0); // z0 offset | |
151 | dataIntern(iroc,1) = mf(iroc,1); // slope in x | |
152 | dataIntern(iroc,2) = mf(iroc,2); // slope in y | |
153 | } | |
154 | ||
155 | f.Close(); | |
156 | ||
157 | fInitLookUp = kFALSE; | |
158 | ||
159 | } | |
160 | ||
161 | void AliTPCROCVoltError3D::GetCorrection(const Float_t x[],const Short_t roc,Float_t dx[]) { | |
162 | // | |
163 | // Calculates the correction due e.g. residual voltage errors on the TPC boundaries | |
164 | // | |
165 | ||
166 | if (!fInitLookUp) { | |
167 | AliInfo("Lookup table was not initialized! Perform the inizialisation now ..."); | |
168 | InitROCVoltError3D(); | |
169 | return; | |
170 | } | |
171 | ||
172 | Int_t order = 1 ; // FIXME: hardcoded? Linear interpolation = 1, Quadratic = 2 | |
173 | ||
174 | Double_t intEr, intEphi, intDeltaEz; | |
175 | Double_t r, phi, z ; | |
176 | Int_t sign; | |
177 | ||
178 | r = TMath::Sqrt( x[0]*x[0] + x[1]*x[1] ) ; | |
179 | phi = TMath::ATan2(x[1],x[0]) ; | |
180 | if ( phi < 0 ) phi += TMath::TwoPi() ; // Table uses phi from 0 to 2*Pi | |
181 | z = x[2] ; // Create temporary copy of x[2] | |
182 | ||
183 | if ( (roc%36) < 18 ) { | |
184 | sign = 1; // (TPC A side) | |
185 | } else { | |
186 | sign = -1; // (TPC C side) | |
187 | } | |
188 | ||
189 | if ( sign==1 && z < fgkZOffSet ) z = fgkZOffSet; // Protect against discontinuity at CE | |
190 | if ( sign==-1 && z > -fgkZOffSet ) z = -fgkZOffSet; // Protect against discontinuity at CE | |
191 | ||
192 | ||
193 | if ( (sign==1 && z<0) || (sign==-1 && z>0) ) // just a consistency check | |
194 | AliError("ROC number does not correspond to z coordinate! Calculation of distortions is most likely wrong!"); | |
195 | ||
196 | // Get the Er and Ephi field integrals plus the integral over DeltaEz | |
197 | intEr = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, | |
198 | fgkRList, fgkZList, fgkPhiList, fLookUpErOverEz ); | |
199 | intEphi = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, | |
200 | fgkRList, fgkZList, fgkPhiList, fLookUpEphiOverEz); | |
201 | intDeltaEz = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, | |
202 | fgkRList, fgkZList, fgkPhiList, fLookUpDeltaEz ); | |
203 | ||
204 | // printf("%lf %lf %lf\n",intEr,intEphi,intDeltaEz); | |
205 | ||
206 | // Calculate distorted position | |
207 | if ( r > 0.0 ) { | |
208 | phi = phi + ( fC0*intEphi - fC1*intEr ) / r; | |
209 | r = r + ( fC0*intEr + fC1*intEphi ); | |
210 | } | |
211 | ||
212 | // Calculate correction in cartesian coordinates | |
213 | dx[0] = r * TMath::Cos(phi) - x[0]; | |
214 | dx[1] = r * TMath::Sin(phi) - x[1]; | |
215 | dx[2] = intDeltaEz; // z distortion - (internally scaled with driftvelocity dependency | |
216 | // on the Ez field plus the actual ROC misalignment (if set TRUE) | |
217 | ||
218 | } | |
219 | ||
220 | void AliTPCROCVoltError3D::InitROCVoltError3D() { | |
221 | // | |
222 | // Initialization of the Lookup table which contains the solutions of the | |
223 | // Dirichlet boundary problem | |
224 | // Calculation of the single 3D-Poisson solver is done just if needed | |
225 | // (see basic lookup tables in header file) | |
226 | // | |
227 | ||
228 | const Int_t order = 1 ; // Linear interpolation = 1, Quadratic = 2 | |
229 | const Float_t gridSizeR = (fgkOFCRadius-fgkIFCRadius) / (kRows-1) ; | |
230 | const Float_t gridSizeZ = fgkTPCZ0 / (kColumns-1) ; | |
231 | const Float_t gridSizePhi = TMath::TwoPi() / ( 18.0 * kPhiSlicesPerSector); | |
232 | ||
233 | // temporary arrays to create the boundary conditions | |
234 | TMatrixD *arrayofArrayV[kPhiSlices], *arrayofCharge[kPhiSlices] ; | |
235 | TMatrixD *arrayofEroverEz[kPhiSlices], *arrayofEphioverEz[kPhiSlices], *arrayofDeltaEz[kPhiSlices] ; | |
236 | ||
237 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
238 | arrayofArrayV[k] = new TMatrixD(kRows,kColumns) ; | |
239 | arrayofCharge[k] = new TMatrixD(kRows,kColumns) ; | |
240 | arrayofEroverEz[k] = new TMatrixD(kRows,kColumns) ; | |
241 | arrayofEphioverEz[k] = new TMatrixD(kRows,kColumns) ; | |
242 | arrayofDeltaEz[k] = new TMatrixD(kRows,kColumns) ; | |
243 | } | |
244 | ||
245 | // list of point as used in the poisson relation and the interpolation (during sum up) | |
246 | Double_t rlist[kRows], zedlist[kColumns] , philist[kPhiSlices]; | |
247 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
248 | philist[k] = gridSizePhi * k; | |
249 | for ( Int_t i = 0 ; i < kRows ; i++ ) { | |
250 | rlist[i] = fgkIFCRadius + i*gridSizeR ; | |
251 | for ( Int_t j = 0 ; j < kColumns ; j++ ) { // Fill Vmatrix with Boundary Conditions | |
252 | zedlist[j] = j * gridSizeZ ; | |
253 | } | |
254 | } | |
255 | } | |
256 | ||
257 | // ========================================================================== | |
258 | // Solve Poisson's equation in 3D cylindrical coordinates by relaxation technique | |
259 | // Allow for different size grid spacing in R and Z directions | |
260 | ||
261 | const Int_t symmetry = 0; | |
262 | ||
263 | // Set bondaries and solve Poisson's equation -------------------------- | |
264 | ||
265 | if ( !fInitLookUp ) { | |
266 | ||
267 | AliInfo(Form("Solving the poisson equation (~ %d sec)",2*10*(int)(kPhiSlices/10))); | |
268 | ||
269 | for ( Int_t side = 0 ; side < 2 ; side++ ) { // Solve Poisson3D twice; once for +Z and once for -Z | |
270 | ||
271 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
272 | TMatrixD &arrayV = *arrayofArrayV[k] ; | |
273 | TMatrixD &charge = *arrayofCharge[k] ; | |
274 | ||
275 | //Fill arrays with initial conditions. V on the boundary and Charge in the volume. | |
276 | for ( Int_t i = 0 ; i < kRows ; i++ ) { | |
277 | for ( Int_t j = 0 ; j < kColumns ; j++ ) { // Fill Vmatrix with Boundary Conditions | |
278 | arrayV(i,j) = 0.0 ; | |
279 | charge(i,j) = 0.0 ; | |
280 | ||
281 | Float_t radius0 = rlist[i] ; | |
282 | Float_t phi0 = gridSizePhi * k ; | |
283 | ||
284 | // To avoid problems at sector boundaries, use an average of +- 1 degree from actual phi location | |
285 | if ( j == (kColumns-1) ) | |
286 | arrayV(i,j) = 0.5* ( GetROCVoltOffset( side, radius0, phi0+0.02 ) + GetROCVoltOffset( side, radius0, phi0-0.02 ) ) ; | |
287 | ||
288 | } | |
289 | } | |
290 | ||
291 | for ( Int_t i = 1 ; i < kRows-1 ; i++ ) { | |
292 | for ( Int_t j = 1 ; j < kColumns-1 ; j++ ) { | |
293 | charge(i,j) = 0.0 ; | |
294 | } | |
295 | } | |
296 | } | |
297 | ||
298 | // Solve Poisson's equation in 3D cylindrical coordinates by relaxation technique | |
299 | // Allow for different size grid spacing in R and Z directions | |
300 | ||
301 | PoissonRelaxation3D( arrayofArrayV, arrayofCharge, | |
302 | arrayofEroverEz, arrayofEphioverEz, arrayofDeltaEz, | |
303 | kRows, kColumns, kPhiSlices, gridSizePhi, kIterations, | |
304 | symmetry, fROCdisplacement) ; | |
305 | ||
306 | ||
307 | //Interpolate results onto a custom grid which is used just for these calculations. | |
308 | Double_t r, phi, z ; | |
309 | for ( Int_t k = 0 ; k < kNPhi ; k++ ) { | |
310 | phi = fgkPhiList[k] ; | |
311 | ||
312 | TMatrixD &erOverEz = *fLookUpErOverEz[k] ; | |
313 | TMatrixD &ephiOverEz = *fLookUpEphiOverEz[k]; | |
314 | TMatrixD &deltaEz = *fLookUpDeltaEz[k] ; | |
315 | ||
316 | for ( Int_t j = 0 ; j < kNZ ; j++ ) { | |
317 | ||
318 | z = TMath::Abs(fgkZList[j]) ; // Symmetric solution in Z that depends only on ABS(Z) | |
319 | ||
320 | if ( side == 0 && fgkZList[j] < 0 ) continue; // Skip rest of this loop if on the wrong side | |
321 | if ( side == 1 && fgkZList[j] > 0 ) continue; // Skip rest of this loop if on the wrong side | |
322 | ||
323 | for ( Int_t i = 0 ; i < kNR ; i++ ) { | |
324 | r = fgkRList[i] ; | |
325 | ||
326 | // Interpolate basicLookup tables; once for each rod, then sum the results | |
327 | erOverEz(i,j) = Interpolate3DTable(order, r, z, phi, kRows, kColumns, kPhiSlices, | |
328 | rlist, zedlist, philist, arrayofEroverEz ); | |
329 | ephiOverEz(i,j) = Interpolate3DTable(order, r, z, phi, kRows, kColumns, kPhiSlices, | |
330 | rlist, zedlist, philist, arrayofEphioverEz); | |
331 | deltaEz(i,j) = Interpolate3DTable(order, r, z, phi, kRows, kColumns, kPhiSlices, | |
332 | rlist, zedlist, philist, arrayofDeltaEz ); | |
333 | ||
334 | if (side == 1) deltaEz(i,j) = - deltaEz(i,j); // negative coordinate system on C side | |
335 | ||
336 | } // end r loop | |
337 | }// end z loop | |
338 | }// end phi loop | |
339 | ||
340 | if ( side == 0 ) AliInfo(" A side done"); | |
341 | if ( side == 1 ) AliInfo(" C side done"); | |
342 | } // end side loop | |
343 | } | |
344 | ||
345 | // clear the temporary arrays lists | |
346 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
347 | delete arrayofArrayV[k]; | |
348 | delete arrayofCharge[k]; | |
349 | delete arrayofEroverEz[k]; | |
350 | delete arrayofEphioverEz[k]; | |
351 | delete arrayofDeltaEz[k]; | |
352 | } | |
353 | ||
354 | ||
355 | fInitLookUp = kTRUE; | |
356 | ||
357 | } | |
358 | ||
359 | ||
360 | Float_t AliTPCROCVoltError3D::GetROCVoltOffset(Int_t side, Float_t r0, Float_t phi0) { | |
361 | // | |
362 | // Returns the dz alignment data (in voltage equivalents) at | |
363 | // the given position | |
364 | // | |
365 | ||
366 | Float_t xp = r0*TMath::Cos(phi0); | |
367 | Float_t yp = r0*TMath::Sin(phi0); | |
368 | ||
369 | // phi0 should be between 0 and 2pi | |
370 | if (phi0<0) phi0+=TMath::TwoPi(); | |
371 | Int_t roc = (Int_t)TMath::Floor((TMath::RadToDeg()*phi0)/20); | |
372 | if (side==1) roc+=18; // C side | |
373 | if (r0>132) roc+=36; // OROC | |
374 | ||
375 | // linear-plane data: z = z0 + kx*x + ky*y | |
376 | TMatrixD &fitData = *fdzDataLinFit; | |
377 | Float_t dz = fitData(roc,0)+fitData(roc,1)*xp + fitData(roc,2)*yp; // value in cm | |
378 | ||
379 | // aproximated Voltage-offset-aquivalent to the z misalignment | |
380 | // (linearly scaled with the z position) | |
381 | Double_t ezField = (fgkCathodeV-fgkGG)/fgkTPCZ0; // = ALICE Electric Field (V/cm) Magnitude ~ -400 V/cm; | |
382 | Float_t voltOff = dz*ezField; // values in "Volt equivalents" | |
383 | ||
384 | return voltOff; | |
385 | } | |
386 | ||
387 | TH2F * AliTPCROCVoltError3D::CreateHistoOfZSurvey(Int_t side, Int_t nx, Int_t ny) { | |
388 | // | |
389 | // return a simple histogramm containing the input to the poisson solver | |
390 | // (z positions of the Read-out chambers, linearly interpolated) | |
391 | ||
392 | char hname[100]; | |
393 | if (side==0) sprintf(hname,"survey_dz_Aside"); | |
394 | if (side==1) sprintf(hname,"survey_dz_Cside"); | |
395 | ||
396 | TH2F *h = new TH2F(hname,hname,nx,-250.,250.,ny,-250.,250.); | |
397 | ||
398 | for (Int_t iy=1;iy<=ny;++iy) { | |
399 | Double_t yp = h->GetYaxis()->GetBinCenter(iy); | |
400 | for (Int_t ix=1;ix<=nx;++ix) { | |
401 | Double_t xp = h->GetXaxis()->GetBinCenter(ix); | |
402 | ||
403 | Float_t r0 = TMath::Sqrt(xp*xp+yp*yp); | |
404 | Float_t phi0 = TMath::ATan2(yp,xp); | |
405 | ||
406 | Float_t dz = GetROCVoltOffset(side,r0,phi0); // in [volt] | |
407 | ||
408 | Double_t ezField = (fgkCathodeV-fgkGG)/fgkTPCZ0; // = ALICE Electric Field (V/cm) Magnitude ~ -400 V/cm; | |
409 | dz = dz/ezField; // in [cm] | |
410 | ||
411 | if (85.<=r0 && r0<=245.) { | |
412 | h->SetBinContent(ix,iy,dz); | |
413 | } else { | |
414 | h->SetBinContent(ix,iy,0.); | |
415 | } | |
416 | } | |
417 | } | |
418 | ||
419 | h->GetXaxis()->SetTitle("x [cm]"); | |
420 | h->GetYaxis()->SetTitle("y [cm]"); | |
421 | h->GetZaxis()->SetTitle("dz [cm]"); | |
422 | h->SetStats(0); | |
423 | // h->DrawCopy("colz"); | |
424 | ||
425 | return h; | |
426 | } | |
427 | ||
428 | void AliTPCROCVoltError3D::Print(const Option_t* option) const { | |
429 | // | |
430 | // Print function to check the settings of the Rod shifts and the rotated clips | |
431 | // option=="a" prints the C0 and C1 coefficents for calibration purposes | |
432 | // | |
433 | ||
434 | TString opt = option; opt.ToLower(); | |
435 | printf("%s\n",GetTitle()); | |
436 | printf(" - Voltage settings on the TPC Read-Out chambers - linearly interpolated\n"); | |
437 | printf(" info: Check the following data-file for more details: %s \n",fROCDataFileName); | |
438 | ||
439 | if (opt.Contains("a")) { // Print all details | |
440 | printf(" - T1: %1.4f, T2: %1.4f \n",fT1,fT2); | |
441 | printf(" - C1: %1.4f, C0: %1.4f \n",fC1,fC0); | |
442 | } | |
443 | ||
444 | if (!fInitLookUp) AliError("Lookup table was not initialized! You should do InitROCVoltError3D() ..."); | |
445 | ||
446 | } |