]>
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 | // AliTPCFCVoltError3D class // | |
19 | // The class calculates the space point distortions due to residual voltage // | |
20 | // errors on the Field Cage boundaries (rods) 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 | // drift velocity changes in dependence of the electric field!!! // | |
26 | // // | |
27 | // date: 08/08/2010 // | |
28 | // Authors: Jim Thomas, Stefan Rossegger // | |
29 | // // | |
30 | // Example usage : // | |
31 | // AliTPCFCVoltError3D fcerror; // | |
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 | ||
41 | #include "TMath.h" | |
42 | #include "AliTPCROC.h" | |
43 | #include "AliTPCFCVoltError3D.h" | |
44 | ||
45 | ClassImp(AliTPCFCVoltError3D) | |
46 | ||
47 | AliTPCFCVoltError3D::AliTPCFCVoltError3D() | |
48 | : AliTPCCorrection("FieldCageVoltErrors","FieldCage (Rods) Voltage Errors"), | |
49 | fC0(0.),fC1(0.), | |
50 | fInitLookUp(kFALSE) | |
51 | { | |
52 | // | |
53 | // default constructor | |
54 | // | |
55 | ||
56 | // flags for filled 'basic lookup tables' | |
57 | for (Int_t i=0; i<5; i++){ | |
58 | fInitLookUpBasic[i]= kFALSE; | |
59 | } | |
60 | ||
61 | // Boundary settings | |
62 | for (Int_t i=0; i<36; i++){ | |
63 | fRodVoltShiftA[i] = 0; | |
64 | fRodVoltShiftC[i] = 0; | |
65 | } | |
66 | for (Int_t i=0; i<2; i++){ | |
67 | fRotatedClipVoltA[i] = 0; | |
68 | fRotatedClipVoltC[i] = 0; | |
69 | } | |
70 | // | |
71 | for (Int_t i=0; i<18; i++){ | |
72 | fOFCRodShiftA[i] = 0; | |
73 | fOFCRodShiftC[i] = 0; | |
74 | } | |
75 | ||
76 | // Array which will contain the solution according to the setted boundary conditions | |
77 | // it represents a sum up of the 4 basic look up tables (created individually) | |
78 | // see InitFCVoltError3D() function | |
79 | for ( Int_t k = 0 ; k < kNPhi ; k++ ) { | |
80 | fLookUpErOverEz[k] = new TMatrixD(kNR,kNZ); | |
81 | fLookUpEphiOverEz[k] = new TMatrixD(kNR,kNZ); | |
82 | fLookUpDeltaEz[k] = new TMatrixD(kNR,kNZ); | |
83 | } | |
84 | ||
85 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
86 | fLookUpBasic1ErOverEz[k] = 0; | |
87 | fLookUpBasic1EphiOverEz[k] = 0; | |
88 | fLookUpBasic1DeltaEz[k] = 0; | |
89 | ||
90 | fLookUpBasic2ErOverEz[k] = 0; | |
91 | fLookUpBasic2EphiOverEz[k] = 0; | |
92 | fLookUpBasic2DeltaEz[k] = 0; | |
93 | ||
94 | fLookUpBasic3ErOverEz[k] = 0; | |
95 | fLookUpBasic3EphiOverEz[k] = 0; | |
96 | fLookUpBasic3DeltaEz[k] = 0; | |
97 | ||
98 | fLookUpBasic4ErOverEz[k] = 0; | |
99 | fLookUpBasic4EphiOverEz[k] = 0; | |
100 | fLookUpBasic4DeltaEz[k] = 0; | |
101 | ||
102 | fLookUpBasic5ErOverEz[k] = 0; | |
103 | fLookUpBasic5EphiOverEz[k] = 0; | |
104 | fLookUpBasic5DeltaEz[k] = 0; | |
105 | } | |
106 | ||
107 | } | |
108 | ||
109 | AliTPCFCVoltError3D::~AliTPCFCVoltError3D() { | |
110 | // | |
111 | // destructor | |
112 | // | |
113 | ||
114 | for ( Int_t k = 0 ; k < kNPhi ; k++ ) { | |
115 | delete fLookUpErOverEz[k]; | |
116 | delete fLookUpEphiOverEz[k]; | |
117 | delete fLookUpDeltaEz[k]; | |
118 | } | |
119 | ||
120 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
121 | delete fLookUpBasic1ErOverEz[k]; // does nothing if pointer is zero! | |
122 | delete fLookUpBasic1EphiOverEz[k]; | |
123 | delete fLookUpBasic1DeltaEz[k]; | |
124 | ||
125 | delete fLookUpBasic2ErOverEz[k]; // does nothing if pointer is zero! | |
126 | delete fLookUpBasic2EphiOverEz[k]; | |
127 | delete fLookUpBasic2DeltaEz[k]; | |
128 | ||
129 | delete fLookUpBasic3ErOverEz[k]; // does nothing if pointer is zero! | |
130 | delete fLookUpBasic3EphiOverEz[k]; | |
131 | delete fLookUpBasic3DeltaEz[k]; | |
132 | ||
133 | delete fLookUpBasic4ErOverEz[k]; // does nothing if pointer is zero! | |
134 | delete fLookUpBasic4EphiOverEz[k]; | |
135 | delete fLookUpBasic4DeltaEz[k]; | |
136 | ||
137 | delete fLookUpBasic5ErOverEz[k]; // does nothing if pointer is zero! | |
138 | delete fLookUpBasic5EphiOverEz[k]; | |
139 | delete fLookUpBasic5DeltaEz[k]; | |
140 | } | |
141 | } | |
142 | ||
143 | void AliTPCFCVoltError3D::Init() { | |
144 | // | |
145 | // Initialization funtion | |
146 | // | |
147 | ||
148 | AliMagF* magF= (AliMagF*)TGeoGlobalMagField::Instance()->GetField(); | |
149 | if (!magF) AliError("Magneticd field - not initialized"); | |
150 | Double_t bzField = magF->SolenoidField()/10.; //field in T | |
151 | AliTPCParam *param= AliTPCcalibDB::Instance()->GetParameters(); | |
152 | if (!param) AliError("Parameters - not initialized"); | |
153 | Double_t vdrift = param->GetDriftV()/1000000.; // [cm/us] // From dataBase: to be updated: per second (ideally) | |
154 | Double_t ezField = 400; // [V/cm] // to be updated: never (hopefully) | |
155 | Double_t wt = -10.0 * (bzField*10) * vdrift / ezField ; | |
156 | // Correction Terms for effective omegaTau; obtained by a laser calibration run | |
157 | SetOmegaTauT1T2(wt,fT1,fT2); | |
158 | ||
159 | InitFCVoltError3D(); | |
160 | } | |
161 | ||
162 | void AliTPCFCVoltError3D::Update(const TTimeStamp &/*timeStamp*/) { | |
163 | // | |
164 | // Update function | |
165 | // | |
166 | AliMagF* magF= (AliMagF*)TGeoGlobalMagField::Instance()->GetField(); | |
167 | if (!magF) AliError("Magneticd field - not initialized"); | |
168 | Double_t bzField = magF->SolenoidField()/10.; //field in T | |
169 | AliTPCParam *param= AliTPCcalibDB::Instance()->GetParameters(); | |
170 | if (!param) AliError("Parameters - not initialized"); | |
171 | Double_t vdrift = param->GetDriftV()/1000000.; // [cm/us] // From dataBase: to be updated: per second (ideally) | |
172 | Double_t ezField = 400; // [V/cm] // to be updated: never (hopefully) | |
173 | Double_t wt = -10.0 * (bzField*10) * vdrift / ezField ; | |
174 | // Correction Terms for effective omegaTau; obtained by a laser calibration run | |
175 | SetOmegaTauT1T2(wt,fT1,fT2); | |
176 | ||
177 | ||
178 | } | |
179 | ||
180 | ||
181 | ||
182 | void AliTPCFCVoltError3D::GetCorrection(const Float_t x[],const Short_t roc,Float_t dx[]) { | |
183 | // | |
184 | // Calculates the correction due e.g. residual voltage errors on the TPC boundaries | |
185 | // | |
186 | ||
187 | if (!fInitLookUp) { | |
188 | AliInfo("Lookup table was not initialized! Perform the inizialisation now ..."); | |
189 | InitFCVoltError3D(); | |
190 | return; | |
191 | } | |
192 | ||
193 | Int_t order = 1 ; // FIXME: hardcoded? Linear interpolation = 1, Quadratic = 2 | |
194 | // note that the poisson solution was linearly mirroed on this grid! | |
195 | Double_t intEr, intEphi, intDeltaEz; | |
196 | Double_t r, phi, z ; | |
197 | Int_t sign; | |
198 | ||
199 | r = TMath::Sqrt( x[0]*x[0] + x[1]*x[1] ) ; | |
200 | phi = TMath::ATan2(x[1],x[0]) ; | |
201 | if ( phi < 0 ) phi += TMath::TwoPi() ; // Table uses phi from 0 to 2*Pi | |
202 | z = x[2] ; // Create temporary copy of x[2] | |
203 | ||
204 | if ( (roc%36) < 18 ) { | |
205 | sign = 1; // (TPC A side) | |
206 | } else { | |
207 | sign = -1; // (TPC C side) | |
208 | } | |
209 | ||
210 | if ( sign==1 && z < fgkZOffSet ) z = fgkZOffSet; // Protect against discontinuity at CE | |
211 | if ( sign==-1 && z > -fgkZOffSet ) z = -fgkZOffSet; // Protect against discontinuity at CE | |
212 | ||
213 | ||
214 | if ( (sign==1 && z<0) || (sign==-1 && z>0) ) // just a consistency check | |
215 | AliError("ROC number does not correspond to z coordinate! Calculation of distortions is most likely wrong!"); | |
216 | ||
217 | // Get the Er and Ephi field integrals plus the integral over DeltaEz | |
218 | intEr = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, | |
219 | fgkRList, fgkZList, fgkPhiList, fLookUpErOverEz ); | |
220 | intEphi = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, | |
221 | fgkRList, fgkZList, fgkPhiList, fLookUpEphiOverEz ); | |
222 | intDeltaEz = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, | |
223 | fgkRList, fgkZList, fgkPhiList, fLookUpDeltaEz ); | |
224 | ||
225 | // printf("%lf %lf %lf\n",intEr,intEphi,intDeltaEz); | |
226 | ||
227 | // Calculate distorted position | |
228 | if ( r > 0.0 ) { | |
229 | phi = phi + ( fC0*intEphi - fC1*intEr ) / r; | |
230 | r = r + ( fC0*intEr + fC1*intEphi ); | |
231 | } | |
232 | ||
233 | // Calculate correction in cartesian coordinates | |
234 | dx[0] = r * TMath::Cos(phi) - x[0]; | |
235 | dx[1] = r * TMath::Sin(phi) - x[1]; | |
236 | dx[2] = intDeltaEz; // z distortion - (internally scaled with driftvelocity dependency | |
237 | // on the Ez field plus the actual ROC misalignment (if set TRUE) | |
238 | ||
239 | } | |
240 | ||
241 | void AliTPCFCVoltError3D::InitFCVoltError3D() { | |
242 | // | |
243 | // Initialization of the Lookup table which contains the solutions of the | |
244 | // Dirichlet boundary problem | |
245 | // Calculation of the single 3D-Poisson solver is done just if needed | |
246 | // (see basic lookup tables in header file) | |
247 | // | |
248 | ||
249 | const Int_t order = 1 ; // Linear interpolation = 1, Quadratic = 2 | |
250 | const Float_t gridSizeR = (fgkOFCRadius-fgkIFCRadius) / (kRows-1) ; | |
251 | const Float_t gridSizeZ = fgkTPCZ0 / (kColumns-1) ; | |
252 | const Float_t gridSizePhi = TMath::TwoPi() / ( 18.0 * kPhiSlicesPerSector); | |
253 | ||
254 | // temporary arrays to create the boundary conditions | |
255 | TMatrixD *arrayofArrayV[kPhiSlices], *arrayofCharge[kPhiSlices] ; | |
256 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
257 | arrayofArrayV[k] = new TMatrixD(kRows,kColumns) ; | |
258 | arrayofCharge[k] = new TMatrixD(kRows,kColumns) ; | |
259 | } | |
260 | Double_t innerList[kPhiSlices], outerList[kPhiSlices] ; | |
261 | ||
262 | // list of point as used in the poisson relation and the interpolation (during sum up) | |
263 | Double_t rlist[kRows], zedlist[kColumns] , philist[kPhiSlices]; | |
264 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
265 | philist[k] = gridSizePhi * k; | |
266 | for ( Int_t i = 0 ; i < kRows ; i++ ) { | |
267 | rlist[i] = fgkIFCRadius + i*gridSizeR ; | |
268 | for ( Int_t j = 0 ; j < kColumns ; j++ ) { // Fill Vmatrix with Boundary Conditions | |
269 | zedlist[j] = j * gridSizeZ ; | |
270 | } | |
271 | } | |
272 | } | |
273 | ||
274 | // ========================================================================== | |
275 | // Solve Poisson's equation in 3D cylindrical coordinates by relaxation technique | |
276 | // Allow for different size grid spacing in R and Z directions | |
277 | ||
278 | const Int_t symmetry[5] = {1,1,-1,-1,1}; // shifted rod (2x), rotated clip (2x), only rod shift on OFC (1x) | |
279 | ||
280 | // check which lookup tables are needed --------------------------------- | |
281 | ||
282 | Bool_t needTable[5] ={kFALSE,kFALSE,kFALSE,kFALSE,kFALSE}; | |
283 | ||
284 | // Shifted rods & strips | |
285 | for ( Int_t rod = 0 ; rod < 18 ; rod++ ) { | |
286 | if (fRodVoltShiftA[rod]!=0) needTable[0]=kTRUE; | |
287 | if (fRodVoltShiftC[rod]!=0) needTable[0]=kTRUE; | |
288 | } | |
289 | for ( Int_t rod = 18 ; rod < 36 ; rod++ ) { | |
290 | if (fRodVoltShiftA[rod]!=0) needTable[1]=kTRUE; | |
291 | if (fRodVoltShiftC[rod]!=0) needTable[1]=kTRUE; | |
292 | } | |
293 | // Rotated clips | |
294 | if (fRotatedClipVoltA[0]!=0 || fRotatedClipVoltC[0]!=0) needTable[2]=kTRUE; | |
295 | if (fRotatedClipVoltA[1]!=0 || fRotatedClipVoltC[1]!=0) needTable[3]=kTRUE; | |
296 | ||
297 | // shifted Copper rods on OFC | |
298 | for ( Int_t rod = 0 ; rod < 18 ; rod++ ) { | |
299 | if (fOFCRodShiftA[rod]!=0) needTable[4]=kTRUE; | |
300 | if (fOFCRodShiftC[rod]!=0) needTable[4]=kTRUE; | |
301 | } | |
302 | ||
303 | ||
304 | // reserve the arrays for the basic lookup tables ---------------------- | |
305 | if (needTable[0] && !fInitLookUpBasic[0]) { | |
306 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { // Possibly make an extra table to be used for 0 == 360 | |
307 | fLookUpBasic1ErOverEz[k] = new TMatrixD(kRows,kColumns); | |
308 | fLookUpBasic1EphiOverEz[k] = new TMatrixD(kRows,kColumns); | |
309 | fLookUpBasic1DeltaEz[k] = new TMatrixD(kRows,kColumns); | |
310 | // will be deleted in destructor | |
311 | } | |
312 | } | |
313 | if (needTable[1] && !fInitLookUpBasic[1]) { | |
314 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { // Possibly make an extra table to be used for 0 == 360 | |
315 | fLookUpBasic2ErOverEz[k] = new TMatrixD(kRows,kColumns); | |
316 | fLookUpBasic2EphiOverEz[k] = new TMatrixD(kRows,kColumns); | |
317 | fLookUpBasic2DeltaEz[k] = new TMatrixD(kRows,kColumns); | |
318 | // will be deleted in destructor | |
319 | } | |
320 | } | |
321 | if (needTable[2] && !fInitLookUpBasic[2]) { | |
322 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { // Possibly make an extra table to be used for 0 == 360 | |
323 | fLookUpBasic3ErOverEz[k] = new TMatrixD(kRows,kColumns); | |
324 | fLookUpBasic3EphiOverEz[k] = new TMatrixD(kRows,kColumns); | |
325 | fLookUpBasic3DeltaEz[k] = new TMatrixD(kRows,kColumns); | |
326 | // will be deleted in destructor | |
327 | } | |
328 | } | |
329 | if (needTable[3] && !fInitLookUpBasic[3]) { | |
330 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { // Possibly make an extra table to be used for 0 == 360 | |
331 | fLookUpBasic4ErOverEz[k] = new TMatrixD(kRows,kColumns); | |
332 | fLookUpBasic4EphiOverEz[k] = new TMatrixD(kRows,kColumns); | |
333 | fLookUpBasic4DeltaEz[k] = new TMatrixD(kRows,kColumns); | |
334 | // will be deleted in destructor | |
335 | } | |
336 | } | |
337 | if (needTable[4] && !fInitLookUpBasic[4]) { | |
338 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { // Possibly make an extra table to be used for 0 == 360 | |
339 | fLookUpBasic5ErOverEz[k] = new TMatrixD(kRows,kColumns); | |
340 | fLookUpBasic5EphiOverEz[k] = new TMatrixD(kRows,kColumns); | |
341 | fLookUpBasic5DeltaEz[k] = new TMatrixD(kRows,kColumns); | |
342 | // will be deleted in destructor | |
343 | } | |
344 | } | |
345 | ||
346 | // Set bondaries and solve Poisson's equation -------------------------- | |
347 | ||
348 | for (Int_t look=0; look<5; look++) { | |
349 | ||
350 | Float_t inner18[18] = { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 } ; | |
351 | Float_t outer18[18] = { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 } ; | |
352 | ||
353 | if (needTable[look] && !fInitLookUpBasic[look]) { | |
354 | ||
355 | // Specify which rod is the reference; other distortions calculated by summing over multiple rotations of refrence | |
356 | // Note: the interpolation later on depends on it!! Do not change if not really needed! | |
357 | if (look==0) { | |
358 | AliInfo(Form("IFC - ROD&Strip shift : Filling table (~ %d sec)",3*(int)(kPhiSlices/5))); | |
359 | inner18[0] = 1; | |
360 | } else if (look==1) { | |
361 | AliInfo(Form("OFC - ROD&Strip shift : Filling table (~ %d sec)",3*(int)(kPhiSlices/5))); | |
362 | outer18[0] = 1; | |
363 | } else if (look==2) { | |
364 | AliInfo(Form("IFC - Clip rot. : Filling table (~ %d sec)",3*(int)(kPhiSlices/5))); | |
365 | inner18[0] = 1; | |
366 | } else if (look==3) { | |
367 | AliInfo(Form("OFC - Clip rot. : Filling table (~ %d sec)",3*(int)(kPhiSlices/5))); | |
368 | outer18[0] = 1; | |
369 | } else if (look==4) { | |
370 | AliInfo(Form("OFC - CopperRod shift : Filling table (~ %d sec)",3*(int)(kPhiSlices/5))); | |
371 | outer18[0] = 1; | |
372 | } | |
373 | // Build potentials on boundary stripes for a rotated clip (SYMMETRY==-1) or a shifted rod (SYMMETRY==1) | |
374 | if (look!=4) { | |
375 | // in these cases, the strips follow the actual rod shift (linear interpolation between the rods) | |
376 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
377 | Int_t slices = kPhiSlicesPerSector ; | |
378 | Int_t ipoint = k/slices ; | |
379 | innerList[k] = (((ipoint+1)*slices-k)*-1*inner18[ipoint]+(k-ipoint*slices)*inner18[ipoint+1])/slices ; | |
380 | outerList[k] = (((ipoint+1)*slices-k)*-1*outer18[ipoint]+(k-ipoint*slices)*outer18[ipoint+1])/slices ; | |
381 | if ( (k%slices) == 0 && symmetry[look] == -1 ) { innerList[k] = 0.0 ; outerList[k] = 0.0 ; } | |
382 | // above, force Zero if Anti-Sym | |
383 | } | |
384 | } else { | |
385 | // in this case, we assume that the strips stay at the correct position, but the copper plated OFC-rods move | |
386 | // the distortion is then much more localized around the rod (indicated by real data) | |
387 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) | |
388 | innerList[k] = outerList[k] = 0; | |
389 | outerList[0] = outer18[0]; // point at rod | |
390 | outerList[1] = outer18[0]/4; // point close to rod (educated-`guessed` scaling) | |
391 | } | |
392 | ||
393 | // Fill arrays with initial conditions. V on the boundary and Charge in the volume. | |
394 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
395 | TMatrixD &arrayV = *arrayofArrayV[k] ; | |
396 | TMatrixD &charge = *arrayofCharge[k] ; | |
397 | for ( Int_t i = 0 ; i < kRows ; i++ ) { | |
398 | for ( Int_t j = 0 ; j < kColumns ; j++ ) { // Fill Vmatrix with Boundary Conditions | |
399 | arrayV(i,j) = 0.0 ; | |
400 | charge(i,j) = 0.0 ; | |
401 | if ( i == 0 ) arrayV(i,j) = innerList[k] ; | |
402 | if ( i == kRows-1 ) arrayV(i,j) = outerList[k] ; | |
403 | } | |
404 | } | |
405 | // no charge in the volume | |
406 | for ( Int_t i = 1 ; i < kRows-1 ; i++ ) { | |
407 | for ( Int_t j = 1 ; j < kColumns-1 ; j++ ) { | |
408 | charge(i,j) = 0.0 ; | |
409 | } | |
410 | } | |
411 | } | |
412 | ||
413 | // Solve Poisson's equation in 3D cylindrical coordinates by relaxation technique | |
414 | // Allow for different size grid spacing in R and Z directions | |
415 | if (look==0) { | |
416 | PoissonRelaxation3D( arrayofArrayV, arrayofCharge, | |
417 | fLookUpBasic1ErOverEz, fLookUpBasic1EphiOverEz, fLookUpBasic1DeltaEz, | |
418 | kRows, kColumns, kPhiSlices, gridSizePhi, kIterations, symmetry[0]) ; | |
419 | AliInfo("IFC - ROD&Strip shift : done "); | |
420 | } else if (look==1) { | |
421 | PoissonRelaxation3D( arrayofArrayV, arrayofCharge, | |
422 | fLookUpBasic2ErOverEz, fLookUpBasic2EphiOverEz, fLookUpBasic2DeltaEz, | |
423 | kRows, kColumns, kPhiSlices, gridSizePhi, kIterations, symmetry[1]) ; | |
424 | ||
425 | AliInfo("OFC - ROD&Strip shift : done "); | |
426 | } else if (look==2) { | |
427 | PoissonRelaxation3D( arrayofArrayV, arrayofCharge, | |
428 | fLookUpBasic3ErOverEz, fLookUpBasic3EphiOverEz, fLookUpBasic3DeltaEz, | |
429 | kRows, kColumns, kPhiSlices, gridSizePhi, kIterations, symmetry[2]) ; | |
430 | AliInfo("IFC - Clip rot. : done "); | |
431 | } else if (look==3) { | |
432 | PoissonRelaxation3D( arrayofArrayV, arrayofCharge, | |
433 | fLookUpBasic4ErOverEz, fLookUpBasic4EphiOverEz, fLookUpBasic4DeltaEz, | |
434 | kRows, kColumns, kPhiSlices, gridSizePhi, kIterations, symmetry[3]) ; | |
435 | AliInfo("OFC - Clip rot. : done "); | |
436 | } else if (look==4) { | |
437 | PoissonRelaxation3D( arrayofArrayV, arrayofCharge, | |
438 | fLookUpBasic5ErOverEz, fLookUpBasic5EphiOverEz, fLookUpBasic5DeltaEz, | |
439 | kRows, kColumns, kPhiSlices, gridSizePhi, kIterations, symmetry[4]) ; | |
440 | AliInfo("OFC - CopperRod shift : done "); | |
441 | } | |
442 | ||
443 | fInitLookUpBasic[look] = kTRUE; | |
444 | } | |
445 | } | |
446 | ||
447 | ||
448 | // ============================================================================= | |
449 | // Create the final lookup table with corresponding ROD shifts or clip rotations | |
450 | ||
451 | Float_t erIntegralSum = 0.0 ; | |
452 | Float_t ephiIntegralSum = 0.0 ; | |
453 | Float_t ezIntegralSum = 0.0 ; | |
454 | ||
455 | Double_t phiPrime = 0. ; | |
456 | Double_t erIntegral = 0. ; | |
457 | Double_t ephiIntegral = 0. ; | |
458 | Double_t ezIntegral = 0. ; | |
459 | ||
460 | ||
461 | AliInfo("Calculation of combined Look-up Table"); | |
462 | ||
463 | // Interpolate and sum the Basic lookup tables onto the standard grid | |
464 | Double_t r, phi, z ; | |
465 | for ( Int_t k = 0 ; k < kNPhi ; k++ ) { | |
466 | phi = fgkPhiList[k] ; | |
467 | ||
468 | TMatrixD &erOverEz = *fLookUpErOverEz[k] ; | |
469 | TMatrixD &ephiOverEz = *fLookUpEphiOverEz[k]; | |
470 | TMatrixD &deltaEz = *fLookUpDeltaEz[k] ; | |
471 | ||
472 | for ( Int_t i = 0 ; i < kNZ ; i++ ) { | |
473 | z = TMath::Abs(fgkZList[i]) ; // Symmetric solution in Z that depends only on ABS(Z) | |
474 | for ( Int_t j = 0 ; j < kNR ; j++ ) { | |
475 | r = fgkRList[j] ; | |
476 | // Interpolate basicLookup tables; once for each rod, then sum the results | |
477 | ||
478 | erIntegralSum = 0.0 ; | |
479 | ephiIntegralSum = 0.0 ; | |
480 | ezIntegralSum = 0.0 ; | |
481 | ||
482 | // SHIFTED RODS ========================================================= | |
483 | ||
484 | // IFC ROD SHIFTS +++++++++++++++++++++++++++++ | |
485 | for ( Int_t rod = 0 ; rod < 18 ; rod++ ) { | |
486 | ||
487 | erIntegral = ephiIntegral = ezIntegral = 0.0 ; | |
488 | ||
489 | if ( fRodVoltShiftA[rod] == 0 && fgkZList[i] > 0) continue ; | |
490 | if ( fRodVoltShiftC[rod] == 0 && fgkZList[i] < 0) continue ; | |
491 | ||
492 | // Rotate to a position where we have maps and prepare to sum | |
493 | phiPrime = phi - rod*kPhiSlicesPerSector*gridSizePhi ; | |
494 | ||
495 | if ( phiPrime < 0 ) phiPrime += TMath::TwoPi() ; // Stay in range from 0 to TwoPi | |
496 | ||
497 | if ( (phiPrime >= 0) && (phiPrime <= kPhiSlices*gridSizePhi) ) { | |
498 | ||
499 | erIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
500 | rlist, zedlist, philist, fLookUpBasic1ErOverEz ); | |
501 | ephiIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
502 | rlist, zedlist, philist, fLookUpBasic1EphiOverEz); | |
503 | ezIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
504 | rlist, zedlist, philist, fLookUpBasic1DeltaEz ); | |
505 | ||
506 | } else if ( (phiPrime <= TMath::TwoPi()) && (phiPrime >= (TMath::TwoPi()-kPhiSlices*gridSizePhi)) ){ | |
507 | ||
508 | phiPrime = TMath::TwoPi() - phiPrime ; | |
509 | ||
510 | erIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
511 | rlist, zedlist, philist, fLookUpBasic1ErOverEz ); | |
512 | ephiIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
513 | rlist, zedlist, philist, fLookUpBasic1EphiOverEz); | |
514 | ezIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
515 | rlist, zedlist, philist, fLookUpBasic1DeltaEz ); | |
516 | ||
517 | // Flip symmetry of phi integral for symmetric boundary conditions | |
518 | if ( symmetry[0] == 1 ) ephiIntegral *= -1 ; | |
519 | // Flip symmetry of r integral if anti-symmetric boundary conditions | |
520 | if ( symmetry[0] == -1 ) erIntegral *= -1 ; | |
521 | ||
522 | } | |
523 | ||
524 | if ( fgkZList[i] > 0 ) { | |
525 | erIntegralSum += fRodVoltShiftA[rod]*erIntegral ; | |
526 | ephiIntegralSum += fRodVoltShiftA[rod]*ephiIntegral ; | |
527 | ezIntegralSum += fRodVoltShiftA[rod]*ezIntegral; | |
528 | } else if ( fgkZList[i] < 0 ) { | |
529 | erIntegralSum += fRodVoltShiftC[rod]*erIntegral ; | |
530 | ephiIntegralSum += fRodVoltShiftC[rod]*ephiIntegral ; | |
531 | ezIntegralSum -= fRodVoltShiftC[rod]*ezIntegral; | |
532 | } | |
533 | } | |
534 | ||
535 | // OFC ROD SHIFTS +++++++++++++++++++++++++++++ | |
536 | for ( Int_t rod = 18 ; rod < 36 ; rod++ ) { | |
537 | ||
538 | erIntegral = ephiIntegral = ezIntegral = 0.0 ; | |
539 | ||
540 | if ( fRodVoltShiftA[rod] == 0 && fgkZList[i] > 0) continue ; | |
541 | if ( fRodVoltShiftC[rod] == 0 && fgkZList[i] < 0) continue ; | |
542 | ||
543 | // Rotate to a position where we have maps and prepare to sum | |
544 | phiPrime = phi - (rod-18)*kPhiSlicesPerSector*gridSizePhi ; | |
545 | ||
546 | ||
547 | if ( phiPrime < 0 ) phiPrime += TMath::TwoPi() ; // Stay in range from 0 to TwoPi | |
548 | ||
549 | if ( (phiPrime >= 0) && (phiPrime <= kPhiSlices*gridSizePhi) ) { | |
550 | ||
551 | erIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
552 | rlist, zedlist, philist, fLookUpBasic2ErOverEz ); | |
553 | ephiIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
554 | rlist, zedlist, philist, fLookUpBasic2EphiOverEz); | |
555 | ezIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
556 | rlist, zedlist, philist, fLookUpBasic2DeltaEz ); | |
557 | ||
558 | } else if ( (phiPrime <= TMath::TwoPi()) && (phiPrime >= (TMath::TwoPi()-kPhiSlices*gridSizePhi)) ){ | |
559 | ||
560 | phiPrime = TMath::TwoPi() - phiPrime ; | |
561 | ||
562 | erIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
563 | rlist, zedlist, philist, fLookUpBasic2ErOverEz ); | |
564 | ephiIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
565 | rlist, zedlist, philist, fLookUpBasic2EphiOverEz); | |
566 | ezIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
567 | rlist, zedlist, philist, fLookUpBasic2DeltaEz ); | |
568 | ||
569 | // Flip symmetry of phi integral for symmetric boundary conditions | |
570 | if ( symmetry[1] == 1 ) ephiIntegral *= -1 ; | |
571 | // Flip symmetry of r integral if anti-symmetric boundary conditions | |
572 | if ( symmetry[1] == -1 ) erIntegral *= -1 ; | |
573 | ||
574 | } | |
575 | ||
576 | if ( fgkZList[i] > 0 ) { | |
577 | erIntegralSum += fRodVoltShiftA[rod]*erIntegral ; | |
578 | ephiIntegralSum += fRodVoltShiftA[rod]*ephiIntegral ; | |
579 | ezIntegralSum += fRodVoltShiftA[rod]*ezIntegral; | |
580 | } else if ( fgkZList[i] < 0 ) { | |
581 | erIntegralSum += fRodVoltShiftC[rod]*erIntegral ; | |
582 | ephiIntegralSum += fRodVoltShiftC[rod]*ephiIntegral ; | |
583 | ezIntegralSum -= fRodVoltShiftC[rod]*ezIntegral; | |
584 | } | |
585 | ||
586 | } // rod loop - shited ROD | |
587 | ||
588 | ||
589 | // Rotated clips ========================================================= | |
590 | ||
591 | Int_t rodIFC = 11; // resistor rod positions, IFC | |
592 | Int_t rodOFC = 3; // resistor rod positions, OFC | |
593 | // just one rod on IFC and OFC | |
594 | ||
595 | // IFC ROTATED CLIP +++++++++++++++++++++++++++++ | |
596 | for ( Int_t rod = rodIFC ; rod < rodIFC+1 ; rod++ ) { // loop over 1 to keep the "ignore"-functionality | |
597 | ||
598 | erIntegral = ephiIntegral = ezIntegral = 0.0 ; | |
599 | if ( fRotatedClipVoltA[0] == 0 && fgkZList[i] > 0) continue ; | |
600 | if ( fRotatedClipVoltC[0] == 0 && fgkZList[i] < 0) continue ; | |
601 | ||
602 | // Rotate to a position where we have maps and prepare to sum | |
603 | phiPrime = phi - rod*kPhiSlicesPerSector*gridSizePhi ; | |
604 | ||
605 | if ( phiPrime < 0 ) phiPrime += TMath::TwoPi() ; // Stay in range from 0 to TwoPi | |
606 | ||
607 | if ( (phiPrime >= 0) && (phiPrime <= kPhiSlices*gridSizePhi) ) { | |
608 | ||
609 | erIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
610 | rlist, zedlist, philist, fLookUpBasic3ErOverEz ); | |
611 | ephiIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
612 | rlist, zedlist, philist, fLookUpBasic3EphiOverEz); | |
613 | ezIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
614 | rlist, zedlist, philist, fLookUpBasic3DeltaEz ); | |
615 | ||
616 | } else if ( (phiPrime <= TMath::TwoPi()) && (phiPrime >= (TMath::TwoPi()-kPhiSlices*gridSizePhi)) ){ | |
617 | ||
618 | phiPrime = TMath::TwoPi() - phiPrime ; | |
619 | ||
620 | erIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
621 | rlist, zedlist, philist, fLookUpBasic3ErOverEz ); | |
622 | ephiIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
623 | rlist, zedlist, philist, fLookUpBasic3EphiOverEz); | |
624 | ezIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
625 | rlist, zedlist, philist, fLookUpBasic3DeltaEz ); | |
626 | ||
627 | // Flip symmetry of phi integral for symmetric boundary conditions | |
628 | if ( symmetry[2] == 1 ) ephiIntegral *= -1 ; | |
629 | // Flip symmetry of r integral if anti-symmetric boundary conditions | |
630 | if ( symmetry[2] == -1 ) erIntegral *= -1 ; | |
631 | ||
632 | } | |
633 | ||
634 | if ( fgkZList[i] > 0 ) { | |
635 | erIntegralSum += fRotatedClipVoltA[0]*erIntegral ; | |
636 | ephiIntegralSum += fRotatedClipVoltA[0]*ephiIntegral ; | |
637 | ezIntegralSum += fRotatedClipVoltA[0]*ezIntegral; | |
638 | } else if ( fgkZList[i] < 0 ) { | |
639 | erIntegralSum += fRotatedClipVoltC[0]*erIntegral ; | |
640 | ephiIntegralSum += fRotatedClipVoltC[0]*ephiIntegral ; | |
641 | ezIntegralSum -= fRotatedClipVoltC[0]*ezIntegral; | |
642 | } | |
643 | } | |
644 | ||
645 | // OFC: ROTATED CLIP +++++++++++++++++++++++++++++ | |
646 | for ( Int_t rod = rodOFC ; rod < rodOFC+1 ; rod++ ) { // loop over 1 to keep the "ignore"-functionality | |
647 | ||
648 | erIntegral = ephiIntegral = ezIntegral = 0.0 ; | |
649 | ||
650 | if ( fRotatedClipVoltA[1] == 0 && fgkZList[i] > 0) continue ; | |
651 | if ( fRotatedClipVoltC[1] == 0 && fgkZList[i] < 0) continue ; | |
652 | ||
653 | // Rotate to a position where we have maps and prepare to sum | |
654 | phiPrime = phi - rod*kPhiSlicesPerSector*gridSizePhi ; | |
655 | ||
656 | ||
657 | if ( phiPrime < 0 ) phiPrime += TMath::TwoPi() ; // Stay in range from 0 to TwoPi | |
658 | ||
659 | if ( (phiPrime >= 0) && (phiPrime <= kPhiSlices*gridSizePhi) ) { | |
660 | ||
661 | erIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
662 | rlist, zedlist, philist, fLookUpBasic4ErOverEz ); | |
663 | ephiIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
664 | rlist, zedlist, philist, fLookUpBasic4EphiOverEz); | |
665 | ezIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
666 | rlist, zedlist, philist, fLookUpBasic4DeltaEz ); | |
667 | ||
668 | } else if ( (phiPrime <= TMath::TwoPi()) && (phiPrime >= (TMath::TwoPi()-kPhiSlices*gridSizePhi)) ){ | |
669 | ||
670 | phiPrime = TMath::TwoPi() - phiPrime ; | |
671 | ||
672 | erIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
673 | rlist, zedlist, philist, fLookUpBasic4ErOverEz ); | |
674 | ephiIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
675 | rlist, zedlist, philist, fLookUpBasic4EphiOverEz); | |
676 | ezIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
677 | rlist, zedlist, philist, fLookUpBasic4DeltaEz ); | |
678 | ||
679 | // Flip symmetry of phi integral for symmetric boundary conditions | |
680 | if ( symmetry[3] == 1 ) ephiIntegral *= -1 ; | |
681 | // Flip symmetry of r integral if anti-symmetric boundary conditions | |
682 | if ( symmetry[3] == -1 ) erIntegral *= -1 ; | |
683 | ||
684 | } | |
685 | ||
686 | if ( fgkZList[i] > 0 ) { | |
687 | erIntegralSum += fRotatedClipVoltA[1]*erIntegral ; | |
688 | ephiIntegralSum += fRotatedClipVoltA[1]*ephiIntegral ; | |
689 | ezIntegralSum += fRotatedClipVoltA[1]*ezIntegral; | |
690 | } else if ( fgkZList[i] < 0 ) { | |
691 | erIntegralSum += fRotatedClipVoltC[1]*erIntegral ; | |
692 | ephiIntegralSum += fRotatedClipVoltC[1]*ephiIntegral ; | |
693 | ezIntegralSum -= fRotatedClipVoltC[1]*ezIntegral; | |
694 | } | |
695 | } | |
696 | ||
697 | // OFC Cooper rod shift +++++++++++++++++++++++++++++ | |
698 | for ( Int_t rod = 0 ; rod < 18 ; rod++ ) { | |
699 | ||
700 | erIntegral = ephiIntegral = ezIntegral = 0.0 ; | |
701 | ||
702 | if ( fOFCRodShiftA[rod] == 0 && fgkZList[i] > 0) continue ; | |
703 | if ( fOFCRodShiftC[rod] == 0 && fgkZList[i] < 0) continue ; | |
704 | ||
705 | // Rotate to a position where we have maps and prepare to sum | |
706 | phiPrime = phi - rod*kPhiSlicesPerSector*gridSizePhi ; | |
707 | ||
708 | if ( phiPrime < 0 ) phiPrime += TMath::TwoPi() ; // Stay in range from 0 to TwoPi | |
709 | ||
710 | if ( (phiPrime >= 0) && (phiPrime <= kPhiSlices*gridSizePhi) ) { | |
711 | ||
712 | erIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
713 | rlist, zedlist, philist, fLookUpBasic5ErOverEz ); | |
714 | ephiIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
715 | rlist, zedlist, philist, fLookUpBasic5EphiOverEz); | |
716 | ezIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
717 | rlist, zedlist, philist, fLookUpBasic5DeltaEz ); | |
718 | ||
719 | } else if ( (phiPrime <= TMath::TwoPi()) && (phiPrime >= (TMath::TwoPi()-kPhiSlices*gridSizePhi)) ){ | |
720 | ||
721 | phiPrime = TMath::TwoPi() - phiPrime ; | |
722 | ||
723 | erIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
724 | rlist, zedlist, philist, fLookUpBasic5ErOverEz ); | |
725 | ephiIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
726 | rlist, zedlist, philist, fLookUpBasic5EphiOverEz); | |
727 | ezIntegral = Interpolate3DTable(order, r, z, phiPrime, kRows, kColumns, kPhiSlices, | |
728 | rlist, zedlist, philist, fLookUpBasic5DeltaEz ); | |
729 | ||
730 | // Flip symmetry of phi integral for symmetric boundary conditions | |
731 | if ( symmetry[0] == 1 ) ephiIntegral *= -1 ; | |
732 | // Flip symmetry of r integral if anti-symmetric boundary conditions | |
733 | if ( symmetry[0] == -1 ) erIntegral *= -1 ; | |
734 | ||
735 | } | |
736 | ||
737 | if ( fgkZList[i] > 0 ) { | |
738 | erIntegralSum += fOFCRodShiftA[rod]*erIntegral ; | |
739 | ephiIntegralSum += fOFCRodShiftA[rod]*ephiIntegral ; | |
740 | ezIntegralSum += fOFCRodShiftA[rod]*ezIntegral; | |
741 | } else if ( fgkZList[i] < 0 ) { | |
742 | erIntegralSum += fOFCRodShiftC[rod]*erIntegral ; | |
743 | ephiIntegralSum += fOFCRodShiftC[rod]*ephiIntegral ; | |
744 | ezIntegralSum -= fOFCRodShiftC[rod]*ezIntegral; | |
745 | } | |
746 | } | |
747 | ||
748 | // put the sum into the final lookup table | |
749 | erOverEz(j,i) = erIntegralSum; | |
750 | ephiOverEz(j,i) = ephiIntegralSum; | |
751 | deltaEz(j,i) = ezIntegralSum; | |
752 | ||
753 | // if (j==1) printf("%lf %lf %lf: %lf %lf %lf\n",r, phi, z, erIntegralSum,ephiIntegralSum,ezIntegralSum); | |
754 | ||
755 | } // endo r loop | |
756 | } // end of z loop | |
757 | } // end of phi loop | |
758 | ||
759 | ||
760 | // clear the temporary arrays lists | |
761 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
762 | delete arrayofArrayV[k]; | |
763 | delete arrayofCharge[k]; | |
764 | } | |
765 | ||
766 | AliInfo(" done"); | |
767 | fInitLookUp = kTRUE; | |
768 | ||
769 | } | |
770 | ||
771 | void AliTPCFCVoltError3D::Print(const Option_t* option) const { | |
772 | // | |
773 | // Print function to check the settings of the Rod shifts and the rotated clips | |
774 | // option=="a" prints the C0 and C1 coefficents for calibration purposes | |
775 | // | |
776 | ||
777 | TString opt = option; opt.ToLower(); | |
778 | printf("%s\n",GetTitle()); | |
779 | printf(" - ROD shifts (residual voltage settings): 40V correspond to 1mm of shift\n"); | |
780 | Int_t count = 0; | |
781 | printf(" : A-side - IFC (Rod & Strips) \n "); | |
782 | for (Int_t i=0; i<18;i++) { | |
783 | if (fRodVoltShiftA[i]!=0) { | |
784 | printf("Rod %2d: %3.1f V \t",i,fRodVoltShiftA[i]); | |
785 | count++; | |
786 | } | |
787 | if (count==0&&i==17) | |
788 | printf("-> all at 0 V\n"); | |
789 | else if (i==17){ | |
790 | printf("\n"); | |
791 | count=0; | |
792 | } | |
793 | } | |
794 | printf(" : C-side - IFC (Rod & Strips) \n "); | |
795 | for (Int_t i=0; i<18;i++) { | |
796 | if (fRodVoltShiftC[i]!=0) { | |
797 | printf("Rod %2d: %3.1f V \t",i,fRodVoltShiftC[i]); | |
798 | count++; | |
799 | } | |
800 | if (count==0&&i==17) | |
801 | printf("-> all at 0 V\n"); | |
802 | else if (i==17){ | |
803 | printf("\n"); | |
804 | count=0; | |
805 | } | |
806 | } | |
807 | printf(" : A-side - OFC (only Rod) \n "); | |
808 | for (Int_t i=18; i<36;i++) { | |
809 | if (fRodVoltShiftA[i]!=0) { | |
810 | printf("Rod %2d: %3.1f V \t",i,fRodVoltShiftA[i]); | |
811 | count++; | |
812 | } | |
813 | if (count==0&&i==35) | |
814 | printf("-> all at 0 V\n"); | |
815 | else if (i==35) { | |
816 | printf("\n"); | |
817 | count=0; | |
818 | } | |
819 | } | |
820 | printf(" : C-side - OFC (only Rod) \n "); | |
821 | for (Int_t i=18; i<36;i++) { | |
822 | if (fRodVoltShiftC[i]!=0) { | |
823 | printf("Rod %2d: %3.1f V \t",i,fRodVoltShiftC[i]); | |
824 | count++; | |
825 | } | |
826 | if (count==0&&i==35) | |
827 | printf("-> all at 0 V\n"); | |
828 | else if (i==35){ | |
829 | printf("\n"); | |
830 | count=0; | |
831 | } | |
832 | } | |
833 | ||
834 | printf(" - Rotated clips (residual voltage settings): 40V correspond to 1mm of 'offset'\n"); | |
835 | if (fRotatedClipVoltA[0]!=0) { printf(" A side (IFC): HV Rod: %3.1f V \n",fRotatedClipVoltA[0]); count++; } | |
836 | if (fRotatedClipVoltA[1]!=0) { printf(" A side (OFC): HV Rod: %3.1f V \n",fRotatedClipVoltA[1]); count++; } | |
837 | if (fRotatedClipVoltC[0]!=0) { printf(" C side (IFC): HV Rod: %3.1f V \n",fRotatedClipVoltC[0]); count++; } | |
838 | if (fRotatedClipVoltC[1]!=0) { printf(" C side (OFC): HV Rod: %3.1f V \n",fRotatedClipVoltC[1]); count++; } | |
839 | if (count==0) | |
840 | printf(" -> no rotated clips \n"); | |
841 | ||
842 | count=0; | |
843 | printf(" - Copper ROD shifts (without strips):\n"); | |
844 | printf(" : A-side - OFC (Copper Rod shift) \n "); | |
845 | for (Int_t i=0; i<18;i++) { | |
846 | if (fOFCRodShiftA[i]!=0) { | |
847 | printf("Rod %2d: %3.1f V \t",i,fOFCRodShiftA[i]); | |
848 | count++; | |
849 | } | |
850 | if (count==0&&i==17) | |
851 | printf("-> all at 0 V\n"); | |
852 | else if (i==17){ | |
853 | printf("\n"); | |
854 | count=0; | |
855 | } | |
856 | } | |
857 | printf(" : C-side - OFC (Copper Rod shift) \n "); | |
858 | for (Int_t i=0; i<18;i++) { | |
859 | if (fOFCRodShiftC[i]!=0) { | |
860 | printf("Rod %2d: %3.1f V \t",i,fOFCRodShiftC[i]); | |
861 | count++; | |
862 | } | |
863 | if (count==0&&i==17) | |
864 | printf("-> all at 0 V\n"); | |
865 | else if (i==17){ | |
866 | printf("\n"); | |
867 | count=0; | |
868 | } | |
869 | } | |
870 | ||
871 | if (opt.Contains("a")) { // Print all details | |
872 | printf(" - T1: %1.4f, T2: %1.4f \n",fT1,fT2); | |
873 | printf(" - C1: %1.4f, C0: %1.4f \n",fC1,fC0); | |
874 | } | |
875 | ||
876 | if (!fInitLookUp) AliError("Lookup table was not initialized! You should do InitFCVoltError3D() ..."); | |
877 | ||
878 | } |