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