]>
Commit | Line | Data |
---|---|---|
b1f0a2a5 | 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 | ||
b1f0a2a5 | 16 | |
0fd1fc0b | 17 | /// \class AliTPCEfield |
18 | /// | |
19 | /// Solution of Laplace equation in cartesian system, with boundary condition. | |
20 | /// | |
21 | /// See details: | |
22 | /// http://web.mit.edu/6.013_book/www/chapter5/5.10.html | |
b1f0a2a5 | 23 | |
24 | #include "TTreeStream.h" | |
25 | #include "TMath.h" | |
26 | #include "TLinearFitter.h" | |
27 | #include "TRandom.h" | |
28 | #include "AliTPCEfield.h" | |
29 | ||
0fd1fc0b | 30 | /// \cond CLASSIMP |
b1f0a2a5 | 31 | ClassImp(AliTPCEfield) |
0fd1fc0b | 32 | /// \endcond |
b1f0a2a5 | 33 | |
34 | AliTPCEfield* AliTPCEfield::fgInstance=0; | |
35 | ||
36 | AliTPCEfield::AliTPCEfield(): | |
37 | TNamed(), | |
38 | fScale(0), | |
39 | fMaxFreq(0), | |
40 | fIs2D(kTRUE), | |
41 | fWorkspace(0) // file with trees, pictures ... | |
42 | { | |
0fd1fc0b | 43 | /// |
44 | ||
b1f0a2a5 | 45 | for (Int_t i=0; i<3; i++){ |
46 | fMin[i]=0; fMax[i]=0; | |
47 | } | |
48 | fgInstance =this; | |
49 | } | |
50 | ||
51 | AliTPCEfield::AliTPCEfield(const char* name, Int_t maxFreq, Bool_t is2D, Bool_t useLinear): | |
52 | TNamed(name,name), | |
53 | fScale(0), | |
54 | fMaxFreq(maxFreq), | |
55 | fIs2D(is2D), | |
56 | fUseLinear(useLinear), | |
57 | fWorkspace(0) // file with trees, pictures ... | |
58 | { | |
0fd1fc0b | 59 | /// |
60 | ||
b1f0a2a5 | 61 | for (Int_t i=0; i<3; i++){ |
62 | fMin[i]=0; fMax[i]=0; | |
63 | } | |
64 | fWorkspace=new TTreeSRedirector(Form("%s.root",name)); | |
65 | MakeFitFunctions(maxFreq); | |
66 | fgInstance =this; | |
67 | } | |
68 | ||
69 | void AliTPCEfield::MakeFitFunctions(Int_t maxFreq){ | |
0fd1fc0b | 70 | /// fit functions = \f$ f(x,y,z) = fx(x)*fy(y)*fz(z) \f$ |
71 | /// function can be of following types: | |
72 | /// 0 - constant | |
73 | /// 1 - linear | |
74 | /// 2 - hypx | |
75 | /// 3 - hypy | |
76 | /// 4 - hypz | |
77 | ||
b1f0a2a5 | 78 | Int_t nfunctions=0; |
79 | if (fIs2D) nfunctions = 1+(maxFreq)*8; | |
80 | if (!fIs2D) nfunctions = 1+(maxFreq)*8; | |
81 | if (fUseLinear) nfunctions+=2; | |
82 | fFitFunctions = new TMatrixD(nfunctions,4); | |
83 | fFitParam = new TVectorD(nfunctions); | |
84 | fFitCovar = new TMatrixD(nfunctions,nfunctions); | |
85 | TMatrixD &fitF = *fFitFunctions; | |
86 | // constant function | |
87 | Int_t counter=1; | |
88 | fitF(0,0)=0; | |
89 | // | |
90 | // linear functions in one cordinate - constan in others | |
91 | if (fUseLinear) | |
92 | for (Int_t ifx=0; ifx<=1; ifx++) | |
93 | for (Int_t ify=0; ify<=1; ify++) | |
94 | for (Int_t ifz=0; ifz<=1; ifz++){ | |
95 | if (fIs2D && ifz>0) continue; | |
96 | if ((ifx+ify+ifz)==0) continue; | |
97 | if ((ifx+ify+ifz)>1) continue; | |
98 | fitF(counter,0)= 1; // function type | |
99 | fitF(counter,1)= ifx; // type of x function | |
100 | fitF(counter,2)= ify; // type of y function | |
101 | fitF(counter,3)= ifz; // type of z function | |
102 | counter++; | |
103 | } | |
104 | // | |
105 | if (fIs2D){ | |
106 | for (Int_t ihyp=0; ihyp<2; ihyp++) | |
107 | for (Int_t ifx=-maxFreq; ifx<=maxFreq; ifx++) | |
108 | for (Int_t ify=-maxFreq; ify<=maxFreq; ify++){ | |
109 | if (TMath::Abs(ify)!=TMath::Abs(ifx)) continue; | |
110 | if (ifx==0) continue; | |
111 | if (ify==0) continue; | |
112 | fitF(counter,0)= 2+ihyp; // function type | |
113 | fitF(counter,1)= ifx; // type of x function - + sinus - cosin | |
114 | fitF(counter,2)= ify; // type of y function | |
115 | fitF(counter,3)= 0; // type of y function | |
116 | counter++; | |
117 | } | |
118 | } | |
119 | ||
120 | } | |
121 | ||
122 | ||
123 | ||
124 | AliTPCEfield::~AliTPCEfield() { | |
0fd1fc0b | 125 | /// Destructor |
126 | ||
b1f0a2a5 | 127 | if (fWorkspace) delete fWorkspace; |
128 | } | |
129 | ||
130 | void AliTPCEfield::SetRange(Double_t x0, Double_t x1, Double_t y0, Double_t y1, Double_t z0,Double_t z1){ | |
0fd1fc0b | 131 | /// Set the ranges - coordinates are rescaled in order to use proper |
132 | /// cos,sin expansion in scaled space | |
133 | ||
b1f0a2a5 | 134 | fMin[0]=x0; fMax[0]=x1; |
135 | fMin[1]=y0; fMax[1]=y1; | |
136 | fMin[2]=z0; fMax[2]=z1; | |
137 | if (fIs2D) fScale=0.5*(TMath::Abs(x1-x0)+TMath::Abs(y1-y0)); | |
138 | if (!fIs2D) fScale=0.5*(TMath::Abs(x1-x0)+TMath::Abs(y1-y0)+TMath::Abs(z1-z0)); | |
139 | } | |
140 | ||
141 | ||
142 | void AliTPCEfield::AddBoundaryLine(Double_t x0,Double_t y0,Double_t z0, Double_t v0, Double_t x1, Double_t y1, Double_t z1,Double_t v1, Int_t id, Int_t npoints){ | |
0fd1fc0b | 143 | /// Add a e field boundary line |
144 | /// From point (x0,y0) to point (x1,y1) | |
145 | /// Linear decrease of potential is assumed | |
146 | /// Boundary can be identified using boundary ID | |
147 | /// The line is written into tree Boundary | |
148 | ||
b1f0a2a5 | 149 | Double_t deltaX = (x1-x0); |
150 | Double_t deltaY = (y1-y0); | |
151 | Double_t deltaZ = (z1-z0); | |
152 | Double_t deltaV = (v1-v0); | |
153 | for (Int_t ipoint=0; ipoint<npoints; ipoint++){ | |
154 | Double_t bpoint=gRandom->Rndm(); | |
155 | Double_t x = x0+deltaX*bpoint; | |
156 | Double_t y = y0+deltaY*bpoint; | |
157 | Double_t z = z0+deltaZ*bpoint; | |
158 | Double_t v = v0+deltaV*bpoint; | |
159 | (*fWorkspace)<<"Boundary"<< | |
160 | "x="<<x<< // x coordinate | |
161 | "y="<<y<< // y coordinate | |
162 | "z="<<z<< // z coordinate | |
163 | "v="<<v<< // potential | |
164 | "id="<<id<< // boundary ID | |
165 | "\n"; | |
166 | } | |
167 | } | |
168 | ||
169 | TTree * AliTPCEfield::GetTree(const char * tname){ | |
0fd1fc0b | 170 | /// |
171 | ||
b1f0a2a5 | 172 | return ((*fWorkspace)<<tname).GetTree(); |
173 | } | |
174 | ||
175 | Double_t AliTPCEfield::Field(Int_t ftype, Double_t ifx, Double_t ify, Double_t ifz, Double_t x, Double_t y, Double_t z){ | |
0fd1fc0b | 176 | /// Field component in |
177 | /// f frequency | |
178 | ||
b1f0a2a5 | 179 | Double_t fx=1,fy=1,fz=1; |
180 | const Double_t kEps=0.01; | |
181 | // | |
182 | if (ftype==0) return 1; | |
183 | if (ftype==1) { | |
184 | if (TMath::Nint(ifx)==1) return x; | |
185 | if (TMath::Nint(ify)==1) return y; | |
186 | if (TMath::Nint(ifz)==1) return z; | |
187 | } | |
188 | Double_t pi = TMath::Pi(); | |
189 | if (ifx>kEps) fx = (ftype==2) ? SinHNorm(ifx*pi*x,TMath::Abs(ifx*pi)): TMath::Sin(ifx*pi*x); | |
190 | if (ifx<-kEps) fx = (ftype==2) ? CosHNorm(ifx*pi*x,TMath::Abs(ifx*pi)): TMath::Cos(ifx*pi*x); | |
191 | // | |
192 | if (ify>kEps) fy = (ftype==3) ? SinHNorm(ify*pi*y,TMath::Abs(ify*pi)): TMath::Sin(ify*pi*y); | |
193 | if (ify<-kEps) fy = (ftype==3) ? CosHNorm(ify*pi*y,TMath::Abs(ify*pi)): TMath::Cos(ify*pi*y); | |
194 | // | |
195 | if (ifz>kEps) fz = (ftype==4) ? SinHNorm(ifz*pi*z,TMath::Abs(ifz*pi)): TMath::Sin(ifz*pi*z); | |
196 | if (ifz<-kEps) fz = (ftype==4) ? CosHNorm(ifz*pi*z,TMath::Abs(ifz*pi)): TMath::Cos(ifz*pi*z); | |
197 | (*fWorkspace)<<"eval"<< | |
198 | "x="<<x<< | |
199 | "y="<<y<< | |
200 | "z="<<z<< | |
201 | "ifx="<<ifx<< | |
202 | "ify="<<ify<< | |
203 | "ifz="<<ifz<< | |
204 | "fx="<<fx<< | |
205 | "fy="<<fy<< | |
206 | "fz="<<fz<< | |
207 | "ftype="<<ftype<< | |
208 | "\n"; | |
209 | return fx*fy*fz; | |
210 | } | |
211 | ||
212 | ||
213 | Double_t AliTPCEfield::FieldDn(Int_t ftype, Double_t ifx, Double_t ify, Double_t ifz, Int_t dn, Double_t x, Double_t y, Double_t z){ | |
0fd1fc0b | 214 | /// Field component in |
215 | /// f frequency | |
216 | ||
b1f0a2a5 | 217 | Double_t fx=1,fy=1,fz=1; |
218 | const Double_t kEps=0.01; | |
219 | // | |
220 | if (ftype==0) return 0.; | |
221 | if (ftype==1) { | |
222 | Double_t value=0; | |
223 | if (TMath::Nint(ifx)==1 &&dn==0) value=1.; | |
224 | if (TMath::Nint(ify)==1 &&dn==1) value=1.; | |
225 | if (TMath::Nint(ifz)==1 &&dn==2) value=1.; | |
226 | return value; | |
227 | } | |
228 | Double_t pi = TMath::Pi(); | |
229 | if (ifx>kEps) fx = (ftype==2) ? SinHNorm(ifx*pi*x,TMath::Abs(ifx*pi)): TMath::Sin(ifx*pi*x); | |
230 | if (ifx<-kEps) fx = (ftype==2) ? CosHNorm(ifx*pi*x,TMath::Abs(ifx*pi)): TMath::Cos(ifx*pi*x); | |
231 | // | |
232 | if (ify>kEps) fy = (ftype==3) ? SinHNorm(ify*pi*y,TMath::Abs(ify*pi)): TMath::Sin(ify*pi*y); | |
233 | if (ify<-kEps) fy = (ftype==3) ? CosHNorm(ify*pi*y,TMath::Abs(ify*pi)): TMath::Cos(ify*pi*y); | |
234 | // | |
235 | if (ifz>kEps) fz = (ftype==4) ? SinHNorm(ifz*pi*z,TMath::Abs(ifz*pi)): TMath::Sin(ifz*pi*z); | |
236 | if (ifz<-kEps) fz = (ftype==4) ? CosHNorm(ifz*pi*z,-TMath::Abs(ifz*pi)): TMath::Cos(ifz*pi*z); | |
237 | ||
238 | if (dn==0){ | |
239 | if (ifx>kEps) fx = (ftype==2) ? CosHNorm(ifx*pi*x,TMath::Abs(ifx*pi)): TMath::Cos(ifx*pi*x); | |
240 | if (ifx<-kEps) fx = (ftype==2) ? SinHNorm(ifx*pi*x,TMath::Abs(ifx*pi)): -TMath::Sin(ifx*pi*x); | |
241 | fx*=ifx*pi; | |
242 | } | |
243 | if (dn==1){ | |
244 | if (ify>kEps) fy = (ftype==3) ? CosHNorm(ify*pi*y,TMath::Abs(ify*pi)): TMath::Cos(ify*pi*y); | |
245 | if (ify<-kEps) fy = (ftype==3) ? SinHNorm(ify*pi*y,TMath::Abs(ify*pi)): -TMath::Sin(ify*pi*y); | |
246 | fy*=ify*pi; | |
247 | } | |
248 | if (dn==2){ | |
249 | if (ifz>kEps) fz = (ftype==4) ? CosHNorm(ifz*pi*z,TMath::Abs(ifz*pi)): TMath::Cos(ifz*pi*z); | |
250 | if (ifz<-kEps) fz = (ftype==4) ? SinHNorm(ifz*pi*z,TMath::Abs(ifz*pi)): -TMath::Sin(ifz*pi*z); | |
251 | fz*=ifz*pi; | |
252 | } | |
253 | ||
254 | return fx*fy*fz; | |
255 | } | |
256 | ||
257 | ||
258 | ||
259 | ||
260 | Double_t AliTPCEfield::EvalField(Int_t ifun, Double_t x, Double_t y, Double_t z, Int_t type){ | |
0fd1fc0b | 261 | /// Evaluate function ifun at position gx amd gy |
262 | /// type == 0 - field | |
263 | /// == 1 - Ex | |
264 | /// == 2 - Ey | |
265 | /// == 3 - Ez | |
266 | ||
b1f0a2a5 | 267 | TMatrixD &mat = *fFitFunctions; |
268 | Int_t fid = TMath::Nint(mat(ifun,0)); | |
269 | Double_t ifx = (mat(ifun,1)); | |
270 | Double_t ify = (mat(ifun,2)); | |
271 | Double_t ifz = (mat(ifun,3)); | |
272 | // | |
273 | if (type==0) return Field(fid,ifx,ify,ifz, x, y,z); | |
274 | if (type>0) return FieldDn(fid,ifx,ify,ifz,type-1, x, y,z); | |
275 | return 0; | |
276 | } | |
277 | ||
278 | Double_t AliTPCEfield::Eval(Double_t x, Double_t y, Double_t z, Int_t type){ | |
0fd1fc0b | 279 | /// Evaluate function ifun at position gx amd gy |
280 | /// type == 0 - field | |
281 | /// == 1 - Ex | |
282 | /// == 2 - Ey | |
283 | /// == 3 - Ez | |
284 | ||
b1f0a2a5 | 285 | Double_t value=0; |
286 | Double_t lx= 2.*(x-(fMin[0]+fMax[0])*0.5)/fScale; | |
287 | Double_t ly= 2.*(y-(fMin[1]+fMax[1])*0.5)/fScale; | |
288 | Double_t lz= 2.*(z-(fMin[2]+fMax[2])*0.5)/fScale; | |
289 | // | |
290 | Int_t nfun=fFitFunctions->GetNrows(); | |
291 | for (Int_t ifun=0; ifun<nfun; ifun++){ | |
292 | if (type==0) value+=(*fFitParam)[ifun]*EvalField(ifun,lx,ly,lz,type); | |
293 | if (type>0) value+=2*(*fFitParam)[ifun]*EvalField(ifun,lx,ly,lz,type)/fScale; | |
294 | } | |
295 | return value; | |
296 | } | |
297 | ||
298 | Double_t AliTPCEfield::EvalS(Double_t x, Double_t y, Double_t z, Int_t type){ | |
0fd1fc0b | 299 | /// static evaluation - possible to use it in the TF1 |
300 | ||
b1f0a2a5 | 301 | return fgInstance->Eval(x,y,z,type); |
302 | } | |
303 | ||
304 | void AliTPCEfield::FitField(){ | |
0fd1fc0b | 305 | /// Fit the e field |
306 | /// Minimize chi2 residuals at the boundary points | |
307 | /// ?Tempoary sollution - integrals can be calculated analytically - | |
308 | ||
b1f0a2a5 | 309 | Int_t nfun=fFitFunctions->GetNrows(); |
310 | Double_t *fun =new Double_t[nfun]; | |
311 | fFitter= new TLinearFitter(nfun, Form("hyp%d", nfun-1)); | |
312 | // | |
313 | TTree * tree = GetTree("Boundary"); | |
314 | Int_t npoints = tree->GetEntries(); | |
315 | Int_t *indexes = new Int_t[npoints]; | |
316 | Double_t *rindex = new Double_t[npoints]; | |
317 | // | |
318 | ||
319 | Double_t x=0, y=0, z=0, v=0; | |
320 | tree->SetBranchAddress("x",&x); | |
321 | tree->SetBranchAddress("y",&y); | |
322 | tree->SetBranchAddress("z",&z); | |
323 | tree->SetBranchAddress("v",&v); | |
324 | TMatrixD valMatrix(npoints,4); | |
325 | for (Int_t ipoint=0; ipoint<npoints; ipoint++){ | |
326 | tree->GetEntry(ipoint); | |
327 | valMatrix(ipoint,0)=2.*(x-(fMin[0]+fMax[0])*0.5)/fScale; | |
328 | valMatrix(ipoint,1)=2.*(y-(fMin[1]+fMax[1])*0.5)/fScale; | |
329 | valMatrix(ipoint,2)=2.*(z-(fMin[2]+fMax[2])*0.5)/fScale; | |
330 | valMatrix(ipoint,3)=v; | |
331 | rindex[ipoint]=gRandom->Rndm(); | |
332 | } | |
333 | TMath::Sort(npoints,rindex,indexes); | |
334 | for (Int_t jpoint=0; jpoint<npoints; jpoint++){ | |
335 | Int_t ipoint=indexes[jpoint]; | |
336 | // | |
337 | Double_t lx= valMatrix(ipoint,0); | |
338 | Double_t ly= valMatrix(ipoint,1); | |
339 | Double_t lz= valMatrix(ipoint,2); | |
340 | Double_t value = valMatrix(ipoint,3); | |
341 | for (Int_t ifun=1; ifun<nfun; ifun++){ | |
342 | Double_t ffun=EvalField(ifun,lx,ly,lz,0); | |
343 | fun[ifun-1]=ffun; | |
344 | } | |
345 | fFitter->AddPoint(fun,value,1); | |
346 | } | |
347 | fFitter->Eval(); | |
348 | fFitter->GetCovarianceMatrix(*fFitCovar); | |
349 | fFitter->GetParameters(*fFitParam); | |
350 | ||
351 | (*fWorkspace)<<"fitGlobal"<< | |
352 | "covar.="<<fFitCovar<< // covariance matrix | |
353 | "param.="<<fFitParam<< // fit parameters | |
354 | "fun.="<<fFitFunctions<< // type of fit functions | |
355 | "\n"; | |
356 | ||
357 | for (Int_t ifun=1; ifun<nfun; ifun++){ | |
358 | ||
359 | } | |
360 | ||
361 | ||
362 | } | |
363 | ||
364 | ||
365 | TMatrixD* AliTPCEfield::MakeCorrelation(TMatrixD &matrix){ | |
0fd1fc0b | 366 | /// |
367 | ||
b1f0a2a5 | 368 | Int_t nrows = matrix.GetNrows(); |
369 | TMatrixD * mat = new TMatrixD(nrows,nrows); | |
370 | for (Int_t irow=0; irow<nrows; irow++) | |
371 | for (Int_t icol=0; icol<nrows; icol++){ | |
372 | (*mat)(irow,icol)= matrix(irow,icol)/TMath::Sqrt(matrix(irow,irow)*matrix(icol,icol)); | |
373 | } | |
374 | return mat; | |
375 | } | |
376 | ||
377 | ||
378 | ||
379 | ||
380 | void AliTPCEfield::DumpField(Double_t gridSize, Double_t step){ | |
0fd1fc0b | 381 | /// |
382 | ||
b1f0a2a5 | 383 | Double_t stepSize=0.001*fScale/fMaxFreq; |
384 | // | |
385 | for (Double_t x = fMin[0]+stepSize; x<=fMax[0]-stepSize; x+=gridSize){ | |
386 | for (Double_t y = fMin[1]+stepSize; y<=fMax[1]-stepSize; y+=gridSize) | |
387 | for (Double_t z = fMin[2]; z<=fMax[2]; z+=gridSize){ | |
388 | // | |
389 | // | |
390 | Double_t v = Eval(x,y,z,0); | |
391 | Double_t ex = Eval(x,y,z,1); | |
392 | Double_t ey = Eval(x,y,z,2); | |
393 | Double_t ez = Eval(x,y,z,3); | |
394 | Double_t dexdx = (Eval(x,y,z,1)-Eval(x-stepSize,y,z,1))/stepSize; // numerical derivative | |
395 | Double_t deydy = (Eval(x,y,z,2)-Eval(x,y-stepSize,z,2))/stepSize; | |
396 | Double_t dezdz = (Eval(x,y,z,3)-Eval(x,y,z-stepSize,3))/stepSize; | |
397 | (*fWorkspace)<<"dumpField"<< | |
398 | "x="<<x<< // position | |
399 | "y="<<y<< | |
400 | "z="<<z<< | |
401 | "v="<<v<< // potential | |
402 | "ex="<<ex<< // Efield | |
403 | "ey="<<ey<< | |
404 | "ez="<<ez<< | |
405 | "dexdx="<<dexdx<< //gradient of e field | |
406 | "deydy="<<deydy<< | |
407 | "dezdz="<<dezdz<< | |
408 | "\n"; | |
409 | } | |
410 | } | |
411 | // | |
412 | // | |
413 | for (Double_t x = fMin[0]+stepSize; x<=fMax[0]-stepSize; x+=gridSize){ | |
414 | Double_t sumEx=0; | |
415 | Double_t sumEy=0; | |
416 | Double_t sumExEy=0; | |
417 | for (Double_t y = fMin[1]+0.2; y<=fMax[1]-stepSize; y+=step) | |
418 | for (Double_t z = fMin[2]; z<=fMax[2]; z+=gridSize){ | |
419 | // | |
420 | // | |
421 | Double_t v = Eval(x,y+step*0.5,z,0); | |
422 | Double_t ex = Eval(x,y+step*0.5,z,1); | |
423 | Double_t ey = Eval(x,y+step*0.5,z,2); | |
424 | Double_t ez = Eval(x,y+step*0.5,z,3); | |
425 | sumEx+=ex*step; | |
426 | sumEy+=ey*step; | |
427 | sumExEy+=(ex/ey)*step; | |
428 | (*fWorkspace)<<"dumpDistortion"<< | |
429 | "x="<<x<< // position | |
430 | "y="<<y<< | |
431 | "z="<<z<< | |
432 | "v="<<v<< // potential | |
433 | "ex="<<ex<< // Efield | |
434 | "ey="<<ey<< | |
435 | "ez="<<ez<< | |
436 | // | |
437 | "sEx="<<sumEx<< // field x integral | |
438 | "sEy="<<sumEy<< // field y integral | |
439 | "sExEy="<<sumExEy<< // tan integral | |
440 | "\n"; | |
441 | } | |
442 | } | |
443 | } | |
444 | ||
445 | ||
446 | void MakeTPC2DExample(AliTPCEfield *field){ | |
447 | // | |
448 | /* | |
449 | .L $ALICE_ROOT/TPC/AliTPCEfield.cxx++ | |
450 | AliTPCEfield *field = new AliTPCEfield("field",20, kTRUE,kTRUE); | |
451 | MakeTPC2DExample(field) | |
452 | field->FitField() | |
453 | sqrt(field->fFitter.GetChisquare()/field->fFitter.GetNpoints()) | |
454 | ||
455 | TF2 f2("f2","AliTPCEfield::EvalS(x,y,0,0)",90,245,0,250); | |
456 | f2->SetNpx(100); f2->SetNpy(100); f2->Draw("colz"); | |
457 | ||
458 | TF2 f2x("f2x","AliTPCEfield::EvalS(x,y,0,1)",90,240,0,240); | |
459 | f2x->SetNpx(100); f2x->SetNpy(100); f2x->Draw("surf2"); | |
460 | ||
461 | TF2 f2y("f2y","AliTPCEfield::EvalS(x,y,0,2)",90,240,0,240); | |
462 | f2y->SetNpx(100); f2y->SetNpy(100); f2y->Draw("surf2"); | |
463 | ||
464 | field->MakeCorrelation(*(field->fFitCovar)).Print() | |
465 | Double_t index[100000]; | |
466 | for (Int_t i=0; i<field->fFitCovar->GetNrows(); i++) index[i]=i; | |
467 | TGraph gr(field->fFitCovar->GetNrows(), index, field->fFitParam->GetMatrixArray()); | |
468 | gr->Draw("alp"); | |
469 | ||
470 | field->GetTree()->Draw("AliTPCEfield::EvalS(x,y,0,0):v"); | |
471 | ||
472 | ||
473 | TF2 f2xdy("f2xdy","AliTPCEfield::EvalS(x,y,0,1)/AliTPCEfield::EvalS(x,y,0,2)",90,240,0,240); | |
474 | f2xdy->SetNpx(100); f2xdy->SetNpy(100); f2xdy->Draw("colz"); | |
475 | ||
476 | */ | |
477 | ||
478 | Double_t p0[4]; | |
479 | Double_t p1[4]; | |
480 | Double_t xmin=85, xmax=245; | |
481 | Double_t ymin=0, ymax=250, deltaY=0.1*ymax; | |
482 | Double_t vup=1; | |
483 | Int_t npoints=1000; | |
484 | field->SetRange(xmin, xmax,ymin,ymax,0,0); | |
485 | // upper part | |
486 | p0[0]=xmin; p0[1]=ymax+deltaY; p0[2]=0; p0[3]=vup; | |
487 | p1[0]=xmax; p1[1]=ymax-deltaY; p1[2]=0; p1[3]=vup; | |
488 | field->AddBoundaryLine(p0[0],p0[1], p0[2],p0[3],p1[0],p1[1], p1[2],p1[3],1,npoints); | |
489 | //left | |
490 | p0[0]=xmin; p0[1]=ymin+deltaY; p0[2]=0; p0[3]=0; | |
491 | p1[0]=xmin; p1[1]=ymax+deltaY; p1[2]=0; p1[3]=vup; | |
492 | field->AddBoundaryLine(p0[0],p0[1], p0[2],p0[3],p1[0],p1[1], p1[2],p1[3],2,npoints); | |
493 | //right | |
494 | p0[0]=xmax; p0[1]=ymin-deltaY; p0[2]=0; p0[3]=0; | |
495 | p1[0]=xmax; p1[1]=ymax-deltaY; p1[2]=0; p1[3]=vup; | |
496 | field->AddBoundaryLine(p0[0],p0[1], p0[2],p0[3],p1[0],p1[1], p1[2],p1[3],3,npoints); | |
497 | // | |
498 | //ROC | |
499 | p0[0]=xmin; p0[1]=deltaY; p0[2]=0; p0[3]=-0; | |
500 | p1[0]=xmax; p1[1]=-deltaY; p1[2]=0; p1[3]=0; | |
501 | field->AddBoundaryLine(p0[0],p0[1], p0[2],p0[3],p1[0],p1[1], p1[2],p1[3],4,npoints); | |
502 | } | |
503 | ||
504 | ||
505 |