]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TPC/AliTPCClusterParam.cxx
Add setter for cluster unfolding (Marian)
[u/mrichter/AliRoot.git] / TPC / AliTPCClusterParam.cxx
CommitLineData
12ca5da1 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// //
19// TPC cluster error and shape parameterization //
20// //
21// //
22///////////////////////////////////////////////////////////////////////////////
23#include "AliTPCClusterParam.h"
24#include "TMath.h"
25#include "TFile.h"
26#include "TTree.h"
27#include <TVectorF.h>
28#include <TLinearFitter.h>
29#include <TH1F.h>
30#include <TProfile2D.h>
0a65832b 31#include <TVectorD.h>
32#include <TObjArray.h>
12ca5da1 33
34ClassImp(AliTPCClusterParam)
35
36
37AliTPCClusterParam* AliTPCClusterParam::fgInstance = 0;
38
39
40/*
41 Example usage fitting parameterization:
42 TFile fres("resol.root"); //tree with resolution and shape
43 TTree * treeRes =(TTree*)fres.Get("Resol");
44
45 AliTPCClusterParam param;
46 param.SetInstance(&param);
47 param.FitResol(treeRes);
48 param.FitRMS(treeRes);
49 TFile fparam("TPCClusterParam.root","recreate");
50 param.Write("Param");
51 //
52 //
53 TFile fparam("TPCClusterParam.root");
54 AliTPCClusterParam *param2 = (AliTPCClusterParam *) fparam.Get("Param");
55 param2->SetInstance(param2);
56 param2->Test(treeRes);
57
58
59 treeRes->Draw("(Resol-AliTPCClusterParam::SGetError0(Dim,Pad,Zm,AngleM))/Resol","Dim==0&&QMean<0")
60
61*/
62
63
64
65
66//_ singleton implementation __________________________________________________
67AliTPCClusterParam* AliTPCClusterParam::Instance()
68{
69 //
70 // Singleton implementation
71 // Returns an instance of this class, it is created if neccessary
72 //
73 if (fgInstance == 0){
74 fgInstance = new AliTPCClusterParam();
75 }
76 return fgInstance;
77}
78
79
80
81
82void AliTPCClusterParam::FitResol0(TTree * tree, Int_t dim, Int_t type, Float_t *param0, Float_t *error){
83 //
84 // Fit z - angular dependence of resolution
85 //
86 // Int_t dim=0, type=0;
87 char varVal[100];
88 sprintf(varVal,"Resol:AngleM:Zm");
89 char varErr[100];
90 sprintf(varErr,"Sigma:AngleS:Zs");
91 char varCut[100];
92 sprintf(varCut,"Dim==%d&&Pad==%d&&QMean<0",dim,type);
93 //
94 Int_t entries = tree->Draw(varVal,varCut);
95 Float_t px[10000], py[10000], pz[10000];
96 Float_t ex[10000], ey[10000], ez[10000];
97 //
98 tree->Draw(varErr,varCut);
99 for (Int_t ipoint=0; ipoint<entries; ipoint++){
100 ex[ipoint]= tree->GetV3()[ipoint];
101 ey[ipoint]= tree->GetV2()[ipoint];
102 ez[ipoint]= tree->GetV1()[ipoint];
103 }
104 tree->Draw(varVal,varCut);
105 for (Int_t ipoint=0; ipoint<entries; ipoint++){
106 px[ipoint]= tree->GetV3()[ipoint];
107 py[ipoint]= tree->GetV2()[ipoint];
108 pz[ipoint]= tree->GetV1()[ipoint];
109 }
110
111 //
112 TLinearFitter fitter(3,"hyp2");
113 for (Int_t ipoint=0; ipoint<entries; ipoint++){
114 Float_t val = pz[ipoint]*pz[ipoint];
115 Float_t err = 2*pz[ipoint]*TMath::Sqrt(ez[ipoint]*ez[ipoint]+fRatio*fRatio*pz[ipoint]*pz[ipoint]);
116 Double_t x[2];
117 x[0] = px[ipoint];
118 x[1] = py[ipoint]*py[ipoint];
119 fitter.AddPoint(x,val,err);
120 }
121 fitter.Eval();
122 TVectorD param(3);
123 fitter.GetParameters(param);
124 param0[0] = param[0];
125 param0[1] = param[1];
126 param0[2] = param[2];
127 Float_t chi2 = fitter.GetChisquare()/entries;
128 param0[3] = chi2;
129 error[0] = (fitter.GetParError(0)*TMath::Sqrt(chi2));
130 error[1] = (fitter.GetParError(1)*TMath::Sqrt(chi2));
131 error[2] = (fitter.GetParError(2)*TMath::Sqrt(chi2));
132}
133
134
135void AliTPCClusterParam::FitResol0Par(TTree * tree, Int_t dim, Int_t type, Float_t *param0, Float_t *error){
136 //
137 // Fit z - angular dependence of resolution
138 //
139 // Int_t dim=0, type=0;
140 char varVal[100];
141 sprintf(varVal,"Resol:AngleM:Zm");
142 char varErr[100];
143 sprintf(varErr,"Sigma:AngleS:Zs");
144 char varCut[100];
145 sprintf(varCut,"Dim==%d&&Pad==%d&&QMean<0",dim,type);
146 //
147 Int_t entries = tree->Draw(varVal,varCut);
148 Float_t px[10000], py[10000], pz[10000];
149 Float_t ex[10000], ey[10000], ez[10000];
150 //
151 tree->Draw(varErr,varCut);
152 for (Int_t ipoint=0; ipoint<entries; ipoint++){
153 ex[ipoint]= tree->GetV3()[ipoint];
154 ey[ipoint]= tree->GetV2()[ipoint];
155 ez[ipoint]= tree->GetV1()[ipoint];
156 }
157 tree->Draw(varVal,varCut);
158 for (Int_t ipoint=0; ipoint<entries; ipoint++){
159 px[ipoint]= tree->GetV3()[ipoint];
160 py[ipoint]= tree->GetV2()[ipoint];
161 pz[ipoint]= tree->GetV1()[ipoint];
162 }
163
164 //
165 TLinearFitter fitter(6,"hyp5");
166 for (Int_t ipoint=0; ipoint<entries; ipoint++){
167 Float_t val = pz[ipoint]*pz[ipoint];
168 Float_t err = 2*pz[ipoint]*TMath::Sqrt(ez[ipoint]*ez[ipoint]+fRatio*fRatio*pz[ipoint]*pz[ipoint]);
169 Double_t x[6];
170 x[0] = px[ipoint];
171 x[1] = py[ipoint]*py[ipoint];
172 x[2] = x[0]*x[0];
173 x[3] = x[1]*x[1];
174 x[4] = x[0]*x[1];
175 fitter.AddPoint(x,val,err);
176 }
177 fitter.Eval();
178 TVectorD param(6);
179 fitter.GetParameters(param);
180 param0[0] = param[0];
181 param0[1] = param[1];
182 param0[2] = param[2];
183 param0[3] = param[3];
184 param0[4] = param[4];
185 param0[5] = param[5];
186 Float_t chi2 = fitter.GetChisquare()/entries;
187 param0[6] = chi2;
188 error[0] = (fitter.GetParError(0)*TMath::Sqrt(chi2));
189 error[1] = (fitter.GetParError(1)*TMath::Sqrt(chi2));
190 error[2] = (fitter.GetParError(2)*TMath::Sqrt(chi2));
191 error[3] = (fitter.GetParError(3)*TMath::Sqrt(chi2));
192 error[4] = (fitter.GetParError(4)*TMath::Sqrt(chi2));
193 error[5] = (fitter.GetParError(5)*TMath::Sqrt(chi2));
194}
195
196
197
198
199
200void AliTPCClusterParam::FitResol1(TTree * tree, Int_t dim, Float_t *param0, Float_t *error){
201 //
202 // Fit z - angular dependence of resolution - pad length scaling
203 //
204 // Int_t dim=0, type=0;
205 char varVal[100];
206 sprintf(varVal,"Resol:AngleM*sqrt(Length):Zm/Length");
207 char varErr[100];
208 sprintf(varErr,"Sigma:AngleS:Zs");
209 char varCut[100];
210 sprintf(varCut,"Dim==%d&&QMean<0",dim);
211 //
212 Int_t entries = tree->Draw(varVal,varCut);
213 Float_t px[10000], py[10000], pz[10000];
214 Float_t ex[10000], ey[10000], ez[10000];
215 //
216 tree->Draw(varErr,varCut);
217 for (Int_t ipoint=0; ipoint<entries; ipoint++){
218 ex[ipoint]= tree->GetV3()[ipoint];
219 ey[ipoint]= tree->GetV2()[ipoint];
220 ez[ipoint]= tree->GetV1()[ipoint];
221 }
222 tree->Draw(varVal,varCut);
223 for (Int_t ipoint=0; ipoint<entries; ipoint++){
224 px[ipoint]= tree->GetV3()[ipoint];
225 py[ipoint]= tree->GetV2()[ipoint];
226 pz[ipoint]= tree->GetV1()[ipoint];
227 }
228
229 //
230 TLinearFitter fitter(3,"hyp2");
231 for (Int_t ipoint=0; ipoint<entries; ipoint++){
232 Float_t val = pz[ipoint]*pz[ipoint];
233 Float_t err = 2*pz[ipoint]*TMath::Sqrt(ez[ipoint]*ez[ipoint]+fRatio*fRatio*pz[ipoint]*pz[ipoint]);
234 Double_t x[2];
235 x[0] = px[ipoint];
236 x[1] = py[ipoint]*py[ipoint];
237 fitter.AddPoint(x,val,err);
238 }
239 fitter.Eval();
240 TVectorD param(3);
241 fitter.GetParameters(param);
242 param0[0] = param[0];
243 param0[1] = param[1];
244 param0[2] = param[2];
245 Float_t chi2 = fitter.GetChisquare()/entries;
246 param0[3] = chi2;
247 error[0] = (fitter.GetParError(0)*TMath::Sqrt(chi2));
248 error[1] = (fitter.GetParError(1)*TMath::Sqrt(chi2));
249 error[2] = (fitter.GetParError(2)*TMath::Sqrt(chi2));
250}
251
252void AliTPCClusterParam::FitResolQ(TTree * tree, Int_t dim, Int_t type, Float_t *param0, Float_t *error){
253 //
254 // Fit z - angular dependence of resolution - Q scaling
255 //
256 // Int_t dim=0, type=0;
257 char varVal[100];
258 sprintf(varVal,"Resol:AngleM/sqrt(QMean):Zm/QMean");
259 char varVal0[100];
260 sprintf(varVal0,"Resol:AngleM:Zm");
261 //
262 char varErr[100];
263 sprintf(varErr,"Sigma:AngleS:Zs");
264 char varCut[100];
265 sprintf(varCut,"Dim==%d&&Pad==%d&&QMean>0",dim,type);
266 //
267 Int_t entries = tree->Draw(varVal,varCut);
268 Float_t px[20000], py[20000], pz[20000], pu[20000], pt[20000];
269 Float_t ex[20000], ey[20000], ez[20000];
270 //
271 tree->Draw(varErr,varCut);
272 for (Int_t ipoint=0; ipoint<entries; ipoint++){
273 ex[ipoint]= tree->GetV3()[ipoint];
274 ey[ipoint]= tree->GetV2()[ipoint];
275 ez[ipoint]= tree->GetV1()[ipoint];
276 }
277 tree->Draw(varVal,varCut);
278 for (Int_t ipoint=0; ipoint<entries; ipoint++){
279 px[ipoint]= tree->GetV3()[ipoint];
280 py[ipoint]= tree->GetV2()[ipoint];
281 pz[ipoint]= tree->GetV1()[ipoint];
282 }
283 tree->Draw(varVal0,varCut);
284 for (Int_t ipoint=0; ipoint<entries; ipoint++){
285 pu[ipoint]= tree->GetV3()[ipoint];
286 pt[ipoint]= tree->GetV2()[ipoint];
287 }
288
289 //
290 TLinearFitter fitter(5,"hyp4");
291 for (Int_t ipoint=0; ipoint<entries; ipoint++){
292 Float_t val = pz[ipoint]*pz[ipoint];
293 Float_t err = 2*pz[ipoint]*TMath::Sqrt(ez[ipoint]*ez[ipoint]+fRatio*fRatio*pz[ipoint]*pz[ipoint]);
294 Double_t x[4];
295 x[0] = pu[ipoint];
296 x[1] = pt[ipoint]*pt[ipoint];
297 x[2] = px[ipoint];
298 x[3] = py[ipoint]*py[ipoint];
299 fitter.AddPoint(x,val,err);
300 }
301
302 fitter.Eval();
303 TVectorD param(5);
304 fitter.GetParameters(param);
305 param0[0] = param[0];
306 param0[1] = param[1];
307 param0[2] = param[2];
308 param0[3] = param[3];
309 param0[4] = param[4];
310 Float_t chi2 = fitter.GetChisquare()/entries;
311 param0[5] = chi2;
312 error[0] = (fitter.GetParError(0)*TMath::Sqrt(chi2));
313 error[1] = (fitter.GetParError(1)*TMath::Sqrt(chi2));
314 error[2] = (fitter.GetParError(2)*TMath::Sqrt(chi2));
315 error[3] = (fitter.GetParError(3)*TMath::Sqrt(chi2));
316 error[4] = (fitter.GetParError(4)*TMath::Sqrt(chi2));
317}
318
319void AliTPCClusterParam::FitResolQPar(TTree * tree, Int_t dim, Int_t type, Float_t *param0, Float_t *error){
320 //
321 // Fit z - angular dependence of resolution - Q scaling - parabolic correction
322 //
323 // Int_t dim=0, type=0;
324 char varVal[100];
325 sprintf(varVal,"Resol:AngleM/sqrt(QMean):Zm/QMean");
326 char varVal0[100];
327 sprintf(varVal0,"Resol:AngleM:Zm");
328 //
329 char varErr[100];
330 sprintf(varErr,"Sigma:AngleS:Zs");
331 char varCut[100];
332 sprintf(varCut,"Dim==%d&&Pad==%d&&QMean>0",dim,type);
333 //
334 Int_t entries = tree->Draw(varVal,varCut);
335 Float_t px[20000], py[20000], pz[20000], pu[20000], pt[20000];
336 Float_t ex[20000], ey[20000], ez[20000];
337 //
338 tree->Draw(varErr,varCut);
339 for (Int_t ipoint=0; ipoint<entries; ipoint++){
340 ex[ipoint]= tree->GetV3()[ipoint];
341 ey[ipoint]= tree->GetV2()[ipoint];
342 ez[ipoint]= tree->GetV1()[ipoint];
343 }
344 tree->Draw(varVal,varCut);
345 for (Int_t ipoint=0; ipoint<entries; ipoint++){
346 px[ipoint]= tree->GetV3()[ipoint];
347 py[ipoint]= tree->GetV2()[ipoint];
348 pz[ipoint]= tree->GetV1()[ipoint];
349 }
350 tree->Draw(varVal0,varCut);
351 for (Int_t ipoint=0; ipoint<entries; ipoint++){
352 pu[ipoint]= tree->GetV3()[ipoint];
353 pt[ipoint]= tree->GetV2()[ipoint];
354 }
355
356 //
357 TLinearFitter fitter(8,"hyp7");
358 for (Int_t ipoint=0; ipoint<entries; ipoint++){
359 Float_t val = pz[ipoint]*pz[ipoint];
360 Float_t err = 2*pz[ipoint]*TMath::Sqrt(ez[ipoint]*ez[ipoint]+fRatio*fRatio*pz[ipoint]*pz[ipoint]);
361 Double_t x[7];
362 x[0] = pu[ipoint];
363 x[1] = pt[ipoint]*pt[ipoint];
364 x[2] = x[0]*x[0];
365 x[3] = x[1]*x[1];
366 x[4] = x[0]*x[1];
367 x[5] = px[ipoint];
368 x[6] = py[ipoint]*py[ipoint];
369 //
370 fitter.AddPoint(x,val,err);
371 }
372
373 fitter.Eval();
374 TVectorD param(8);
375 fitter.GetParameters(param);
376 param0[0] = param[0];
377 param0[1] = param[1];
378 param0[2] = param[2];
379 param0[3] = param[3];
380 param0[4] = param[4];
381 param0[5] = param[5];
382 param0[6] = param[6];
383 param0[7] = param[7];
384
385 Float_t chi2 = fitter.GetChisquare()/entries;
386 param0[8] = chi2;
387 error[0] = (fitter.GetParError(0)*TMath::Sqrt(chi2));
388 error[1] = (fitter.GetParError(1)*TMath::Sqrt(chi2));
389 error[2] = (fitter.GetParError(2)*TMath::Sqrt(chi2));
390 error[3] = (fitter.GetParError(3)*TMath::Sqrt(chi2));
391 error[4] = (fitter.GetParError(4)*TMath::Sqrt(chi2));
392 error[5] = (fitter.GetParError(5)*TMath::Sqrt(chi2));
393 error[6] = (fitter.GetParError(6)*TMath::Sqrt(chi2));
394 error[7] = (fitter.GetParError(7)*TMath::Sqrt(chi2));
395}
396
397
398
399void AliTPCClusterParam::FitRMS0(TTree * tree, Int_t dim, Int_t type, Float_t *param0, Float_t *error){
400 //
401 // Fit z - angular dependence of resolution
402 //
403 // Int_t dim=0, type=0;
404 char varVal[100];
405 sprintf(varVal,"RMSm:AngleM:Zm");
406 char varErr[100];
407 sprintf(varErr,"sqrt((1./(100.*sqrt(12.))^2)+RMSe0^2):AngleS:Zs");
408 char varCut[100];
409 sprintf(varCut,"Dim==%d&&Pad==%d&&QMean<0",dim,type);
410 //
411 Int_t entries = tree->Draw(varVal,varCut);
412 Float_t px[10000], py[10000], pz[10000];
413 Float_t ex[10000], ey[10000], ez[10000];
414 //
415 tree->Draw(varErr,varCut);
416 for (Int_t ipoint=0; ipoint<entries; ipoint++){
417 ex[ipoint]= tree->GetV3()[ipoint];
418 ey[ipoint]= tree->GetV2()[ipoint];
419 ez[ipoint]= tree->GetV1()[ipoint];
420 }
421 tree->Draw(varVal,varCut);
422 for (Int_t ipoint=0; ipoint<entries; ipoint++){
423 px[ipoint]= tree->GetV3()[ipoint];
424 py[ipoint]= tree->GetV2()[ipoint];
425 pz[ipoint]= tree->GetV1()[ipoint];
426 }
427
428 //
429 TLinearFitter fitter(3,"hyp2");
430 for (Int_t ipoint=0; ipoint<entries; ipoint++){
431 Float_t val = pz[ipoint]*pz[ipoint];
432 Float_t err = 2*pz[ipoint]*TMath::Sqrt(ez[ipoint]*ez[ipoint]+fRatio*fRatio*pz[ipoint]*pz[ipoint]);
433 Double_t x[2];
434 x[0] = px[ipoint];
435 x[1] = py[ipoint]*py[ipoint];
436 fitter.AddPoint(x,val,err);
437 }
438 fitter.Eval();
439 TVectorD param(3);
440 fitter.GetParameters(param);
441 param0[0] = param[0];
442 param0[1] = param[1];
443 param0[2] = param[2];
444 Float_t chi2 = fitter.GetChisquare()/entries;
445 param0[3] = chi2;
446 error[0] = (fitter.GetParError(0)*TMath::Sqrt(chi2));
447 error[1] = (fitter.GetParError(1)*TMath::Sqrt(chi2));
448 error[2] = (fitter.GetParError(2)*TMath::Sqrt(chi2));
449}
450
451void AliTPCClusterParam::FitRMS1(TTree * tree, Int_t dim, Float_t *param0, Float_t *error){
452 //
453 // Fit z - angular dependence of resolution - pad length scaling
454 //
455 // Int_t dim=0, type=0;
456 char varVal[100];
457 sprintf(varVal,"RMSm:AngleM*Length:Zm");
458 char varErr[100];
459 sprintf(varErr,"sqrt((1./(100.*sqrt(12.))^2)+RMSe0^2):AngleS:Pad");
460 char varCut[100];
461 sprintf(varCut,"Dim==%d&&QMean<0",dim);
462 //
463 Int_t entries = tree->Draw(varVal,varCut);
464 Float_t px[10000], py[10000], pz[10000];
465 Float_t type[10000], ey[10000], ez[10000];
466 //
467 tree->Draw(varErr,varCut);
468 for (Int_t ipoint=0; ipoint<entries; ipoint++){
469 type[ipoint] = tree->GetV3()[ipoint];
470 ey[ipoint] = tree->GetV2()[ipoint];
471 ez[ipoint] = tree->GetV1()[ipoint];
472 }
473 tree->Draw(varVal,varCut);
474 for (Int_t ipoint=0; ipoint<entries; ipoint++){
475 px[ipoint]= tree->GetV3()[ipoint];
476 py[ipoint]= tree->GetV2()[ipoint];
477 pz[ipoint]= tree->GetV1()[ipoint];
478 }
479
480 //
481 TLinearFitter fitter(4,"hyp3");
482 for (Int_t ipoint=0; ipoint<entries; ipoint++){
483 Float_t val = pz[ipoint]*pz[ipoint];
484 Float_t err = 2*pz[ipoint]*TMath::Sqrt(ez[ipoint]*ez[ipoint]+fRatio*fRatio*pz[ipoint]*pz[ipoint]);
485 Double_t x[3];
486 x[0] = (type[ipoint]<0.5)? 0.:1.;
487 x[1] = px[ipoint];
488 x[2] = py[ipoint]*py[ipoint];
489 fitter.AddPoint(x,val,err);
490 }
491 fitter.Eval();
492 TVectorD param(4);
493 fitter.GetParameters(param);
494 param0[0] = param[0];
495 param0[1] = param[0]+param[1];
496 param0[2] = param[2];
497 param0[3] = param[3];
498 Float_t chi2 = fitter.GetChisquare()/entries;
499 param0[4] = chi2;
500 error[0] = (fitter.GetParError(0)*TMath::Sqrt(chi2));
501 error[1] = (fitter.GetParError(1)*TMath::Sqrt(chi2));
502 error[2] = (fitter.GetParError(2)*TMath::Sqrt(chi2));
503 error[3] = (fitter.GetParError(3)*TMath::Sqrt(chi2));
504}
505
506void AliTPCClusterParam::FitRMSQ(TTree * tree, Int_t dim, Int_t type, Float_t *param0, Float_t *error){
507 //
508 // Fit z - angular dependence of resolution - Q scaling
509 //
510 // Int_t dim=0, type=0;
511 char varVal[100];
512 sprintf(varVal,"RMSm:AngleM/sqrt(QMean):Zm/QMean");
513 char varVal0[100];
514 sprintf(varVal0,"RMSm:AngleM:Zm");
515 //
516 char varErr[100];
517 sprintf(varErr,"sqrt((1./(100.*sqrt(12.))^2)+RMSe0^2):AngleS:Zs");
518 char varCut[100];
519 sprintf(varCut,"Dim==%d&&Pad==%d&&QMean>0",dim,type);
520 //
521 Int_t entries = tree->Draw(varVal,varCut);
522 Float_t px[20000], py[20000], pz[20000], pu[20000], pt[20000];
523 Float_t ex[20000], ey[20000], ez[20000];
524 //
525 tree->Draw(varErr,varCut);
526 for (Int_t ipoint=0; ipoint<entries; ipoint++){
527 ex[ipoint]= tree->GetV3()[ipoint];
528 ey[ipoint]= tree->GetV2()[ipoint];
529 ez[ipoint]= tree->GetV1()[ipoint];
530 }
531 tree->Draw(varVal,varCut);
532 for (Int_t ipoint=0; ipoint<entries; ipoint++){
533 px[ipoint]= tree->GetV3()[ipoint];
534 py[ipoint]= tree->GetV2()[ipoint];
535 pz[ipoint]= tree->GetV1()[ipoint];
536 }
537 tree->Draw(varVal0,varCut);
538 for (Int_t ipoint=0; ipoint<entries; ipoint++){
539 pu[ipoint]= tree->GetV3()[ipoint];
540 pt[ipoint]= tree->GetV2()[ipoint];
541 }
542
543 //
544 TLinearFitter fitter(5,"hyp4");
545 for (Int_t ipoint=0; ipoint<entries; ipoint++){
546 Float_t val = pz[ipoint]*pz[ipoint];
547 Float_t err = 2*pz[ipoint]*TMath::Sqrt(ez[ipoint]*ez[ipoint]+fRatio*fRatio*pz[ipoint]*pz[ipoint]);
548 Double_t x[4];
549 x[0] = pu[ipoint];
550 x[1] = pt[ipoint]*pt[ipoint];
551 x[2] = px[ipoint];
552 x[3] = py[ipoint]*py[ipoint];
553 fitter.AddPoint(x,val,err);
554 }
555
556 fitter.Eval();
557 TVectorD param(5);
558 fitter.GetParameters(param);
559 param0[0] = param[0];
560 param0[1] = param[1];
561 param0[2] = param[2];
562 param0[3] = param[3];
563 param0[4] = param[4];
564 Float_t chi2 = fitter.GetChisquare()/entries;
565 param0[5] = chi2;
566 error[0] = (fitter.GetParError(0)*TMath::Sqrt(chi2));
567 error[1] = (fitter.GetParError(1)*TMath::Sqrt(chi2));
568 error[2] = (fitter.GetParError(2)*TMath::Sqrt(chi2));
569 error[3] = (fitter.GetParError(3)*TMath::Sqrt(chi2));
570 error[4] = (fitter.GetParError(4)*TMath::Sqrt(chi2));
571}
572
573
574void AliTPCClusterParam::FitRMSSigma(TTree * tree, Int_t dim, Int_t type, Float_t *param0, Float_t */*error*/){
575 //
576 // Fit z - angular dependence of resolution - Q scaling
577 //
578 // Int_t dim=0, type=0;
579 char varVal[100];
580 sprintf(varVal,"RMSs:RMSm");
581 //
582 char varCut[100];
583 sprintf(varCut,"Dim==%d&&Pad==%d&&QMean<0",dim,type);
584 //
585 Int_t entries = tree->Draw(varVal,varCut);
586 Float_t px[20000], py[20000];
587 //
588 tree->Draw(varVal,varCut);
589 for (Int_t ipoint=0; ipoint<entries; ipoint++){
590 px[ipoint]= tree->GetV2()[ipoint];
591 py[ipoint]= tree->GetV1()[ipoint];
592 }
593 TLinearFitter fitter(2,"pol1");
594 for (Int_t ipoint=0; ipoint<entries; ipoint++){
595 Float_t val = py[ipoint];
596 Float_t err = fRatio*px[ipoint];
597 Double_t x[4];
598 x[0] = px[ipoint];
599 fitter.AddPoint(x,val,err);
600 }
601 fitter.Eval();
602 param0[0]= fitter.GetParameter(0);
603 param0[1]= fitter.GetParameter(1);
604}
605
606
607
608Float_t AliTPCClusterParam::GetError0(Int_t dim, Int_t type, Float_t z, Float_t angle){
609 //
610 //
611 //
612 Float_t value=0;
613 value += fParamS0[dim][type][0];
614 value += fParamS0[dim][type][1]*z;
615 value += fParamS0[dim][type][2]*angle*angle;
616 value = TMath::Sqrt(TMath::Abs(value));
617 return value;
618}
619
620
621Float_t AliTPCClusterParam::GetError0Par(Int_t dim, Int_t type, Float_t z, Float_t angle){
622 //
623 //
624 //
625 Float_t value=0;
626 value += fParamS0Par[dim][type][0];
627 value += fParamS0Par[dim][type][1]*z;
628 value += fParamS0Par[dim][type][2]*angle*angle;
629 value += fParamS0Par[dim][type][3]*z*z;
630 value += fParamS0Par[dim][type][4]*angle*angle*angle*angle;
631 value += fParamS0Par[dim][type][5]*z*angle*angle;
632 value = TMath::Sqrt(TMath::Abs(value));
633 return value;
634}
635
636
637
638Float_t AliTPCClusterParam::GetError1(Int_t dim, Int_t type, Float_t z, Float_t angle){
639 //
640 //
641 //
642 Float_t value=0;
643 Float_t length=0.75;
644 if (type==1) length=1;
645 if (type==2) length=1.5;
646 value += fParamS1[dim][0];
647 value += fParamS1[dim][1]*z/length;
648 value += fParamS1[dim][2]*angle*angle*length;
649 value = TMath::Sqrt(TMath::Abs(value));
650 return value;
651}
652
653Float_t AliTPCClusterParam::GetErrorQ(Int_t dim, Int_t type, Float_t z, Float_t angle, Float_t Qmean){
654 //
655 //
656 //
657 Float_t value=0;
658 value += fParamSQ[dim][type][0];
659 value += fParamSQ[dim][type][1]*z;
660 value += fParamSQ[dim][type][2]*angle*angle;
661 value += fParamSQ[dim][type][3]*z/Qmean;
662 value += fParamSQ[dim][type][4]*angle*angle/Qmean;
663 value = TMath::Sqrt(TMath::Abs(value));
664 return value;
665
666
667}
668
669Float_t AliTPCClusterParam::GetErrorQPar(Int_t dim, Int_t type, Float_t z, Float_t angle, Float_t Qmean){
670 //
671 //
672 //
673 Float_t value=0;
674 value += fParamSQPar[dim][type][0];
675 value += fParamSQPar[dim][type][1]*z;
676 value += fParamSQPar[dim][type][2]*angle*angle;
677 value += fParamSQPar[dim][type][3]*z*z;
678 value += fParamSQPar[dim][type][4]*angle*angle*angle*angle;
679 value += fParamSQPar[dim][type][5]*z*angle*angle;
680 value += fParamSQPar[dim][type][6]*z/Qmean;
681 value += fParamSQPar[dim][type][7]*angle*angle/Qmean;
682 value = TMath::Sqrt(TMath::Abs(value));
683 return value;
684
685
686}
687
688Float_t AliTPCClusterParam::GetErrorQParScaled(Int_t dim, Int_t type, Float_t z, Float_t angle, Float_t Qmean){
689 //
690 //
691 //
692 Float_t value=0;
693 value += fParamSQPar[dim][type][0];
694 value += fParamSQPar[dim][type][1]*z;
695 value += fParamSQPar[dim][type][2]*angle*angle;
696 value += fParamSQPar[dim][type][3]*z*z;
697 value += fParamSQPar[dim][type][4]*angle*angle*angle*angle;
698 value += fParamSQPar[dim][type][5]*z*angle*angle;
699 value += fParamSQPar[dim][type][6]*z/Qmean;
700 value += fParamSQPar[dim][type][7]*angle*angle/Qmean;
701 Float_t valueMean = GetError0Par(dim,type,z,angle);
702 value -= 0.35*0.35*valueMean*valueMean;
703 value = TMath::Sqrt(TMath::Abs(value));
704 return value;
705
706
707}
708
709Float_t AliTPCClusterParam::GetRMS0(Int_t dim, Int_t type, Float_t z, Float_t angle){
710 //
711 // calculate mean RMS of cluster - z,angle - parameters for each pad and dimension separatelly
712 //
713 Float_t value=0;
714 value += fParamRMS0[dim][type][0];
715 value += fParamRMS0[dim][type][1]*z;
716 value += fParamRMS0[dim][type][2]*angle*angle;
717 value = TMath::Sqrt(TMath::Abs(value));
718 return value;
719}
720
721Float_t AliTPCClusterParam::GetRMS1(Int_t dim, Int_t type, Float_t z, Float_t angle){
722 //
723 // calculate mean RMS of cluster - z,angle - pad length scalling
724 //
725 Float_t value=0;
726 Float_t length=0.75;
727 if (type==1) length=1;
728 if (type==2) length=1.5;
729 if (type==0){
730 value += fParamRMS1[dim][0];
731 }else{
732 value += fParamRMS1[dim][1];
733 }
734 value += fParamRMS1[dim][2]*z;
735 value += fParamRMS1[dim][3]*angle*angle*length*length;
736 value = TMath::Sqrt(TMath::Abs(value));
737 return value;
738}
739
740Float_t AliTPCClusterParam::GetRMSQ(Int_t dim, Int_t type, Float_t z, Float_t angle, Float_t Qmean){
741 //
742 // calculate mean RMS of cluster - z,angle, Q dependence
743 //
744 Float_t value=0;
745 value += fParamRMSQ[dim][type][0];
746 value += fParamRMSQ[dim][type][1]*z;
747 value += fParamRMSQ[dim][type][2]*angle*angle;
748 value += fParamRMSQ[dim][type][3]*z/Qmean;
749 value += fParamRMSQ[dim][type][4]*angle*angle/Qmean;
750 value = TMath::Sqrt(TMath::Abs(value));
751 return value;
752}
753
754Float_t AliTPCClusterParam::GetRMSSigma(Int_t dim, Int_t type, Float_t z, Float_t angle, Float_t Qmean){
755 //
756 // calculates RMS of signal shape fluctuation
757 //
758 Float_t mean = GetRMSQ(dim,type,z,angle,Qmean);
759 Float_t value = fRMSSigmaFit[dim][type][0];
760 value+= fRMSSigmaFit[dim][type][1]*mean;
761 return value;
762}
763
764Float_t AliTPCClusterParam::GetShapeFactor(Int_t dim, Int_t type, Float_t z, Float_t angle, Float_t Qmean, Float_t rmsL, Float_t rmsM){
765 //
766 // calculates vallue - sigma distortion contribution
767 //
768 Double_t value =0;
769 //
770 Float_t rmsMeanQ = GetRMSQ(dim,type,z,angle,Qmean);
771 if (rmsL<rmsMeanQ) return value;
772 //
773 Float_t rmsSigma = GetRMSSigma(dim,type,z,angle,Qmean);
774 //
775 if ((rmsM-rmsMeanQ)>2.0*(rmsSigma+fErrorRMSSys[dim])){
776 //1.5 sigma cut on mean
777 value+= rmsL*rmsL+2*rmsM*rmsM-3*rmsMeanQ*rmsMeanQ;
778 }else{
779 if ((rmsL-rmsMeanQ)>3.*(rmsSigma+fErrorRMSSys[dim])){
780 //3 sigma cut on local
781 value+= rmsL*rmsL-rmsMeanQ*rmsMeanQ;
782 }
783 }
784 return TMath::Sqrt(value);
785}
786
787
788
789void AliTPCClusterParam::FitData(TTree * tree){
790 //
791 // make fits for error param and shape param
792 //
793 FitResol(tree);
794 FitRMS(tree);
795
796}
797
798void AliTPCClusterParam::FitResol(TTree * tree){
799 //
800 SetInstance(this);
801 for (Int_t idir=0;idir<2; idir++){
802 for (Int_t itype=0; itype<3; itype++){
803 Float_t param0[10];
804 Float_t error0[10];
805 // model error param
806 FitResol0(tree, idir, itype,param0,error0);
807 printf("\nResol\t%d\t%d\tchi2=%f\n",idir,itype,param0[3]);
808 printf("%f\t%f\t%f\n", param0[0],param0[1],param0[2]);
809 printf("%f\t%f\t%f\n", error0[0],error0[1],error0[2]);
810 for (Int_t ipar=0;ipar<4; ipar++){
811 fParamS0[idir][itype][ipar] = param0[ipar];
812 fErrorS0[idir][itype][ipar] = param0[ipar];
813 }
814 // error param with parabolic correction
815 FitResol0Par(tree, idir, itype,param0,error0);
816 printf("\nResolPar\t%d\t%d\tchi2=%f\n",idir,itype,param0[6]);
817 printf("%f\t%f\t%f\t%f\t%f\t%f\n", param0[0],param0[1],param0[2],param0[3],param0[4],param0[5]);
818 printf("%f\t%f\t%f\t%f\t%f\t%f\n", error0[0],error0[1],error0[2],error0[3],error0[4],error0[5]);
819 for (Int_t ipar=0;ipar<7; ipar++){
820 fParamS0Par[idir][itype][ipar] = param0[ipar];
821 fErrorS0Par[idir][itype][ipar] = param0[ipar];
822 }
823 //
824 FitResolQ(tree, idir, itype,param0,error0);
825 printf("\nResolQ\t%d\t%d\tchi2=%f\n",idir,itype,param0[5]);
826 printf("%f\t%f\t%f\t%f\t%f\n", param0[0],param0[1],param0[2],param0[3],param0[4]);
827 printf("%f\t%f\t%f\t%f\t%f\n", error0[0],error0[1],error0[2],error0[3],error0[4]);
828 for (Int_t ipar=0;ipar<6; ipar++){
829 fParamSQ[idir][itype][ipar] = param0[ipar];
830 fErrorSQ[idir][itype][ipar] = param0[ipar];
831 }
832 //
833 FitResolQPar(tree, idir, itype,param0,error0);
834 printf("\nResolQ\t%d\t%d\tchi2=%f\n",idir,itype,param0[8]);
835 printf("%f\t%f\t%f\t%f\t%f\t%f\t%f\t%f\n", param0[0],param0[1],param0[2],param0[3],param0[4],param0[5],param0[6],param0[7]);
836 printf("%f\t%f\t%f\t%f\t%f\t%f\t%f\t%f\n", error0[0],error0[1],error0[2],error0[3],error0[4],error0[5],error0[6],error0[7]);
837 for (Int_t ipar=0;ipar<9; ipar++){
838 fParamSQPar[idir][itype][ipar] = param0[ipar];
839 fErrorSQPar[idir][itype][ipar] = param0[ipar];
840 }
841 }
842 }
843 //
844 printf("Resol z-scaled\n");
845 for (Int_t idir=0;idir<2; idir++){
846 Float_t param0[4];
847 Float_t error0[4];
848 FitResol1(tree, idir,param0,error0);
849 printf("\nResol\t%d\tchi2=%f\n",idir,param0[3]);
850 printf("%f\t%f\t%f\n", param0[0],param0[1],param0[2]);
851 printf("%f\t%f\t%f\n", error0[0],error0[1],error0[2]);
852 for (Int_t ipar=0;ipar<4; ipar++){
853 fParamS1[idir][ipar] = param0[ipar];
854 fErrorS1[idir][ipar] = param0[ipar];
855 }
856 }
857
858 for (Int_t idir=0;idir<2; idir++){
859 printf("\nDirection %d\n",idir);
860 printf("%d\t%f\t%f\t%f\n", -1,fParamS1[idir][0],fParamS1[idir][1],fParamS1[idir][2]);
861 for (Int_t itype=0; itype<3; itype++){
862 Float_t length=0.75;
863 if (itype==1) length=1;
864 if (itype==2) length=1.5;
865 printf("%d\t%f\t%f\t%f\n", itype,fParamS0[idir][itype][0],fParamS0[idir][itype][1]*TMath::Sqrt(length),fParamS0[idir][itype][2]/TMath::Sqrt(length));
866 }
867 }
868}
869
870
871
872void AliTPCClusterParam::FitRMS(TTree * tree){
873 //
874 SetInstance(this);
875 for (Int_t idir=0;idir<2; idir++){
876 for (Int_t itype=0; itype<3; itype++){
877 Float_t param0[6];
878 Float_t error0[6];
879 FitRMS0(tree, idir, itype,param0,error0);
880 printf("\nRMS\t%d\t%d\tchi2=%f\n",idir,itype,param0[3]);
881 printf("%f\t%f\t%f\n", param0[0],param0[1],param0[2]);
882 printf("%f\t%f\t%f\n", error0[0],error0[1],error0[2]);
883 for (Int_t ipar=0;ipar<4; ipar++){
884 fParamRMS0[idir][itype][ipar] = param0[ipar];
885 fErrorRMS0[idir][itype][ipar] = param0[ipar];
886 }
887 FitRMSQ(tree, idir, itype,param0,error0);
888 printf("\nRMSQ\t%d\t%d\tchi2=%f\n",idir,itype,param0[5]);
889 printf("%f\t%f\t%f\t%f\t%f\n", param0[0],param0[1],param0[2],param0[3],param0[4]);
890 printf("%f\t%f\t%f\t%f\t%f\n", error0[0],error0[1],error0[2],error0[3],error0[4]);
891 for (Int_t ipar=0;ipar<6; ipar++){
892 fParamRMSQ[idir][itype][ipar] = param0[ipar];
893 fErrorRMSQ[idir][itype][ipar] = param0[ipar];
894 }
895 }
896 }
897 //
898 printf("RMS z-scaled\n");
899 for (Int_t idir=0;idir<2; idir++){
900 Float_t param0[5];
901 Float_t error0[5];
902 FitRMS1(tree, idir,param0,error0);
903 printf("\nRMS\t%d\tchi2=%f\n",idir,param0[4]);
904 printf("%f\t%f\t%f\t%f\n", param0[0],param0[1],param0[2], param0[3]);
905 printf("%f\t%f\t%f\t%f\n", error0[0],error0[1],error0[2], error0[3]);
906 for (Int_t ipar=0;ipar<5; ipar++){
907 fParamRMS1[idir][ipar] = param0[ipar];
908 fErrorRMS1[idir][ipar] = param0[ipar];
909 }
910 }
911
912 for (Int_t idir=0;idir<2; idir++){
913 printf("\nDirection %d\n",idir);
914 printf("%d\t%f\t%f\t%f\t%f\n", -1,fParamRMS1[idir][0],fParamRMS1[idir][1],fParamRMS1[idir][2], fParamRMS1[idir][3]);
915 for (Int_t itype=0; itype<3; itype++){
916 Float_t length=0.75;
917 if (itype==1) length=1;
918 if (itype==2) length=1.5;
919 if (itype==0) printf("%d\t%f\t\t\t%f\t%f\n", itype,fParamRMS0[idir][itype][0],fParamRMS0[idir][itype][1],fParamRMS0[idir][itype][2]/length);
920 if (itype>0) printf("%d\t\t\t%f\t%f\t%f\n", itype,fParamRMS0[idir][itype][0],fParamRMS0[idir][itype][1],fParamRMS0[idir][itype][2]/length);
921 }
922 }
923 //
924 // Fit RMS sigma
925 //
926 printf("RMS fluctuation parameterization \n");
927 for (Int_t idir=0;idir<2; idir++){
928 for (Int_t itype=0; itype<3; itype++){
929 Float_t param0[5];
930 Float_t error0[5];
931 FitRMSSigma(tree, idir,itype,param0,error0);
932 printf("\t%d\t%d\t%f\t%f\n", idir, itype, param0[0],param0[1]);
933 for (Int_t ipar=0;ipar<2; ipar++){
934 fRMSSigmaFit[idir][itype][ipar] = param0[ipar];
935 }
936 }
937 }
938 //
939 // store systematic error end RMS fluctuation parameterization
940 //
941 TH1F hratio("hratio","hratio",100,-0.1,0.1);
942 tree->Draw("(RMSm-AliTPCClusterParam::SGetRMSQ(Dim,Pad,Zm,AngleM,QMean))/RMSm>>hratio","Dim==0&&QMean>0");
943 fErrorRMSSys[0] = hratio.GetRMS();
944 tree->Draw("(RMSm-AliTPCClusterParam::SGetRMSQ(Dim,Pad,Zm,AngleM,QMean))/RMSm>>hratio","Dim==1&&QMean>0");
945 fErrorRMSSys[1] = hratio.GetRMS();
946 TH1F hratioR("hratioR","hratioR",100,0,0.2);
947 tree->Draw("RMSs/RMSm>>hratioR","Dim==0&&QMean>0");
948 fRMSSigmaRatio[0][0]=hratioR.GetMean();
949 fRMSSigmaRatio[0][1]=hratioR.GetRMS();
950 tree->Draw("RMSs/RMSm>>hratioR","Dim==1&&QMean>0");
951 fRMSSigmaRatio[1][0]=hratioR.GetMean();
952 fRMSSigmaRatio[1][1]=hratioR.GetRMS();
953}
954
955void AliTPCClusterParam::Test(TTree * tree, const char *output){
956 //
957 // Draw standard quality histograms to output file
958 //
959 SetInstance(this);
960 TFile f(output,"recreate");
961 f.cd();
962 //
963 // 1D histograms - resolution
964 //
965 for (Int_t idim=0; idim<2; idim++){
966 for (Int_t ipad=0; ipad<3; ipad++){
967 char hname1[300];
968 char hcut1[300];
969 char hexp1[300];
970 //
971 sprintf(hname1,"Delta0 Dir %d Pad %d",idim,ipad);
972 sprintf(hcut1,"Dim==%d&&QMean<0&&Pad==%d",idim,ipad);
973 sprintf(hexp1,"(Resol-AliTPCClusterParam::SGetError0(Dim,Pad,Zm,AngleM))/Resol>>%s",hname1);
974 TH1F his1DRel0(hname1, hname1, 100,-0.2, 0.2);
975 sprintf(hname1,"Dim==%d&&QMean<0&&Pad=%d",idim,ipad);
976 tree->Draw(hexp1,hcut1,"");
977 his1DRel0.Write();
978 //
979 sprintf(hname1,"Delta0Par Dir %d Pad %d",idim,ipad);
980 sprintf(hcut1,"Dim==%d&&QMean<0&&Pad==%d",idim,ipad);
981 sprintf(hexp1,"(Resol-AliTPCClusterParam::SGetError0Par(Dim,Pad,Zm,AngleM))/Resol>>%s",hname1);
982 TH1F his1DRel0Par(hname1, hname1, 100,-0.2, 0.2);
983 sprintf(hname1,"Dim==%d&&QMean<0&&Pad=%d",idim,ipad);
984 tree->Draw(hexp1,hcut1,"");
985 his1DRel0Par.Write();
986 //
987 }
988 }
989 //
990 // 2D histograms - resolution
991 //
992 for (Int_t idim=0; idim<2; idim++){
993 for (Int_t ipad=0; ipad<3; ipad++){
994 char hname1[300];
995 char hcut1[300];
996 char hexp1[300];
997 //
998 sprintf(hname1,"2DDelta0 Dir %d Pad %d",idim,ipad);
999 sprintf(hcut1,"Dim==%d&&QMean<0&&Pad==%d",idim,ipad);
1000 sprintf(hexp1,"(Resol-AliTPCClusterParam::SGetError0(Dim,Pad,Zm,AngleM))/Resol:AngleM:Zm>>%s",hname1);
1001 TProfile2D profDRel0(hname1, hname1, 6,0,250,6,0,1);
1002 sprintf(hname1,"Dim==%d&&QMean<0&&Pad=%d",idim,ipad);
1003 tree->Draw(hexp1,hcut1,"");
1004 profDRel0.Write();
1005 //
1006 sprintf(hname1,"2DDelta0Par Dir %d Pad %d",idim,ipad);
1007 sprintf(hcut1,"Dim==%d&&QMean<0&&Pad==%d",idim,ipad);
1008 sprintf(hexp1,"(Resol-AliTPCClusterParam::SGetError0Par(Dim,Pad,Zm,AngleM))/Resol:AngleM:Zm>>%s",hname1);
1009 TProfile2D profDRel0Par(hname1, hname1,6,0,250,6,0,1);
1010 sprintf(hname1,"Dim==%d&&QMean<0&&Pad=%d",idim,ipad);
1011 tree->Draw(hexp1,hcut1,"");
1012 profDRel0Par.Write();
1013 //
1014 }
1015 }
1016}
1017
1018
1019
1020void AliTPCClusterParam::Print(Option_t* /*option*/) const{
1021 //
1022 // Print param Information
1023 //
1024
1025 //
1026 // Error parameterization
1027 //
1028 printf("\nResolution Scaled factors\n");
1029 printf("Dir\tPad\tP0\t\tP1\t\tP2\t\tchi2\n");
1030 printf("Y\tall\t%f\t%f\t%f\t%f\n", TMath::Sqrt(TMath::Abs(fParamS1[0][0])),TMath::Sqrt(fParamS1[0][1]),
1031 TMath::Sqrt(fParamS1[0][2]),TMath::Sqrt(fParamS1[0][3]));
1032 for (Int_t ipad=0; ipad<3; ipad++){
1033 Float_t length=0.75;
1034 if (ipad==1) length=1;
1035 if (ipad==2) length=1.5;
1036 printf("\t%d\t%f\t%f\t%f\t%f\n", ipad,
1037 TMath::Sqrt(TMath::Abs(fParamS0[0][ipad][0])),
1038 TMath::Sqrt(fParamS0[0][ipad][1]*length),
1039 TMath::Sqrt(fParamS0[0][ipad][2]/length),
1040 TMath::Sqrt(fParamS0[0][ipad][3]));
1041 }
1042 for (Int_t ipad=0; ipad<3; ipad++){
1043 Float_t length=0.75;
1044 if (ipad==1) length=1;
1045 if (ipad==2) length=1.5;
1046 printf("\t%dPar\t%f\t%f\t%f\t%f\n", ipad,
1047 TMath::Sqrt(TMath::Abs(fParamS0Par[0][ipad][0])),
1048 TMath::Sqrt(fParamS0Par[0][ipad][1]*length),
1049 TMath::Sqrt(fParamS0Par[0][ipad][2]/length),
1050 TMath::Sqrt(fParamS0Par[0][ipad][6]));
1051 }
1052 printf("Z\tall\t%f\t%f\t%f\t%f\n", TMath::Sqrt(TMath::Abs(fParamS1[1][0])),TMath::Sqrt(fParamS1[1][1]),
1053 TMath::Sqrt(fParamS1[1][2]), TMath::Sqrt(fParamS1[1][3]));
1054
1055 for (Int_t ipad=0; ipad<3; ipad++){
1056 Float_t length=0.75;
1057 if (ipad==1) length=1;
1058 if (ipad==2) length=1.5;
1059 printf("\t%d\t%f\t%f\t%f\t%f\n", ipad,
1060 TMath::Sqrt(TMath::Abs(fParamS0[1][ipad][0])),
1061 TMath::Sqrt(fParamS0[1][ipad][1]*length),
1062 TMath::Sqrt(fParamS0[1][ipad][2]/length),
1063 TMath::Sqrt(fParamS0[1][ipad][3]));
1064 }
1065 for (Int_t ipad=0; ipad<3; ipad++){
1066 Float_t length=0.75;
1067 if (ipad==1) length=1;
1068 if (ipad==2) length=1.5;
1069 printf("\t%dPar\t%f\t%f\t%f\t%f\n", ipad,
1070 TMath::Sqrt(TMath::Abs(fParamS0Par[1][ipad][0])),
1071 TMath::Sqrt(fParamS0Par[1][ipad][1]*length),
1072 TMath::Sqrt(fParamS0Par[1][ipad][2]/length),
1073 TMath::Sqrt(fParamS0Par[1][ipad][6]));
1074 }
1075
1076 //
1077 // RMS scaling
1078 //
1079 printf("\n");
1080 printf("\nRMS Scaled factors\n");
1081 printf("Dir\tPad\tP00\t\tP01\t\tP1\t\tP2\t\tchi2\n");
1082 printf("Y\tall\t%f\t%f\t%f\t%f\t%f\n", TMath::Sqrt(TMath::Abs(fParamRMS1[0][0])),TMath::Sqrt(fParamRMS1[0][1]),
1083 TMath::Sqrt(fParamRMS1[0][2]),TMath::Sqrt(fParamRMS1[0][3]),TMath::Sqrt(fParamRMS1[0][4]));
1084 for (Int_t ipad=0; ipad<3; ipad++){
1085 Float_t length=0.75;
1086 if (ipad==1) length=1;
1087 if (ipad==2) length=1.5;
1088 if (ipad==0){
1089 printf("\t%d\t%f\t%f\t%f\t%f\t%f\n", ipad,
1090 TMath::Sqrt(TMath::Abs(fParamRMS0[0][ipad][0])),
1091 0.,
1092 TMath::Sqrt(fParamRMS0[0][ipad][1]),
1093 TMath::Sqrt(fParamRMS0[0][ipad][2]/(length*length)),
1094 TMath::Sqrt(fParamRMS0[0][ipad][3]));
1095 }else{
1096 printf("\t%d\t%f\t%f\t%f\t%f\t%f\n", ipad,
1097 0.,
1098 TMath::Sqrt(TMath::Abs(fParamRMS0[0][ipad][0])),
1099 TMath::Sqrt(fParamRMS0[0][ipad][1]),
1100 TMath::Sqrt(fParamRMS0[0][ipad][2]/(length*length)),
1101 TMath::Sqrt(fParamRMS0[0][ipad][3]));
1102 }
1103 }
1104 printf("\n");
1105 printf("Z\tall\t%f\t%f\t%f\t%f\t%f\n", TMath::Sqrt(TMath::Abs(fParamRMS1[1][0])),TMath::Sqrt(fParamRMS1[1][1]),
1106 TMath::Sqrt(fParamRMS1[1][2]),TMath::Sqrt(fParamRMS1[1][3]),TMath::Sqrt(fParamRMS1[1][4]));
1107 for (Int_t ipad=0; ipad<3; ipad++){
1108 Float_t length=0.75;
1109 if (ipad==1) length=1;
1110 if (ipad==2) length=1.5;
1111 if (ipad==0){
1112 printf("\t%d\t%f\t%f\t%f\t%f\t%f\n", ipad,
1113 TMath::Sqrt(TMath::Abs(fParamRMS0[1][ipad][0])),
1114 0.,
1115 TMath::Sqrt(fParamRMS0[1][ipad][1]),
1116 TMath::Sqrt(fParamRMS0[1][ipad][2]/(length*length)),
1117 TMath::Sqrt(fParamRMS0[1][ipad][3]));
1118 }else{
1119 printf("\t%d\t%f\t%f\t%f\t%f\t%f\n", ipad,
1120 0.,
1121 TMath::Sqrt(TMath::Abs(fParamRMS0[1][ipad][0])),
1122 TMath::Sqrt(fParamRMS0[1][ipad][1]),
1123 TMath::Sqrt(fParamRMS0[1][ipad][2]/(length*length)),
1124 TMath::Sqrt(fParamRMS0[1][ipad][3]));
1125 }
1126 }
1127}
1128
1129
1130
1131
1132
0a65832b 1133Float_t AliTPCClusterParam::Qnorm(Int_t ipad, Int_t itype, Float_t dr, Float_t ty, Float_t tz){
1134 // get Q normalization
1135 // type - 0 Qtot 1 Qmax
1136 // ipad - 0 (0.75 cm) ,1 (1 cm), 2 (1.5 cm)
1137 //
1138 //
1139 //formula= dr++tz++ty++dr*tz++dr*ty++ty*tz++dr**2++ty**2++tz**2
1140 if (!fQNorm) return 0;
1141 TVectorD * norm = (TVectorD*)fQNorm->At(3*itype+ipad);
1142 if (!norm) return 0;
1143 TVectorD &no = *norm;
1144 Float_t res= no[0]+no[1]*dr+no[2]*tz+no[3]*ty+no[4]*dr*tz+no[5]*dr*ty+no[6]*ty*tz
1145 +no[7]*dr*dr+no[8]*ty*ty+no[9]*tz*tz;
1146
1147 return res;
1148}
1149
1150
1151
1152void AliTPCClusterParam::SetQnorm(Int_t ipad, Int_t itype, TVectorD * norm){
1153 //
1154 // set normalization
1155 //
1156 // type - 0 Qtot 1 Qmax
1157 // ipad - 0 (0.75 cm) ,1 (1 cm), 2 (1.5 cm)
1158 //
1159
1160 if (fQNorm==0) fQNorm = new TObjArray(6);
1161 fQNorm->AddAt(new TVectorD(*norm), itype*3+ipad);
1162}