Correct sign for calculated b_y (outside measured region).
[u/mrichter/AliRoot.git] / TPC / AliTPCCalPad.cxx
CommitLineData
07627591 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/* $Id$ */
17
18///////////////////////////////////////////////////////////////////////////////
19// //
a6d2bd0c 20// TPC calibration class for parameters which are saved per pad //
21// Each AliTPCCalPad consists of 72 AliTPCCalROC-objects //
07627591 22// //
23///////////////////////////////////////////////////////////////////////////////
24
25#include "AliTPCCalPad.h"
26#include "AliTPCCalROC.h"
184bcc16 27#include <TObjArray.h>
28#include <TAxis.h>
29#include <TGraph.h>
200be8a6 30#include <TGraph2D.h>
31#include <TH2F.h>
586007f3 32#include "TTreeStream.h"
ca5dca67 33#include "TFile.h"
34#include "TKey.h"
586007f3 35
07627591 36ClassImp(AliTPCCalPad)
37
38//_____________________________________________________________________________
39AliTPCCalPad::AliTPCCalPad():TNamed()
40{
41 //
42 // AliTPCCalPad default constructor
43 //
44
45 for (Int_t isec = 0; isec < kNsec; isec++) {
46 fROC[isec] = 0;
47 }
48
49}
50
51//_____________________________________________________________________________
52AliTPCCalPad::AliTPCCalPad(const Text_t *name, const Text_t *title)
53 :TNamed(name,title)
54{
55 //
56 // AliTPCCalPad constructor
57 //
58 for (Int_t isec = 0; isec < kNsec; isec++) {
59 fROC[isec] = new AliTPCCalROC(isec);
60 }
61}
62
63
64//_____________________________________________________________________________
65AliTPCCalPad::AliTPCCalPad(const AliTPCCalPad &c):TNamed(c)
66{
67 //
68 // AliTPCCalPad copy constructor
69 //
70
e6c51e6e 71 for (Int_t isec = 0; isec < kNsec; isec++) {
7fb8d357 72 fROC[isec] = 0;
73 if (c.fROC[isec])
74 fROC[isec] = new AliTPCCalROC(*(c.fROC[isec]));
e6c51e6e 75 }
07627591 76}
77
184bcc16 78//_____________________________________________________________________________
79AliTPCCalPad::AliTPCCalPad(TObjArray * array):TNamed()
80{
81 //
82 // AliTPCCalPad default constructor
83 //
84
85 for (Int_t isec = 0; isec < kNsec; isec++) {
86 fROC[isec] = (AliTPCCalROC *)array->At(isec);
87 }
88
89}
90
91
07627591 92///_____________________________________________________________________________
93AliTPCCalPad::~AliTPCCalPad()
94{
95 //
96 // AliTPCCalPad destructor
97 //
98
99 for (Int_t isec = 0; isec < kNsec; isec++) {
100 if (fROC[isec]) {
101 delete fROC[isec];
102 fROC[isec] = 0;
103 }
104 }
105
106}
107
108//_____________________________________________________________________________
109AliTPCCalPad &AliTPCCalPad::operator=(const AliTPCCalPad &c)
110{
111 //
112 // Assignment operator
113 //
114
115 if (this != &c) ((AliTPCCalPad &) c).Copy(*this);
116 return *this;
117
118}
119
120//_____________________________________________________________________________
121void AliTPCCalPad::Copy(TObject &c) const
122{
123 //
124 // Copy function
125 //
126
127 for (Int_t isec = 0; isec < kNsec; isec++) {
128 if (fROC[isec]) {
129 fROC[isec]->Copy(*((AliTPCCalPad &) c).fROC[isec]);
130 }
131 }
132 TObject::Copy(c);
133}
184bcc16 134
a6d2bd0c 135
136void AliTPCCalPad::SetCalROC(AliTPCCalROC* roc, Int_t sector){
137 //
138 // Set AliTPCCalROC copies values from 'roc'
139 // if sector == -1 the sector specified in 'roc' is used
140 // else sector specified in 'roc' is ignored and specified sector is filled
141 //
142 if (sector == -1) sector = roc->GetSector();
143 for (UInt_t ichannel = 0; ichannel < roc->GetNchannels(); ichannel++)
144 fROC[sector]->SetValue(ichannel, roc->GetValue(ichannel));
145}
146
147
148
90127643 149//_____________________________________________________________________________
150void AliTPCCalPad::Add(Float_t c1)
151{
152 //
a6d2bd0c 153 // add constant c1 to all channels of all ROCs
90127643 154 //
155
156 for (Int_t isec = 0; isec < kNsec; isec++) {
157 if (fROC[isec]){
158 fROC[isec]->Add(c1);
159 }
160 }
161}
162
163//_____________________________________________________________________________
164void AliTPCCalPad::Multiply(Float_t c1)
165{
a6d2bd0c 166 //
167 // multiply each channel of all ROCs with c1
168 //
90127643 169 for (Int_t isec = 0; isec < kNsec; isec++) {
170 if (fROC[isec]){
171 fROC[isec]->Multiply(c1);
172 }
173 }
174}
175
176//_____________________________________________________________________________
177void AliTPCCalPad::Add(const AliTPCCalPad * pad, Double_t c1)
178{
a6d2bd0c 179 //
180 // multiply AliTPCCalPad 'pad' by c1 and add each channel to the coresponing channel in all ROCs
181 // - pad by pad -
182 //
90127643 183 for (Int_t isec = 0; isec < kNsec; isec++) {
184 if (fROC[isec]){
185 fROC[isec]->Add(pad->GetCalROC(isec),c1);
186 }
187 }
188}
189
190//_____________________________________________________________________________
191void AliTPCCalPad::Multiply(const AliTPCCalPad * pad)
192{
a6d2bd0c 193 //
194 // multiply each channel of all ROCs with the coresponding channel of 'pad'
195 // - pad by pad -
196 //
197 for (Int_t isec = 0; isec < kNsec; isec++) {
90127643 198 if (fROC[isec]){
199 fROC[isec]->Multiply(pad->GetCalROC(isec));
200 }
201 }
202}
203
204//_____________________________________________________________________________
205void AliTPCCalPad::Divide(const AliTPCCalPad * pad)
206{
a6d2bd0c 207 //
208 // divide each channel of all ROCs by the coresponding channel of 'pad'
209 // - pad by pad -
210 //
90127643 211 for (Int_t isec = 0; isec < kNsec; isec++) {
212 if (fROC[isec]){
213 fROC[isec]->Divide(pad->GetCalROC(isec));
214 }
215 }
216}
217
218//_____________________________________________________________________________
184bcc16 219TGraph * AliTPCCalPad::MakeGraph(Int_t type, Float_t ratio){
220 //
221 // type=1 - mean
222 // 2 - median
223 // 3 - LTM
a6d2bd0c 224 //
184bcc16 225 Int_t npoints = 0;
226 for (Int_t i=0;i<72;i++) if (fROC[i]) npoints++;
227 TGraph * graph = new TGraph(npoints);
228 npoints=0;
229 for (Int_t isec=0;isec<72;isec++){
230 if (!fROC[isec]) continue;
231 if (type==0) graph->SetPoint(npoints,isec,fROC[isec]->GetMean());
232 if (type==1) graph->SetPoint(npoints,isec,fROC[isec]->GetMedian());
233 if (type==2) graph->SetPoint(npoints,isec,fROC[isec]->GetLTM(0,ratio));
234 npoints++;
235 }
236
237 graph->GetXaxis()->SetTitle("Sector");
238 if (type==0) {
239 graph->GetYaxis()->SetTitle("Mean");
240 graph->SetMarkerStyle(22);
241 }
242 if (type==1) {
243 graph->GetYaxis()->SetTitle("Median");
244 graph->SetMarkerStyle(22);
245 }
246 if (type==2) {
247 graph->GetYaxis()->SetTitle(Form("Mean%f",ratio));
248 graph->SetMarkerStyle(24);
249 }
250
251 return graph;
252}
200be8a6 253
90127643 254//_____________________________________________________________________________
255Double_t AliTPCCalPad::GetMeanRMS(Double_t &rms)
256{
257 //
a6d2bd0c 258 // Calculates mean and RMS of all ROCs
90127643 259 //
260 Double_t sum = 0, sum2 = 0, n=0, val=0;
261 for (Int_t isec = 0; isec < kNsec; isec++) {
262 AliTPCCalROC *calRoc = fROC[isec];
263 if ( calRoc ){
264 for (UInt_t irow=0; irow<calRoc->GetNrows(); irow++){
265 for (UInt_t ipad=0; ipad<calRoc->GetNPads(irow); ipad++){
266 val = calRoc->GetValue(irow,ipad);
267 sum+=val;
268 sum2+=val*val;
269 n++;
270 }
271 }
272
273 }
274 }
275 Double_t n1 = 1./n;
276 Double_t mean = sum*n1;
277 rms = TMath::Sqrt(TMath::Abs(sum2*n1-mean*mean));
278 return mean;
279}
280
281
282//_____________________________________________________________________________
ca5dca67 283Double_t AliTPCCalPad::GetMean(AliTPCCalPad* outlierPad)
90127643 284{
285 //
286 // return mean of the mean of all ROCs
287 //
288 Double_t arr[kNsec];
289 Int_t n=0;
290 for (Int_t isec = 0; isec < kNsec; isec++) {
ca5dca67 291 AliTPCCalROC *calRoc = fROC[isec];
292 if ( calRoc ){
293 AliTPCCalROC* outlierROC = 0;
294 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
295 arr[n] = calRoc->GetMean(outlierROC);
296 n++;
297 }
90127643 298 }
299 return TMath::Mean(n,arr);
300}
301
302//_____________________________________________________________________________
ca5dca67 303Double_t AliTPCCalPad::GetRMS(AliTPCCalPad* outlierPad)
90127643 304{
305 //
306 // return mean of the RMS of all ROCs
307 //
308 Double_t arr[kNsec];
309 Int_t n=0;
310 for (Int_t isec = 0; isec < kNsec; isec++) {
ca5dca67 311 AliTPCCalROC *calRoc = fROC[isec];
312 if ( calRoc ){
313 AliTPCCalROC* outlierROC = 0;
314 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
315 arr[n] = calRoc->GetRMS(outlierROC);
316 n++;
317 }
90127643 318 }
319 return TMath::Mean(n,arr);
320}
321
322//_____________________________________________________________________________
ca5dca67 323Double_t AliTPCCalPad::GetMedian(AliTPCCalPad* outlierPad)
90127643 324{
325 //
326 // return mean of the median of all ROCs
327 //
328 Double_t arr[kNsec];
329 Int_t n=0;
330 for (Int_t isec = 0; isec < kNsec; isec++) {
ca5dca67 331 AliTPCCalROC *calRoc = fROC[isec];
332 if ( calRoc ){
333 AliTPCCalROC* outlierROC = 0;
334 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
335 arr[n] = calRoc->GetMedian(outlierROC);
336 n++;
337 }
90127643 338 }
339 return TMath::Mean(n,arr);
340}
341
342//_____________________________________________________________________________
ca5dca67 343Double_t AliTPCCalPad::GetLTM(Double_t *sigma, Double_t fraction, AliTPCCalPad* outlierPad)
90127643 344{
345 //
346 // return mean of the LTM and sigma of all ROCs
347 //
348 Double_t arrm[kNsec];
349 Double_t arrs[kNsec];
350 Double_t *sTemp=0x0;
351 Int_t n=0;
352
353 for (Int_t isec = 0; isec < kNsec; isec++) {
354 AliTPCCalROC *calRoc = fROC[isec];
355 if ( calRoc ){
356 if ( sigma ) sTemp=arrs+n;
ca5dca67 357 AliTPCCalROC* outlierROC = 0;
358 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
359 arrm[n] = calRoc->GetLTM(sTemp,fraction, outlierROC);
90127643 360 n++;
361 }
362 }
363 if ( sigma ) *sigma = TMath::Mean(n,arrs);
364 return TMath::Mean(n,arrm);
365}
366
367//_____________________________________________________________________________
368TH1F * AliTPCCalPad::MakeHisto1D(Float_t min, Float_t max,Int_t type){
369 //
370 // make 1D histo
371 // type -1 = user defined range
372 // 0 = nsigma cut nsigma=min
a6d2bd0c 373 //
90127643 374 if (type>=0){
375 if (type==0){
376 // nsigma range
377 Float_t mean = GetMean();
378 Float_t sigma = GetRMS();
379 Float_t nsigma = TMath::Abs(min);
380 min = mean-nsigma*sigma;
381 max = mean+nsigma*sigma;
382 }
383 if (type==1){
384 // fixed range
385 Float_t mean = GetMedian();
386 Float_t delta = min;
387 min = mean-delta;
388 max = mean+delta;
389 }
390 if (type==2){
391 //
392 // LTM mean +- nsigma
393 //
394 Double_t sigma;
395 Float_t mean = GetLTM(&sigma,max);
396 sigma*=min;
397 min = mean-sigma;
398 max = mean+sigma;
399 }
400 }
401 char name[1000];
402 sprintf(name,"%s Pad 1D",GetTitle());
403 TH1F * his = new TH1F(name,name,100, min,max);
404 for (Int_t isec = 0; isec < kNsec; isec++) {
405 if (fROC[isec]){
406 for (UInt_t irow=0; irow<fROC[isec]->GetNrows(); irow++){
407 UInt_t npads = (Int_t)fROC[isec]->GetNPads(irow);
408 for (UInt_t ipad=0; ipad<npads; ipad++){
409 his->Fill(fROC[isec]->GetValue(irow,ipad));
410 }
411 }
412 }
413 }
414 return his;
415}
416
417//_____________________________________________________________________________
200be8a6 418TH2F *AliTPCCalPad::MakeHisto2D(Int_t side){
419 //
420 // Make 2D graph
421 // side - specify the side A = 0 C = 1
422 // type - used types of determination of boundaries in z
a6d2bd0c 423 //
200be8a6 424 Float_t kEpsilon = 0.000000000001;
425 TH2F * his = new TH2F(GetName(), GetName(), 250,-250,250,250,-250,250);
426 AliTPCROC * roc = AliTPCROC::Instance();
427 for (Int_t isec=0; isec<72; isec++){
428 if (side==0 && isec%36>=18) continue;
429 if (side>0 && isec%36<18) continue;
430 if (fROC[isec]){
431 AliTPCCalROC * calRoc = fROC[isec];
432 for (UInt_t irow=0; irow<calRoc->GetNrows(); irow++)
433 for (UInt_t ipad=0; ipad<calRoc->GetNPads(irow); ipad++)
434 if (TMath::Abs(calRoc->GetValue(irow,ipad))>kEpsilon){
435 Float_t xyz[3];
436 roc->GetPositionGlobal(isec,irow,ipad,xyz);
437 Int_t binx = 1+TMath::Nint((xyz[0]+250.)*0.5);
438 Int_t biny = 1+TMath::Nint((xyz[1]+250.)*0.5);
439 Float_t value = calRoc->GetValue(irow,ipad);
440 his->SetBinContent(binx,biny,value);
441 }
442 }
443 }
444 his->SetXTitle("x (cm)");
445 his->SetYTitle("y (cm)");
446 return his;
447}
448
449
a6d2bd0c 450AliTPCCalPad* AliTPCCalPad::LocalFit(const char* padName, Int_t rowRadius, Int_t padRadius, AliTPCCalPad* PadOutliers, Bool_t robust, Double_t chi2Threshold, Double_t robustFraction){
451 //
452 // Loops over all AliTPCCalROCs and performs a localFit in each ROC
453 // AliTPCCalPad with fit-data is returned
454 // rowRadius and padRadius specifies a window around a given pad in one sector.
455 // The data of this window are fitted with a parabolic function.
456 // This function is evaluated at the pad's position.
457 // At the edges the window is shifted, so that the specified pad is not anymore in the center of the window.
458 // rowRadius - radius - rows to be used for smoothing
459 // padradius - radius - pads to be used for smoothing
460 // ROCoutlier - map of outliers - pads not to be used for local smoothing
461 // robust - robust method of fitting - (much slower)
462 // chi2Threshold: Threshold for chi2 when EvalRobust is called
463 // robustFraction: Fraction of data that will be used in EvalRobust
464 //
465 //
466 AliTPCCalPad* pad = new AliTPCCalPad(padName, padName);
467 for (Int_t isec = 0; isec < 72; isec++){
468 if (PadOutliers)
469 SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, PadOutliers->GetCalROC(isec), robust, chi2Threshold, robustFraction));
470 else
471 SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, 0, robust, chi2Threshold, robustFraction));
472 }
473 return pad;
474}
475
476
477AliTPCCalPad* AliTPCCalPad::GlobalFit(const char* padName, AliTPCCalPad* PadOutliers, Bool_t robust, Int_t fitType, Double_t chi2Threshold, Double_t robustFraction){
478 //
479 // Loops over all AliTPCCalROCs and performs a globalFit in each ROC
480 // AliTPCCalPad with fit-data is returned
481 // chi2Threshold: Threshold for chi2 when EvalRobust is called
482 // robustFraction: Fraction of data that will be used in EvalRobust
483 // chi2Threshold: Threshold for chi2 when EvalRobust is called
484 // robustFraction: Fraction of data that will be used in EvalRobust
485 //
486 AliTPCCalPad* pad = new AliTPCCalPad(padName, padName);
487 TVectorD fitParam(0);
488 TMatrixD covMatrix(0,0);
489 Float_t chi2 = 0;
490 for (Int_t isec = 0; isec < 72; isec++){
491 if (PadOutliers)
492 GetCalROC(isec)->GlobalFit(PadOutliers->GetCalROC(isec), robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction);
493 else
494 GetCalROC(isec)->GlobalFit(0, robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction);
495 pad->SetCalROC(AliTPCCalROC::CreateGlobalFitCalROC(fitParam, isec));
496 }
497 return pad;
498}
499
500
501void AliTPCCalPad::GlobalSidesFit(const AliTPCCalPad* PadOutliers, TVectorD &fitParamSideA, TVectorD &fitParamSideC,TMatrixD &covMatrixSideA, TMatrixD &covMatrixSideC, Float_t & chi2SideA, Float_t & chi2SideC, Int_t fitType, Bool_t robust, Double_t chi2Threshold, Double_t robustFraction){
502 //
503 // Makes a GlobalFit over each side and return fit-parameters, covariance and chi2 for each side
504 // fitType == 0: fit plane function
505 // fitType == 1: fit parabolic function
506 // PadOutliers - pads with value !=0 are not used in fitting procedure
507 // chi2Threshold: Threshold for chi2 when EvalRobust is called
508 // robustFraction: Fraction of data that will be used in EvalRobust
509 //
510 TLinearFitter* fitterGA = 0;
511 TLinearFitter* fitterGC = 0;
512
513 if (fitType == 1) {
514 fitterGA = new TLinearFitter (6,"x0++x1++x2++x3++x4++x5");
515 fitterGC = new TLinearFitter (6,"x0++x1++x2++x3++x4++x5");
516 }
517 else {
518 fitterGA = new TLinearFitter(3,"x0++x1++x2");
519 fitterGC = new TLinearFitter(3,"x0++x1++x2");
520 }
521 fitterGA->StoreData(kTRUE);
522 fitterGC->StoreData(kTRUE);
523 fitterGA->ClearPoints();
524 fitterGC->ClearPoints();
525 Double_t xx[6];
526 Int_t npointsA=0;
527 Int_t npointsC=0;
528
529 Float_t localXY[3] = {0}; // pad's position, needed to get the pad's position
530 Float_t lx, ly; // pads position
531
532 AliTPCROC* tpcROCinstance = AliTPCROC::Instance(); // to calculate the pad's position
533
534 // loop over all sectors and pads and read data into fitterGA and fitterGC
535 if (fitType == 1) {
536 // parabolic fit
537 fitParamSideA.ResizeTo(6);
538 fitParamSideC.ResizeTo(6);
539 covMatrixSideA.ResizeTo(6,6);
540 covMatrixSideC.ResizeTo(6,6);
541 for (UInt_t isec = 0; isec<72; isec++){
542 for (UInt_t irow = 0; irow < GetCalROC(isec)->GetNrows(); irow++) {
543 for (UInt_t ipad = 0; ipad < GetCalROC(isec)->GetNPads(irow); ipad++) {
544 // fill fitterG
545 tpcROCinstance->GetPositionLocal(isec, irow, ipad, localXY); // calculate position localXY by sector, pad and row number
546 lx = localXY[0];
547 ly = localXY[1];
548 xx[0] = 1;
549 xx[1] = lx;
550 xx[2] = ly;
551 xx[3] = lx*lx;
552 xx[4] = ly*ly;
553 xx[5] = lx*ly;
554 if (!PadOutliers || PadOutliers->GetCalROC(isec)->GetValue(irow, ipad) != 1) {
555 // if given pad is no outlier, add it to TLinearFitter, decide to which of both
556// sector 0 - 17: IROC, A
557// sector 18 - 35: IROC, C
558// sector 36 - 53: OROC, A
559// sector 54 - 71: CROC, C
560 if (isec <= 17 || (isec >= 36 && isec <= 53)) { // Side A
561 npointsA++;
562 fitterGA->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);
563 }
564 else { // side C
565 npointsC++;
566 fitterGC->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);
567 }
568 }
569 }
570 }
571 }
572 }
573 else {
574 // linear fit
575 fitParamSideA.ResizeTo(3);
576 fitParamSideC.ResizeTo(3);
577 covMatrixSideA.ResizeTo(3,3);
578 covMatrixSideC.ResizeTo(3,3);
579
580 for (UInt_t isec = 0; isec<72; isec++){
581 for (UInt_t irow = 0; irow < GetCalROC(isec)->GetNrows(); irow++) {
582 for (UInt_t ipad = 0; ipad < GetCalROC(isec)->GetNPads(irow); ipad++) {
583 // fill fitterG
584 tpcROCinstance->GetPositionLocal(isec, irow, ipad, localXY); // calculate position localXY by sector, pad and row number
585 lx = localXY[0];
586 ly = localXY[1];
587 xx[0] = 1;
588 xx[1] = lx;
589 xx[2] = ly;
590 if (!PadOutliers || PadOutliers->GetCalROC(isec)->GetValue(irow, ipad) != 1) {
591 // if given pad is no outlier, add it to TLinearFitter, decide to which of both
592// sector 0 - 17: IROC, A
593// sector 18 - 35: IROC, C
594// sector 36 - 53: OROC, A
595// sector 54 - 71: CROC, C
596 if (isec <= 17 || (isec >= 36 && isec <= 53)) {
597 // Side A
598 npointsA++;
599 fitterGA->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);
600 }
601 else {
602 // side C
603 npointsC++;
604 fitterGC->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);
605 }
606 }
607 }
608 }
609 }
610 }
611
612 fitterGA->Eval();
613 fitterGC->Eval();
614 fitterGA->GetParameters(fitParamSideA);
615 fitterGC->GetParameters(fitParamSideC);
616 fitterGA->GetCovarianceMatrix(covMatrixSideA);
617 fitterGC->GetCovarianceMatrix(covMatrixSideC);
618 if (fitType == 1){
619 chi2SideA = fitterGA->GetChisquare()/(npointsA-6.);
620 chi2SideC = fitterGC->GetChisquare()/(npointsC-6.);
621 }
622 else {
623 chi2SideA = fitterGA->GetChisquare()/(npointsA-3.);
624 chi2SideC = fitterGC->GetChisquare()/(npointsC-3.);
625 }
626 if (robust && chi2SideA > chi2Threshold) {
627 // std::cout << "robust fitter called... " << std::endl;
628 fitterGA->EvalRobust(robustFraction);
629 fitterGA->GetParameters(fitParamSideA);
630 }
631 if (robust && chi2SideC > chi2Threshold) {
632 // std::cout << "robust fitter called... " << std::endl;
633 fitterGC->EvalRobust(robustFraction);
634 fitterGC->GetParameters(fitParamSideC);
635 }
636 delete fitterGA;
637 delete fitterGC;
638}
200be8a6 639
586007f3 640
ca5dca67 641void AliTPCCalPad::MakeTree(const char * fileName, TObjArray * array, const char * mapFileName, AliTPCCalPad* outlierPad, Float_t ltmFraction) {
586007f3 642 //
a6d2bd0c 643 // Write a tree with all available information
644 // im mapFileName is speciefied, the Map information are also written to the tree
645 // pads specified in outlierPad are not used for calculating statistics
646 // - this function will be moved to AliTPCCalibViewer -
647 // - DO NOT USE THIS FUNCTION ANY MORE -
586007f3 648 //
586007f3 649 AliTPCROC* tpcROCinstance = AliTPCROC::Instance();
ca5dca67 650
651 TObjArray* mapIROCs = 0;
652 TObjArray* mapOROCs = 0;
653 TVectorF *mapIROCArray = 0;
654 TVectorF *mapOROCArray = 0;
655 Int_t mapEntries = 0;
656 TString* mapNames = 0;
657
658 if (mapFileName) {
659 TFile mapFile(mapFileName, "read");
660
661 TList* listOfROCs = mapFile.GetListOfKeys();
662 mapEntries = listOfROCs->GetEntries()/2;
663 mapIROCs = new TObjArray(mapEntries*2);
664 mapOROCs = new TObjArray(mapEntries*2);
665 mapIROCArray = new TVectorF[mapEntries];
666 mapOROCArray = new TVectorF[mapEntries];
667
668 mapNames = new TString[mapEntries];
669 for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
670 TString ROCname(((TKey*)(listOfROCs->At(ivalue*2)))->GetName());
671 ROCname.Remove(ROCname.Length()-4, 4);
672 mapIROCs->AddAt((AliTPCCalROC*)mapFile.Get((ROCname + "IROC").Data()), ivalue);
673 mapOROCs->AddAt((AliTPCCalROC*)mapFile.Get((ROCname + "OROC").Data()), ivalue);
674 mapNames[ivalue].Append(ROCname);
675 }
676
677 for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
678 mapIROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(0));
679 mapOROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(36));
680
681 for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(0); ichannel++)
682 (mapIROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapIROCs->At(ivalue)))->GetValue(ichannel);
683 for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(36); ichannel++)
684 (mapOROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapOROCs->At(ivalue)))->GetValue(ichannel);
685 }
686
687 } // if (mapFileName)
688
689 TTreeSRedirector cstream(fileName);
7fb8d357 690 Int_t arrayEntries = array->GetEntries();
586007f3 691
7fb8d357 692 TString* names = new TString[arrayEntries];
693 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
586007f3 694 names[ivalue].Append(((AliTPCCalPad*)array->At(ivalue))->GetName());
695
696 for (UInt_t isector = 0; isector < tpcROCinstance->GetNSectors(); isector++) {
697 //
698 // get statistic for given sector
699 //
7fb8d357 700 TVectorF median(arrayEntries);
701 TVectorF mean(arrayEntries);
702 TVectorF rms(arrayEntries);
703 TVectorF ltm(arrayEntries);
ca5dca67 704 TVectorF ltmrms(arrayEntries);
705 TVectorF medianWithOut(arrayEntries);
706 TVectorF meanWithOut(arrayEntries);
707 TVectorF rmsWithOut(arrayEntries);
708 TVectorF ltmWithOut(arrayEntries);
709 TVectorF ltmrmsWithOut(arrayEntries);
586007f3 710
7fb8d357 711 TVectorF *vectorArray = new TVectorF[arrayEntries];
712 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
586007f3 713 vectorArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
714
7fb8d357 715 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
586007f3 716 AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
717 AliTPCCalROC* calROC = calPad->GetCalROC(isector);
ca5dca67 718 AliTPCCalROC* outlierROC = 0;
719 if (outlierPad) outlierROC = outlierPad->GetCalROC(isector);
7fb8d357 720 if (calROC) {
721 median[ivalue] = calROC->GetMedian();
722 mean[ivalue] = calROC->GetMean();
723 rms[ivalue] = calROC->GetRMS();
ca5dca67 724 Double_t ltmrmsValue = 0;
725 ltm[ivalue] = calROC->GetLTM(&ltmrmsValue, ltmFraction);
726 ltmrms[ivalue] = ltmrmsValue;
727 if (outlierROC) {
728 medianWithOut[ivalue] = calROC->GetMedian(outlierROC);
729 meanWithOut[ivalue] = calROC->GetMean(outlierROC);
730 rmsWithOut[ivalue] = calROC->GetRMS(outlierROC);
731 ltmrmsValue = 0;
732 ltmWithOut[ivalue] = calROC->GetLTM(&ltmrmsValue, ltmFraction, outlierROC);
733 ltmrmsWithOut[ivalue] = ltmrmsValue;
734 }
7fb8d357 735 }
736 else {
737 median[ivalue] = 0.;
738 mean[ivalue] = 0.;
739 rms[ivalue] = 0.;
740 ltm[ivalue] = 0.;
ca5dca67 741 ltmrms[ivalue] = 0.;
742 medianWithOut[ivalue] = 0.;
743 meanWithOut[ivalue] = 0.;
744 rmsWithOut[ivalue] = 0.;
745 ltmWithOut[ivalue] = 0.;
746 ltmrmsWithOut[ivalue] = 0.;
7fb8d357 747 }
586007f3 748 }
749
750 //
751 // fill vectors of variable per pad
752 //
f62f5cfe 753 TVectorF *posArray = new TVectorF[8];
754 for (Int_t ivalue = 0; ivalue < 8; ivalue++)
586007f3 755 posArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
756
757 Float_t posG[3] = {0};
758 Float_t posL[3] = {0};
759 Int_t ichannel = 0;
760 for (UInt_t irow = 0; irow < tpcROCinstance->GetNRows(isector); irow++) {
761 for (UInt_t ipad = 0; ipad < tpcROCinstance->GetNPads(isector, irow); ipad++) {
762 tpcROCinstance->GetPositionLocal(isector, irow, ipad, posL);
763 tpcROCinstance->GetPositionGlobal(isector, irow, ipad, posG);
764 posArray[0][ichannel] = irow;
765 posArray[1][ichannel] = ipad;
766 posArray[2][ichannel] = posL[0];
767 posArray[3][ichannel] = posL[1];
768 posArray[4][ichannel] = posG[0];
769 posArray[5][ichannel] = posG[1];
f62f5cfe 770 posArray[6][ichannel] = (Int_t)(ipad - (Double_t)(tpcROCinstance->GetNPads(isector, irow))/2);
771 posArray[7][ichannel] = ichannel;
586007f3 772
773 // loop over array containing AliTPCCalPads
7fb8d357 774 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
586007f3 775 AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
7fb8d357 776 AliTPCCalROC* calROC = calPad->GetCalROC(isector);
777 if (calROC)
778 (vectorArray[ivalue])[ichannel] = calROC->GetValue(irow, ipad);
779 else
780 (vectorArray[ivalue])[ichannel] = 0;
586007f3 781 }
782 ichannel++;
783 }
784 }
785
786 cstream << "calPads" <<
787 "sector=" << isector;
788
7fb8d357 789 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
586007f3 790 cstream << "calPads" <<
7fb8d357 791 (Char_t*)((names[ivalue] + "_Median=").Data()) << median[ivalue] <<
792 (Char_t*)((names[ivalue] + "_Mean=").Data()) << mean[ivalue] <<
793 (Char_t*)((names[ivalue] + "_RMS=").Data()) << rms[ivalue] <<
ca5dca67 794 (Char_t*)((names[ivalue] + "_LTM=").Data()) << ltm[ivalue] <<
795 (Char_t*)((names[ivalue] + "_RMS_LTM=").Data()) << ltmrms[ivalue];
796 if (outlierPad) {
797 cstream << "calPads" <<
798 (Char_t*)((names[ivalue] + "_Median_OutlierCutted=").Data()) << medianWithOut[ivalue] <<
799 (Char_t*)((names[ivalue] + "_Mean_OutlierCutted=").Data()) << meanWithOut[ivalue] <<
800 (Char_t*)((names[ivalue] + "_RMS_OutlierCutted=").Data()) << rmsWithOut[ivalue] <<
801 (Char_t*)((names[ivalue] + "_LTM_OutlierCutted=").Data()) << ltmWithOut[ivalue] <<
802 (Char_t*)((names[ivalue] + "_RMS_LTM_OutlierCutted=").Data()) << ltmrmsWithOut[ivalue];
803 }
586007f3 804 }
805
7fb8d357 806 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
586007f3 807 cstream << "calPads" <<
808 (Char_t*)((names[ivalue] + ".=").Data()) << &vectorArray[ivalue];
809 }
810
ca5dca67 811 if (mapFileName) {
812 for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
813 if (isector < 36)
814 cstream << "calPads" <<
815 (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapIROCArray[ivalue];
816 else
817 cstream << "calPads" <<
818 (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapOROCArray[ivalue];
819 }
820 }
821
586007f3 822 cstream << "calPads" <<
823 "row.=" << &posArray[0] <<
824 "pad.=" << &posArray[1] <<
825 "lx.=" << &posArray[2] <<
826 "ly.=" << &posArray[3] <<
827 "gx.=" << &posArray[4] <<
f62f5cfe 828 "gy.=" << &posArray[5] <<
829 "rpad.=" << &posArray[6] <<
830 "channel.=" << &posArray[7];
586007f3 831
832 cstream << "calPads" <<
833 "\n";
834
835 delete[] posArray;
836 delete[] vectorArray;
837 }
838 delete[] names;
ca5dca67 839 if (mapFileName) {
840 delete mapIROCs;
841 delete mapOROCs;
842 delete[] mapIROCArray;
843 delete[] mapOROCArray;
844 delete[] mapNames;
845 }
586007f3 846}
847
848