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