]>
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 | ||
b4caed64 | 16 | // _________________________________________________________________ |
17 | // | |
18 | // Begin_Html | |
19 | // <h2> AliTPCROCVoltError3D class </h2> | |
20 | // The class calculates the space point distortions due to z offsets of the | |
21 | // ROCs via the residual voltage technique (Poisson relaxation) in 3D. | |
22 | // Since the GG (part of the ROCs) represents the closure of the FC in z direction, | |
23 | // every misalignment in z produces not only dz distortions but also electrical | |
24 | // field inhomogeneities throughout the volume, which produces additional dr and rd$\phi$ distortions. | |
25 | // <p> | |
26 | // Each ROC can be misaligned (in z direction) in three ways. A general z0 offset, | |
27 | // an inclination along the x and an inclination along the y axis. The z-misalignment's | |
28 | // can be set via the function SetROCData(TMatrixD *mat) for each single chamber (ROC). | |
29 | // The array size has to be (72,3) representing the 72 chambers according to the | |
30 | // offline numbering scheme (IROC: roc$<$36; OROC: roc$\geq$36) and the three misalignment's | |
31 | // (see the source code for further details). | |
32 | // <p> | |
33 | // Internally, these z offsets (unit is cm) are recalculated into residual voltage | |
34 | // equivalents in order to make use of the relaxation technique. | |
35 | // <p> | |
36 | // One has two possibilities when calculating the $dz$ distortions. The resulting | |
37 | // distortions are purely due to the change of the drift velocity (along with the | |
38 | // change of the drift field) when the SetROCDisplacement is FALSE. <br> | |
39 | // For this class, this is a rather unphysical setting and should be avoided. <br> | |
40 | // When the flag is set to TRUE, the corresponding offset in z is added to the dz | |
41 | // calculation of the outer ROCs. <br> | |
42 | // For the Alice TPC gas, both effects are of similar magnitude. This means, if the | |
43 | // drift length is sufficiently large, a z-offset of a chamber appears to have (approx.) | |
44 | // twice the magnitude when one looks at the actual dz distortions. | |
45 | // <p> | |
46 | // In addition, this class allows a correction regarding the different arrival times | |
47 | // of the electrons due to the geometrical difference of the inner and outer chambers. | |
48 | // The impact was simulated via Garfield. If the flag is set via the | |
49 | // function SetElectronArrivalCorrection, the electron-arrival correction is added to the dz calculation. | |
50 | // End_Html | |
51 | // | |
6a1caa6b | 52 | // Begin_Macro(source) |
b4caed64 | 53 | // { |
54 | // gROOT->SetStyle("Plain"); gStyle->SetPalette(1); | |
6a1caa6b | 55 | // TCanvas *c2 = new TCanvas("cAliTPCROCVoltError3D","cAliTPCROCVoltError3D",500,400); |
b4caed64 | 56 | // AliTPCROCVoltError3D roc; |
919db3ce | 57 | // roc.SetROCDataFileName("$ALICE_ROOT/TPC/Calib/maps/TPCROCdzSurvey.root"); |
b4caed64 | 58 | // roc.SetElectronArrivalCorrection(kFALSE); // Correction for electron arrival offset, IROC vs OROC |
59 | // roc.SetROCDisplacement(kTRUE); // include the chamber offset in z when calculating the dz | |
60 | // roc.SetOmegaTauT1T2(0,1,1); // B=0 | |
61 | // roc.CreateHistoDZinXY(1.,300,300)->Draw("colz"); | |
62 | // return c2; | |
63 | // } | |
64 | // End_Macro | |
65 | // | |
66 | // Begin_Html | |
67 | // <p> | |
68 | // Date: 08/08/2010 <br> | |
69 | // Authors: Jim Thomas, Stefan Rossegger | |
70 | // End_Html | |
71 | // _________________________________________________________________ | |
72 | ||
c9cbd2f2 | 73 | |
74 | #include "AliMagF.h" | |
75 | #include "TGeoGlobalMagField.h" | |
76 | #include "AliTPCcalibDB.h" | |
77 | #include "AliTPCParam.h" | |
78 | #include "AliLog.h" | |
79 | #include "TMatrixD.h" | |
80 | #include "TFile.h" | |
81 | ||
82 | #include "TMath.h" | |
83 | #include "AliTPCROC.h" | |
84 | #include "AliTPCROCVoltError3D.h" | |
85 | ||
86 | ClassImp(AliTPCROCVoltError3D) | |
87 | ||
88 | AliTPCROCVoltError3D::AliTPCROCVoltError3D() | |
89 | : AliTPCCorrection("ROCVoltErrors","ROC z alignment Errors"), | |
90 | fC0(0.),fC1(0.), | |
91 | fROCdisplacement(kTRUE), | |
7a348589 | 92 | fElectronArrivalCorrection(kTRUE), |
c9cbd2f2 | 93 | fInitLookUp(kFALSE), |
acf5907b | 94 | fROCDataFileName(""), |
c9cbd2f2 | 95 | fdzDataLinFit(0) |
96 | { | |
97 | // | |
98 | // default constructor | |
99 | // | |
100 | ||
101 | // Array which will contain the solution according to the setted boundary conditions | |
102 | // main input: z alignment of the Read Out chambers | |
103 | // see InitROCVoltError3D() function | |
104 | for ( Int_t k = 0 ; k < kNPhi ; k++ ) { | |
2bf29b72 | 105 | fLookUpErOverEz[k] = new TMatrixF(kNR,kNZ); |
106 | fLookUpEphiOverEz[k] = new TMatrixF(kNR,kNZ); | |
107 | fLookUpDeltaEz[k] = new TMatrixF(kNR,kNZ); | |
c9cbd2f2 | 108 | } |
41afd7d3 | 109 | // fROCDataFileName="$ALICE_ROOT/TPC/Calib/maps/TPCROCdzSurvey.root"; |
110 | // SetROCDataFileName(fROCDataFileName.Data()); // initialization of fdzDataLinFit is included | |
c9cbd2f2 | 111 | |
112 | } | |
113 | ||
114 | AliTPCROCVoltError3D::~AliTPCROCVoltError3D() { | |
115 | // | |
116 | // destructor | |
117 | // | |
118 | ||
119 | for ( Int_t k = 0 ; k < kNPhi ; k++ ) { | |
120 | delete fLookUpErOverEz[k]; | |
121 | delete fLookUpEphiOverEz[k]; | |
122 | delete fLookUpDeltaEz[k]; | |
123 | } | |
124 | ||
125 | delete fdzDataLinFit; | |
126 | } | |
127 | ||
acf5907b | 128 | void AliTPCROCVoltError3D::SetROCData(TMatrixD * matrix){ |
2bbac918 | 129 | // |
c756f562 | 130 | // Set a z alignment map of the chambers not via a file, but |
7a348589 | 131 | // directly via a TMatrix(72,3), where dz = p0 + p1*(lx-133.4) + p2*ly (all in cm) |
2bbac918 | 132 | // |
133 | if (!fdzDataLinFit) fdzDataLinFit=new TMatrixD(*matrix); | |
134 | else *fdzDataLinFit = *matrix; | |
135 | } | |
136 | ||
137 | ||
c9cbd2f2 | 138 | void AliTPCROCVoltError3D::Init() { |
139 | // | |
140 | // Initialization funtion | |
141 | // | |
142 | ||
143 | AliMagF* magF= (AliMagF*)TGeoGlobalMagField::Instance()->GetField(); | |
144 | if (!magF) AliError("Magneticd field - not initialized"); | |
145 | Double_t bzField = magF->SolenoidField()/10.; //field in T | |
146 | AliTPCParam *param= AliTPCcalibDB::Instance()->GetParameters(); | |
147 | if (!param) AliError("Parameters - not initialized"); | |
148 | Double_t vdrift = param->GetDriftV()/1000000.; // [cm/us] // From dataBase: to be updated: per second (ideally) | |
149 | Double_t ezField = 400; // [V/cm] // to be updated: never (hopefully) | |
150 | Double_t wt = -10.0 * (bzField*10) * vdrift / ezField ; | |
151 | // Correction Terms for effective omegaTau; obtained by a laser calibration run | |
152 | SetOmegaTauT1T2(wt,fT1,fT2); | |
153 | ||
35ae345f | 154 | if (!fInitLookUp) InitROCVoltError3D(); |
c9cbd2f2 | 155 | } |
156 | ||
157 | void AliTPCROCVoltError3D::Update(const TTimeStamp &/*timeStamp*/) { | |
158 | // | |
159 | // Update function | |
160 | // | |
161 | AliMagF* magF= (AliMagF*)TGeoGlobalMagField::Instance()->GetField(); | |
162 | if (!magF) AliError("Magneticd field - not initialized"); | |
163 | Double_t bzField = magF->SolenoidField()/10.; //field in T | |
164 | AliTPCParam *param= AliTPCcalibDB::Instance()->GetParameters(); | |
165 | if (!param) AliError("Parameters - not initialized"); | |
166 | Double_t vdrift = param->GetDriftV()/1000000.; // [cm/us] // From dataBase: to be updated: per second (ideally) | |
167 | Double_t ezField = 400; // [V/cm] // to be updated: never (hopefully) | |
168 | Double_t wt = -10.0 * (bzField*10) * vdrift / ezField ; | |
169 | // Correction Terms for effective omegaTau; obtained by a laser calibration run | |
170 | SetOmegaTauT1T2(wt,fT1,fT2); | |
171 | ||
172 | } | |
173 | ||
acf5907b | 174 | void AliTPCROCVoltError3D::SetROCDataFileName(const char * fname) { |
c9cbd2f2 | 175 | // |
176 | // Set / load the ROC data (linear fit of ROC misalignments) | |
177 | // | |
178 | ||
179 | fROCDataFileName = fname; | |
180 | ||
acf5907b | 181 | TFile f(fROCDataFileName.Data(),"READ"); |
c9cbd2f2 | 182 | TMatrixD *m = (TMatrixD*) f.Get("dzSurveyLinFitData"); |
183 | TMatrixD &mf = *m; | |
184 | ||
185 | // prepare some space | |
186 | ||
187 | if (fdzDataLinFit) delete fdzDataLinFit; | |
188 | fdzDataLinFit = new TMatrixD(72,3); | |
189 | TMatrixD &dataIntern = *fdzDataLinFit; | |
190 | ||
191 | for (Int_t iroc=0;iroc<72;iroc++) { | |
192 | dataIntern(iroc,0) = mf(iroc,0); // z0 offset | |
193 | dataIntern(iroc,1) = mf(iroc,1); // slope in x | |
194 | dataIntern(iroc,2) = mf(iroc,2); // slope in y | |
195 | } | |
196 | ||
197 | f.Close(); | |
198 | ||
199 | fInitLookUp = kFALSE; | |
200 | ||
201 | } | |
202 | ||
203 | void AliTPCROCVoltError3D::GetCorrection(const Float_t x[],const Short_t roc,Float_t dx[]) { | |
204 | // | |
205 | // Calculates the correction due e.g. residual voltage errors on the TPC boundaries | |
206 | // | |
2bf29b72 | 207 | const Double_t kEpsilon=Double_t(FLT_MIN); |
c9cbd2f2 | 208 | if (!fInitLookUp) { |
209 | AliInfo("Lookup table was not initialized! Perform the inizialisation now ..."); | |
210 | InitROCVoltError3D(); | |
2bf29b72 | 211 | } |
212 | static Bool_t forceInit=kTRUE; //temporary needed for back compatibility with old OCDB entries | |
213 | if (forceInit&&fLookUpErOverEz[0]){ | |
32150d2c | 214 | if (TMath::Abs(fLookUpErOverEz[0]->Sum())<kEpsilon){//temporary needed for back compatibility with old OCDB entries |
2bf29b72 | 215 | ForceInitROCVoltError3D(); |
216 | } | |
217 | forceInit=kFALSE; | |
c9cbd2f2 | 218 | } |
219 | ||
2bf29b72 | 220 | |
7a348589 | 221 | Int_t order = 1 ; // FIXME: hardcoded? Linear interpolation = 1, Quadratic = 2 |
c9cbd2f2 | 222 | |
2bf29b72 | 223 | Float_t intEr, intEphi, intDeltaEz; |
c9cbd2f2 | 224 | Double_t r, phi, z ; |
225 | Int_t sign; | |
226 | ||
227 | r = TMath::Sqrt( x[0]*x[0] + x[1]*x[1] ) ; | |
228 | phi = TMath::ATan2(x[1],x[0]) ; | |
229 | if ( phi < 0 ) phi += TMath::TwoPi() ; // Table uses phi from 0 to 2*Pi | |
230 | z = x[2] ; // Create temporary copy of x[2] | |
231 | ||
232 | if ( (roc%36) < 18 ) { | |
233 | sign = 1; // (TPC A side) | |
234 | } else { | |
235 | sign = -1; // (TPC C side) | |
236 | } | |
237 | ||
238 | if ( sign==1 && z < fgkZOffSet ) z = fgkZOffSet; // Protect against discontinuity at CE | |
239 | if ( sign==-1 && z > -fgkZOffSet ) z = -fgkZOffSet; // Protect against discontinuity at CE | |
240 | ||
241 | ||
242 | if ( (sign==1 && z<0) || (sign==-1 && z>0) ) // just a consistency check | |
243 | AliError("ROC number does not correspond to z coordinate! Calculation of distortions is most likely wrong!"); | |
244 | ||
245 | // Get the Er and Ephi field integrals plus the integral over DeltaEz | |
246 | intEr = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, | |
247 | fgkRList, fgkZList, fgkPhiList, fLookUpErOverEz ); | |
248 | intEphi = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, | |
249 | fgkRList, fgkZList, fgkPhiList, fLookUpEphiOverEz); | |
250 | intDeltaEz = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, | |
251 | fgkRList, fgkZList, fgkPhiList, fLookUpDeltaEz ); | |
252 | ||
253 | // printf("%lf %lf %lf\n",intEr,intEphi,intDeltaEz); | |
254 | ||
255 | // Calculate distorted position | |
256 | if ( r > 0.0 ) { | |
257 | phi = phi + ( fC0*intEphi - fC1*intEr ) / r; | |
258 | r = r + ( fC0*intEr + fC1*intEphi ); | |
259 | } | |
260 | ||
261 | // Calculate correction in cartesian coordinates | |
262 | dx[0] = r * TMath::Cos(phi) - x[0]; | |
263 | dx[1] = r * TMath::Sin(phi) - x[1]; | |
264 | dx[2] = intDeltaEz; // z distortion - (internally scaled with driftvelocity dependency | |
265 | // on the Ez field plus the actual ROC misalignment (if set TRUE) | |
266 | ||
7a348589 | 267 | |
268 | if (fElectronArrivalCorrection) { | |
269 | ||
270 | // correction for the OROC (in average, a 0.014usec longer drift time | |
271 | // due to different position of the anode wires) -> vd*dt -> 2.64*0.014 = 0.0369 cm | |
272 | // FIXME: correction are token from Magboltz simulations | |
273 | // should be loaded from a database | |
274 | ||
275 | AliTPCROC * rocInfo = AliTPCROC::Instance(); | |
276 | Double_t rCrossingROC = (rocInfo->GetPadRowRadii(0,62)+rocInfo->GetPadRowRadii(36,0))/2; | |
277 | ||
278 | if (r>rCrossingROC) { | |
279 | if (sign==1) | |
280 | dx[2] += 0.0369; // A side - negative correction | |
281 | else | |
282 | dx[2] -= 0.0369; // C side - positive correction | |
283 | } | |
284 | ||
285 | } | |
286 | ||
c9cbd2f2 | 287 | } |
288 | ||
289 | void AliTPCROCVoltError3D::InitROCVoltError3D() { | |
290 | // | |
291 | // Initialization of the Lookup table which contains the solutions of the | |
292 | // Dirichlet boundary problem | |
293 | // Calculation of the single 3D-Poisson solver is done just if needed | |
294 | // (see basic lookup tables in header file) | |
295 | // | |
296 | ||
297 | const Int_t order = 1 ; // Linear interpolation = 1, Quadratic = 2 | |
298 | const Float_t gridSizeR = (fgkOFCRadius-fgkIFCRadius) / (kRows-1) ; | |
299 | const Float_t gridSizeZ = fgkTPCZ0 / (kColumns-1) ; | |
300 | const Float_t gridSizePhi = TMath::TwoPi() / ( 18.0 * kPhiSlicesPerSector); | |
301 | ||
302 | // temporary arrays to create the boundary conditions | |
303 | TMatrixD *arrayofArrayV[kPhiSlices], *arrayofCharge[kPhiSlices] ; | |
304 | TMatrixD *arrayofEroverEz[kPhiSlices], *arrayofEphioverEz[kPhiSlices], *arrayofDeltaEz[kPhiSlices] ; | |
305 | ||
306 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
307 | arrayofArrayV[k] = new TMatrixD(kRows,kColumns) ; | |
308 | arrayofCharge[k] = new TMatrixD(kRows,kColumns) ; | |
309 | arrayofEroverEz[k] = new TMatrixD(kRows,kColumns) ; | |
310 | arrayofEphioverEz[k] = new TMatrixD(kRows,kColumns) ; | |
311 | arrayofDeltaEz[k] = new TMatrixD(kRows,kColumns) ; | |
312 | } | |
313 | ||
314 | // list of point as used in the poisson relation and the interpolation (during sum up) | |
315 | Double_t rlist[kRows], zedlist[kColumns] , philist[kPhiSlices]; | |
316 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
317 | philist[k] = gridSizePhi * k; | |
318 | for ( Int_t i = 0 ; i < kRows ; i++ ) { | |
319 | rlist[i] = fgkIFCRadius + i*gridSizeR ; | |
320 | for ( Int_t j = 0 ; j < kColumns ; j++ ) { // Fill Vmatrix with Boundary Conditions | |
321 | zedlist[j] = j * gridSizeZ ; | |
322 | } | |
323 | } | |
324 | } | |
325 | ||
326 | // ========================================================================== | |
327 | // Solve Poisson's equation in 3D cylindrical coordinates by relaxation technique | |
328 | // Allow for different size grid spacing in R and Z directions | |
329 | ||
330 | const Int_t symmetry = 0; | |
331 | ||
332 | // Set bondaries and solve Poisson's equation -------------------------- | |
333 | ||
334 | if ( !fInitLookUp ) { | |
335 | ||
336 | AliInfo(Form("Solving the poisson equation (~ %d sec)",2*10*(int)(kPhiSlices/10))); | |
337 | ||
338 | for ( Int_t side = 0 ; side < 2 ; side++ ) { // Solve Poisson3D twice; once for +Z and once for -Z | |
339 | ||
340 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
341 | TMatrixD &arrayV = *arrayofArrayV[k] ; | |
342 | TMatrixD &charge = *arrayofCharge[k] ; | |
343 | ||
344 | //Fill arrays with initial conditions. V on the boundary and Charge in the volume. | |
345 | for ( Int_t i = 0 ; i < kRows ; i++ ) { | |
346 | for ( Int_t j = 0 ; j < kColumns ; j++ ) { // Fill Vmatrix with Boundary Conditions | |
347 | arrayV(i,j) = 0.0 ; | |
348 | charge(i,j) = 0.0 ; | |
349 | ||
350 | Float_t radius0 = rlist[i] ; | |
351 | Float_t phi0 = gridSizePhi * k ; | |
352 | ||
353 | // To avoid problems at sector boundaries, use an average of +- 1 degree from actual phi location | |
35108d57 | 354 | if ( j == (kColumns-1) ) { |
c9cbd2f2 | 355 | arrayV(i,j) = 0.5* ( GetROCVoltOffset( side, radius0, phi0+0.02 ) + GetROCVoltOffset( side, radius0, phi0-0.02 ) ) ; |
356 | ||
35108d57 | 357 | if (side==1) // C side |
358 | arrayV(i,j) = -arrayV(i,j); // minus sign on the C side to allow a consistent usage of global z when setting the boundaries | |
359 | } | |
c9cbd2f2 | 360 | } |
361 | } | |
362 | ||
363 | for ( Int_t i = 1 ; i < kRows-1 ; i++ ) { | |
364 | for ( Int_t j = 1 ; j < kColumns-1 ; j++ ) { | |
365 | charge(i,j) = 0.0 ; | |
366 | } | |
367 | } | |
368 | } | |
369 | ||
370 | // Solve Poisson's equation in 3D cylindrical coordinates by relaxation technique | |
371 | // Allow for different size grid spacing in R and Z directions | |
372 | ||
373 | PoissonRelaxation3D( arrayofArrayV, arrayofCharge, | |
374 | arrayofEroverEz, arrayofEphioverEz, arrayofDeltaEz, | |
375 | kRows, kColumns, kPhiSlices, gridSizePhi, kIterations, | |
376 | symmetry, fROCdisplacement) ; | |
377 | ||
378 | ||
379 | //Interpolate results onto a custom grid which is used just for these calculations. | |
380 | Double_t r, phi, z ; | |
381 | for ( Int_t k = 0 ; k < kNPhi ; k++ ) { | |
382 | phi = fgkPhiList[k] ; | |
383 | ||
2bf29b72 | 384 | TMatrixF &erOverEz = *fLookUpErOverEz[k] ; |
385 | TMatrixF &ephiOverEz = *fLookUpEphiOverEz[k]; | |
386 | TMatrixF &deltaEz = *fLookUpDeltaEz[k] ; | |
c9cbd2f2 | 387 | |
388 | for ( Int_t j = 0 ; j < kNZ ; j++ ) { | |
389 | ||
390 | z = TMath::Abs(fgkZList[j]) ; // Symmetric solution in Z that depends only on ABS(Z) | |
391 | ||
392 | if ( side == 0 && fgkZList[j] < 0 ) continue; // Skip rest of this loop if on the wrong side | |
393 | if ( side == 1 && fgkZList[j] > 0 ) continue; // Skip rest of this loop if on the wrong side | |
394 | ||
395 | for ( Int_t i = 0 ; i < kNR ; i++ ) { | |
396 | r = fgkRList[i] ; | |
397 | ||
398 | // Interpolate basicLookup tables; once for each rod, then sum the results | |
399 | erOverEz(i,j) = Interpolate3DTable(order, r, z, phi, kRows, kColumns, kPhiSlices, | |
400 | rlist, zedlist, philist, arrayofEroverEz ); | |
401 | ephiOverEz(i,j) = Interpolate3DTable(order, r, z, phi, kRows, kColumns, kPhiSlices, | |
402 | rlist, zedlist, philist, arrayofEphioverEz); | |
403 | deltaEz(i,j) = Interpolate3DTable(order, r, z, phi, kRows, kColumns, kPhiSlices, | |
404 | rlist, zedlist, philist, arrayofDeltaEz ); | |
405 | ||
406 | if (side == 1) deltaEz(i,j) = - deltaEz(i,j); // negative coordinate system on C side | |
407 | ||
408 | } // end r loop | |
409 | }// end z loop | |
410 | }// end phi loop | |
411 | ||
412 | if ( side == 0 ) AliInfo(" A side done"); | |
413 | if ( side == 1 ) AliInfo(" C side done"); | |
414 | } // end side loop | |
415 | } | |
416 | ||
417 | // clear the temporary arrays lists | |
418 | for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) { | |
419 | delete arrayofArrayV[k]; | |
420 | delete arrayofCharge[k]; | |
421 | delete arrayofEroverEz[k]; | |
422 | delete arrayofEphioverEz[k]; | |
423 | delete arrayofDeltaEz[k]; | |
424 | } | |
425 | ||
426 | ||
427 | fInitLookUp = kTRUE; | |
428 | ||
429 | } | |
430 | ||
431 | ||
432 | Float_t AliTPCROCVoltError3D::GetROCVoltOffset(Int_t side, Float_t r0, Float_t phi0) { | |
433 | // | |
434 | // Returns the dz alignment data (in voltage equivalents) at | |
435 | // the given position | |
436 | // | |
437 | ||
438 | Float_t xp = r0*TMath::Cos(phi0); | |
439 | Float_t yp = r0*TMath::Sin(phi0); | |
440 | ||
441 | // phi0 should be between 0 and 2pi | |
442 | if (phi0<0) phi0+=TMath::TwoPi(); | |
443 | Int_t roc = (Int_t)TMath::Floor((TMath::RadToDeg()*phi0)/20); | |
444 | if (side==1) roc+=18; // C side | |
445 | if (r0>132) roc+=36; // OROC | |
446 | ||
c756f562 | 447 | // linear-plane data: z = z0 + kx*lx + ky*ly (rotation in local coordinates) |
c9cbd2f2 | 448 | TMatrixD &fitData = *fdzDataLinFit; |
c756f562 | 449 | |
450 | // local coordinates | |
451 | Double_t secAlpha = TMath::DegToRad()*(10.+20.*(((Int_t)roc)%18)); | |
452 | Float_t lx = xp*TMath::Cos(secAlpha)+yp*TMath::Sin(secAlpha); | |
453 | Float_t ly = yp*TMath::Cos(secAlpha)-xp*TMath::Sin(secAlpha); | |
454 | ||
7a348589 | 455 | // reference of rotation in lx is at the intersection between OROC and IROC |
456 | // necessary, since the Fitprozedure is otherwise useless | |
457 | ||
458 | AliTPCROC * rocInfo = AliTPCROC::Instance(); | |
459 | Double_t lxRef = (rocInfo->GetPadRowRadii(0,62)+rocInfo->GetPadRowRadii(36,0))/2; | |
460 | ||
461 | Float_t dz = fitData(roc,0)+fitData(roc,1)*(lx-lxRef) + fitData(roc,2)*ly; // value in cm | |
c9cbd2f2 | 462 | |
463 | // aproximated Voltage-offset-aquivalent to the z misalignment | |
464 | // (linearly scaled with the z position) | |
465 | Double_t ezField = (fgkCathodeV-fgkGG)/fgkTPCZ0; // = ALICE Electric Field (V/cm) Magnitude ~ -400 V/cm; | |
466 | Float_t voltOff = dz*ezField; // values in "Volt equivalents" | |
467 | ||
468 | return voltOff; | |
469 | } | |
470 | ||
7a348589 | 471 | TH2F * AliTPCROCVoltError3D::CreateHistoOfZAlignment(Int_t side, Int_t nx, Int_t ny) { |
c9cbd2f2 | 472 | // |
473 | // return a simple histogramm containing the input to the poisson solver | |
474 | // (z positions of the Read-out chambers, linearly interpolated) | |
475 | ||
476 | char hname[100]; | |
7ee86790 | 477 | if (side==0) snprintf(hname,100,"survey_dz_Aside"); |
478 | if (side==1) snprintf(hname,100,"survey_dz_Cside"); | |
c9cbd2f2 | 479 | |
480 | TH2F *h = new TH2F(hname,hname,nx,-250.,250.,ny,-250.,250.); | |
481 | ||
482 | for (Int_t iy=1;iy<=ny;++iy) { | |
483 | Double_t yp = h->GetYaxis()->GetBinCenter(iy); | |
484 | for (Int_t ix=1;ix<=nx;++ix) { | |
485 | Double_t xp = h->GetXaxis()->GetBinCenter(ix); | |
486 | ||
487 | Float_t r0 = TMath::Sqrt(xp*xp+yp*yp); | |
488 | Float_t phi0 = TMath::ATan2(yp,xp); | |
489 | ||
490 | Float_t dz = GetROCVoltOffset(side,r0,phi0); // in [volt] | |
491 | ||
492 | Double_t ezField = (fgkCathodeV-fgkGG)/fgkTPCZ0; // = ALICE Electric Field (V/cm) Magnitude ~ -400 V/cm; | |
493 | dz = dz/ezField; // in [cm] | |
494 | ||
495 | if (85.<=r0 && r0<=245.) { | |
496 | h->SetBinContent(ix,iy,dz); | |
497 | } else { | |
498 | h->SetBinContent(ix,iy,0.); | |
499 | } | |
500 | } | |
501 | } | |
502 | ||
503 | h->GetXaxis()->SetTitle("x [cm]"); | |
504 | h->GetYaxis()->SetTitle("y [cm]"); | |
505 | h->GetZaxis()->SetTitle("dz [cm]"); | |
506 | h->SetStats(0); | |
507 | // h->DrawCopy("colz"); | |
508 | ||
509 | return h; | |
510 | } | |
511 | ||
512 | void AliTPCROCVoltError3D::Print(const Option_t* option) const { | |
513 | // | |
514 | // Print function to check the settings of the Rod shifts and the rotated clips | |
515 | // option=="a" prints the C0 and C1 coefficents for calibration purposes | |
516 | // | |
517 | ||
518 | TString opt = option; opt.ToLower(); | |
519 | printf("%s\n",GetTitle()); | |
7a348589 | 520 | printf(" - z aligmnet of the TPC Read-Out chambers \n"); |
521 | printf(" (linear interpolation within the chamber: dz = z0 + kx*(lx-133) + ky*ly [cm] ) \n"); | |
522 | printf(" Info: Check the following data-file for more details: %s \n",fROCDataFileName.Data()); | |
c9cbd2f2 | 523 | |
524 | if (opt.Contains("a")) { // Print all details | |
7a348589 | 525 | TMatrixD &fitData = *fdzDataLinFit; |
526 | printf(" A side: IROC ROCX=(z0,kx,ky): \n"); | |
527 | for (Int_t roc = 0; roc<18; roc++) | |
528 | printf("ROC%d:(%.2e,%.2e,%.2e) ",roc,fitData(roc,0),fitData(roc,1),fitData(roc,2)); | |
529 | printf("\n A side: OROC ROCX=(z0,kx,ky): \n"); | |
530 | for (Int_t roc = 36; roc<54; roc++) | |
531 | printf("ROC%d:(%.2e,%.2e,%.2e) ",roc,fitData(roc,0),fitData(roc,1),fitData(roc,2)); | |
532 | printf("\n C side: IROC ROCX=(z0,kx,ky): \n"); | |
533 | for (Int_t roc = 18; roc<36; roc++) | |
534 | printf("ROC%d:(%.2e,%.2e,%.2e) ",roc,fitData(roc,0),fitData(roc,1),fitData(roc,2)); | |
535 | printf("\n C side: OROC ROCX=(z0,kx,ky): \n"); | |
536 | for (Int_t roc = 54; roc<72; roc++) | |
537 | printf("ROC%d:(%.2e,%.2e,%.2e) ",roc,fitData(roc,0),fitData(roc,1),fitData(roc,2)); | |
538 | printf("\n\n"); | |
c9cbd2f2 | 539 | printf(" - T1: %1.4f, T2: %1.4f \n",fT1,fT2); |
540 | printf(" - C1: %1.4f, C0: %1.4f \n",fC1,fC0); | |
541 | } | |
542 | ||
543 | if (!fInitLookUp) AliError("Lookup table was not initialized! You should do InitROCVoltError3D() ..."); | |
544 | ||
545 | } |