1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // TPC calibration class for parameters which are saved per pad //
21 // Each AliTPCCalPad consists of 72 AliTPCCalROC-objects //
23 ///////////////////////////////////////////////////////////////////////////////
25 #include "AliTPCCalPad.h"
26 #include "AliTPCCalROC.h"
27 #include <TObjArray.h>
32 #include "TTreeStream.h"
37 #include <TObjString.h>
41 ClassImp(AliTPCCalPad)
43 //_____________________________________________________________________________
44 AliTPCCalPad::AliTPCCalPad():TNamed()
47 // AliTPCCalPad default constructor
50 for (Int_t isec = 0; isec < kNsec; isec++) {
56 //_____________________________________________________________________________
57 AliTPCCalPad::AliTPCCalPad(const Text_t *name, const Text_t *title)
61 // AliTPCCalPad constructor
63 for (Int_t isec = 0; isec < kNsec; isec++) {
64 fROC[isec] = new AliTPCCalROC(isec);
69 //_____________________________________________________________________________
70 AliTPCCalPad::AliTPCCalPad(const AliTPCCalPad &c):TNamed(c)
73 // AliTPCCalPad copy constructor
76 for (Int_t isec = 0; isec < kNsec; isec++) {
79 fROC[isec] = new AliTPCCalROC(*(c.fROC[isec]));
83 //_____________________________________________________________________________
84 AliTPCCalPad::AliTPCCalPad(TObjArray * array):TNamed()
87 // AliTPCCalPad default constructor
90 for (Int_t isec = 0; isec < kNsec; isec++) {
91 fROC[isec] = (AliTPCCalROC *)array->At(isec);
97 ///_____________________________________________________________________________
98 AliTPCCalPad::~AliTPCCalPad()
101 // AliTPCCalPad destructor
104 for (Int_t isec = 0; isec < kNsec; isec++) {
113 //_____________________________________________________________________________
114 AliTPCCalPad &AliTPCCalPad::operator=(const AliTPCCalPad &c)
117 // Assignment operator
120 if (this != &c) ((AliTPCCalPad &) c).Copy(*this);
125 //_____________________________________________________________________________
126 void AliTPCCalPad::Copy(TObject &c) const
132 for (Int_t isec = 0; isec < kNsec; isec++) {
134 fROC[isec]->Copy(*((AliTPCCalPad &) c).fROC[isec]);
141 void AliTPCCalPad::SetCalROC(AliTPCCalROC* roc, Int_t sector){
143 // Set AliTPCCalROC copies values from 'roc'
144 // if sector == -1 the sector specified in 'roc' is used
145 // else sector specified in 'roc' is ignored and specified sector is filled
147 if (sector == -1) sector = roc->GetSector();
148 if (!fROC[sector]) fROC[sector] = new AliTPCCalROC(sector);
149 for (UInt_t ichannel = 0; ichannel < roc->GetNchannels(); ichannel++)
150 fROC[sector]->SetValue(ichannel, roc->GetValue(ichannel));
155 //_____________________________________________________________________________
156 void AliTPCCalPad::Add(Float_t c1)
159 // add constant c1 to all channels of all ROCs
162 for (Int_t isec = 0; isec < kNsec; isec++) {
169 //_____________________________________________________________________________
170 void AliTPCCalPad::Multiply(Float_t c1)
173 // multiply each channel of all ROCs with c1
175 for (Int_t isec = 0; isec < kNsec; isec++) {
177 fROC[isec]->Multiply(c1);
182 //_____________________________________________________________________________
183 void AliTPCCalPad::Add(const AliTPCCalPad * pad, Double_t c1)
186 // multiply AliTPCCalPad 'pad' by c1 and add each channel to the coresponing channel in all ROCs
189 for (Int_t isec = 0; isec < kNsec; isec++) {
190 if (fROC[isec] && pad->GetCalROC(isec)){
191 fROC[isec]->Add(pad->GetCalROC(isec),c1);
196 //_____________________________________________________________________________
197 void AliTPCCalPad::Multiply(const AliTPCCalPad * pad)
200 // multiply each channel of all ROCs with the coresponding channel of 'pad'
203 for (Int_t isec = 0; isec < kNsec; isec++) {
205 fROC[isec]->Multiply(pad->GetCalROC(isec));
210 //_____________________________________________________________________________
211 void AliTPCCalPad::Divide(const AliTPCCalPad * pad)
214 // divide each channel of all ROCs by the coresponding channel of 'pad'
217 for (Int_t isec = 0; isec < kNsec; isec++) {
219 fROC[isec]->Divide(pad->GetCalROC(isec));
224 //_____________________________________________________________________________
225 TGraph * AliTPCCalPad::MakeGraph(Int_t type, Float_t ratio){
232 for (Int_t i=0;i<72;i++) if (fROC[i]) npoints++;
233 TGraph * graph = new TGraph(npoints);
235 for (Int_t isec=0;isec<72;isec++){
236 if (!fROC[isec]) continue;
237 if (type==0) graph->SetPoint(npoints,isec,fROC[isec]->GetMean());
238 if (type==1) graph->SetPoint(npoints,isec,fROC[isec]->GetMedian());
239 if (type==2) graph->SetPoint(npoints,isec,fROC[isec]->GetLTM(0,ratio));
243 graph->GetXaxis()->SetTitle("Sector");
245 graph->GetYaxis()->SetTitle("Mean");
246 graph->SetMarkerStyle(22);
249 graph->GetYaxis()->SetTitle("Median");
250 graph->SetMarkerStyle(22);
253 graph->GetYaxis()->SetTitle(Form("Mean%f",ratio));
254 graph->SetMarkerStyle(24);
260 //_____________________________________________________________________________
261 Double_t AliTPCCalPad::GetMeanRMS(Double_t &rms)
264 // Calculates mean and RMS of all ROCs
266 Double_t sum = 0, sum2 = 0, n=0, val=0;
267 for (Int_t isec = 0; isec < kNsec; isec++) {
268 AliTPCCalROC *calRoc = fROC[isec];
270 for (UInt_t irow=0; irow<calRoc->GetNrows(); irow++){
271 for (UInt_t ipad=0; ipad<calRoc->GetNPads(irow); ipad++){
272 val = calRoc->GetValue(irow,ipad);
282 Double_t mean = sum*n1;
283 rms = TMath::Sqrt(TMath::Abs(sum2*n1-mean*mean));
288 //_____________________________________________________________________________
289 Double_t AliTPCCalPad::GetMean(AliTPCCalPad* outlierPad)
292 // return mean of the mean of all ROCs
296 for (Int_t isec = 0; isec < kNsec; isec++) {
297 AliTPCCalROC *calRoc = fROC[isec];
299 AliTPCCalROC* outlierROC = 0;
300 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
301 arr[n] = calRoc->GetMean(outlierROC);
305 return TMath::Mean(n,arr);
308 //_____________________________________________________________________________
309 Double_t AliTPCCalPad::GetRMS(AliTPCCalPad* outlierPad)
312 // return mean of the RMS of all ROCs
316 for (Int_t isec = 0; isec < kNsec; isec++) {
317 AliTPCCalROC *calRoc = fROC[isec];
319 AliTPCCalROC* outlierROC = 0;
320 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
321 arr[n] = calRoc->GetRMS(outlierROC);
325 return TMath::Mean(n,arr);
328 //_____________________________________________________________________________
329 Double_t AliTPCCalPad::GetMedian(AliTPCCalPad* outlierPad)
332 // return mean of the median of all ROCs
336 for (Int_t isec = 0; isec < kNsec; isec++) {
337 AliTPCCalROC *calRoc = fROC[isec];
339 AliTPCCalROC* outlierROC = 0;
340 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
341 arr[n] = calRoc->GetMedian(outlierROC);
345 return TMath::Mean(n,arr);
348 //_____________________________________________________________________________
349 Double_t AliTPCCalPad::GetLTM(Double_t *sigma, Double_t fraction, AliTPCCalPad* outlierPad)
352 // return mean of the LTM and sigma of all ROCs
354 Double_t arrm[kNsec];
355 Double_t arrs[kNsec];
359 for (Int_t isec = 0; isec < kNsec; isec++) {
360 AliTPCCalROC *calRoc = fROC[isec];
362 if ( sigma ) sTemp=arrs+n;
363 AliTPCCalROC* outlierROC = 0;
364 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
365 arrm[n] = calRoc->GetLTM(sTemp,fraction, outlierROC);
369 if ( sigma ) *sigma = TMath::Mean(n,arrs);
370 return TMath::Mean(n,arrm);
373 //_____________________________________________________________________________
374 TH1F * AliTPCCalPad::MakeHisto1D(Float_t min, Float_t max,Int_t type){
377 // type -1 = user defined range
378 // 0 = nsigma cut nsigma=min
383 Float_t mean = GetMean();
384 Float_t sigma = GetRMS();
385 Float_t nsigma = TMath::Abs(min);
386 min = mean-nsigma*sigma;
387 max = mean+nsigma*sigma;
391 Float_t mean = GetMedian();
398 // LTM mean +- nsigma
401 Float_t mean = GetLTM(&sigma,max);
408 sprintf(name,"%s Pad 1D",GetTitle());
409 TH1F * his = new TH1F(name,name,100, min,max);
410 for (Int_t isec = 0; isec < kNsec; isec++) {
412 for (UInt_t irow=0; irow<fROC[isec]->GetNrows(); irow++){
413 UInt_t npads = (Int_t)fROC[isec]->GetNPads(irow);
414 for (UInt_t ipad=0; ipad<npads; ipad++){
415 his->Fill(fROC[isec]->GetValue(irow,ipad));
423 //_____________________________________________________________________________
424 TH2F *AliTPCCalPad::MakeHisto2D(Int_t side){
427 // side - specify the side A = 0 C = 1
428 // type - used types of determination of boundaries in z
430 Float_t kEpsilon = 0.000000000001;
431 TH2F * his = new TH2F(GetName(), GetName(), 250,-250,250,250,-250,250);
432 AliTPCROC * roc = AliTPCROC::Instance();
433 for (Int_t isec=0; isec<72; isec++){
434 if (side==0 && isec%36>=18) continue;
435 if (side>0 && isec%36<18) continue;
437 AliTPCCalROC * calRoc = fROC[isec];
438 for (UInt_t irow=0; irow<calRoc->GetNrows(); irow++)
439 for (UInt_t ipad=0; ipad<calRoc->GetNPads(irow); ipad++)
440 if (TMath::Abs(calRoc->GetValue(irow,ipad))>kEpsilon){
442 roc->GetPositionGlobal(isec,irow,ipad,xyz);
443 Int_t binx = 1+TMath::Nint((xyz[0]+250.)*0.5);
444 Int_t biny = 1+TMath::Nint((xyz[1]+250.)*0.5);
445 Float_t value = calRoc->GetValue(irow,ipad);
446 his->SetBinContent(binx,biny,value);
450 his->SetXTitle("x (cm)");
451 his->SetYTitle("y (cm)");
456 AliTPCCalPad* AliTPCCalPad::LocalFit(const char* padName, Int_t rowRadius, Int_t padRadius, AliTPCCalPad* PadOutliers, Bool_t robust, Double_t chi2Threshold, Double_t robustFraction, Bool_t printCurrentSector) const {
458 // Loops over all AliTPCCalROCs and performs a localFit in each ROC
459 // AliTPCCalPad with fit-data is returned
460 // rowRadius and padRadius specifies a window around a given pad in one sector.
461 // The data of this window are fitted with a parabolic function.
462 // This function is evaluated at the pad's position.
463 // At the edges the window is shifted, so that the specified pad is not anymore in the center of the window.
464 // rowRadius - radius - rows to be used for smoothing
465 // padradius - radius - pads to be used for smoothing
466 // ROCoutlier - map of outliers - pads not to be used for local smoothing
467 // robust - robust method of fitting - (much slower)
468 // chi2Threshold: Threshold for chi2 when EvalRobust is called
469 // robustFraction: Fraction of data that will be used in EvalRobust
472 AliTPCCalPad* pad = new AliTPCCalPad(padName, padName);
473 for (Int_t isec = 0; isec < 72; isec++){
474 if (printCurrentSector) std::cout << "LocalFit in sector " << isec << "\r" << std::flush;
476 pad->SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, PadOutliers->GetCalROC(isec), robust, chi2Threshold, robustFraction));
478 pad->SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, 0, robust, chi2Threshold, robustFraction));
484 AliTPCCalPad* AliTPCCalPad::GlobalFit(const char* padName, AliTPCCalPad* PadOutliers, Bool_t robust, Int_t fitType, Double_t chi2Threshold, Double_t robustFraction, Double_t err, TObjArray *fitParArr, TObjArray *fitCovArr){
486 // Loops over all AliTPCCalROCs and performs a globalFit in each ROC
487 // AliTPCCalPad with fit-data is returned
488 // chi2Threshold: Threshold for chi2 when EvalRobust is called
489 // robustFraction: Fraction of data that will be used in EvalRobust
490 // chi2Threshold: Threshold for chi2 when EvalRobust is called
491 // robustFraction: Fraction of data that will be used in EvalRobust
492 // err: error of the data points
493 // if fitParArr and/or fitCovArr is given, write fitParameters and/or covariance Matrices into the array
495 AliTPCCalPad* pad = new AliTPCCalPad(padName, padName);
496 TVectorD fitParam(0);
497 TMatrixD covMatrix(0,0);
499 for (Int_t isec = 0; isec < 72; isec++){
501 GetCalROC(isec)->GlobalFit(PadOutliers->GetCalROC(isec), robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction, err);
503 GetCalROC(isec)->GlobalFit(0, robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction, err);
505 AliTPCCalROC *roc=AliTPCCalROC::CreateGlobalFitCalROC(fitParam, isec);
508 if ( fitParArr ) fitParArr->AddAtAndExpand(new TVectorD(fitParam), isec);
509 if ( fitCovArr ) fitCovArr->AddAtAndExpand(new TMatrixD(covMatrix), isec);
513 //_____________________________________________________________________________
514 TObjArray* AliTPCCalPad::CreateFormulaArray(const char *fitFormula)
517 // create an array of TFormulas for the each parameter of the fit function
520 // split fit string in single parameters
521 // find dimension of the fit:
522 TString fitString(fitFormula);
523 fitString.ReplaceAll("++","#");
524 fitString.ReplaceAll(" ","");
525 TObjArray *arrFitParams = fitString.Tokenize("#");
526 Int_t ndim = arrFitParams->GetEntries();
527 //create array of TFormulas to evaluate the parameters
528 TObjArray *arrFitFormulas = new TObjArray(ndim);
529 arrFitFormulas->SetOwner(kTRUE);
530 for (Int_t idim=0;idim<ndim;++idim){
531 TString s=((TObjString*)arrFitParams->At(idim))->GetString();
532 s.ReplaceAll("gx","[0]");
533 s.ReplaceAll("gy","[1]");
534 s.ReplaceAll("lx","[2]");
535 s.ReplaceAll("ly","[3]");
536 s.ReplaceAll("sector","[4]");
537 arrFitFormulas->AddAt(new TFormula(Form("param%02d",idim),s.Data()),idim);
541 return arrFitFormulas;
543 //_____________________________________________________________________________
544 void AliTPCCalPad::EvalFormulaArray(const TObjArray &arrFitFormulas, TVectorD &results,
545 const Int_t sec, const Int_t row, const Int_t pad)
548 // evaluate the fit formulas
550 Int_t ndim=arrFitFormulas.GetEntries();
551 results.ResizeTo(ndim);
553 AliTPCROC* tpcROCinstance = AliTPCROC::Instance(); // to calculate the pad's position
555 Float_t globalXYZ[3];
556 tpcROCinstance->GetPositionLocal(sec, row, pad, localXYZ);
557 tpcROCinstance->GetPositionGlobal(sec, row, pad, globalXYZ);
558 //calculate parameter values
559 for (Int_t idim=0;idim<ndim;++idim){
560 TFormula *f=(TFormula*)arrFitFormulas.At(idim);
561 f->SetParameters(globalXYZ[0],globalXYZ[1],localXYZ[0],localXYZ[1],sec);
562 results[idim]=f->Eval(0);
565 //_____________________________________________________________________________
566 void AliTPCCalPad::GlobalSidesFit(const AliTPCCalPad* PadOutliers, const char* fitFormula, TVectorD &fitParamSideA, TVectorD &fitParamSideC,TMatrixD &covMatrixSideA, TMatrixD &covMatrixSideC, Float_t & chi2SideA, Float_t & chi2SideC, AliTPCCalPad *pointError, Bool_t robust, Double_t robustFraction){
568 // Performs a fit on both sides.
569 // Valid information for the fitFormula are the variables
570 // - gx, gy, lx ,ly: meaning global x, global y, local x, local y value of the padName
571 // - sector: the sector number.
572 // eg. a formula might look 'gy' or '(sector<36) ++ gy' or 'gx ++ gy' or 'gx ++ gy ++ lx ++ lx^2' and so on
574 // PadOutliers - pads with value !=0 are not used in fitting procedure
575 // chi2Threshold: Threshold for chi2 when EvalRobust is called
576 // robustFraction: Fraction of data that will be used in EvalRobust
579 TObjArray* arrFitFormulas=CreateFormulaArray(fitFormula);
580 Int_t ndim = arrFitFormulas->GetEntries();
581 //resize output data arrays
582 fitParamSideA.ResizeTo(ndim+1);
583 fitParamSideC.ResizeTo(ndim+1);
584 covMatrixSideA.ResizeTo(ndim+1,ndim+1);
585 covMatrixSideC.ResizeTo(ndim+1,ndim+1);
586 // create linear fitter for A- and C- Side
587 TLinearFitter* fitterGA = new TLinearFitter(ndim+1,Form("hyp%d",ndim));
588 TLinearFitter* fitterGC = new TLinearFitter(ndim+1,Form("hyp%d",ndim));
589 fitterGA->StoreData(kTRUE);
590 fitterGC->StoreData(kTRUE);
592 TVectorD parValues(ndim);
594 AliTPCCalROC *rocErr=0x0;
596 for (UInt_t isec = 0; isec<kNsec; ++isec){
597 AliTPCCalROC *rocOut=PadOutliers->GetCalROC(isec);
598 AliTPCCalROC *rocData=GetCalROC(isec);
599 if (pointError) rocErr=pointError->GetCalROC(isec);
600 if (!rocData) continue;
601 for (UInt_t irow = 0; irow < GetCalROC(isec)->GetNrows(); irow++) {
602 for (UInt_t ipad = 0; ipad < GetCalROC(isec)->GetNPads(irow); ipad++) {
604 if (rocOut && rocOut->GetValue(irow,ipad)) continue;
605 //calculate parameter values
606 EvalFormulaArray(*arrFitFormulas,parValues,isec,irow,ipad);
608 Float_t value=rocData->GetValue(irow,ipad);
612 err=rocErr->GetValue(irow,ipad);
615 //add points to the fitters
617 fitterGA->AddPoint(parValues.GetMatrixArray(),value,err);
619 fitterGC->AddPoint(parValues.GetMatrixArray(),value,err);
625 fitterGA->EvalRobust(robustFraction);
626 fitterGC->EvalRobust(robustFraction);
631 chi2SideA=fitterGA->GetChisquare()/(fitterGA->GetNpoints()-(ndim+1));
632 chi2SideC=fitterGC->GetChisquare()/(fitterGC->GetNpoints()-(ndim+1));
633 fitterGA->GetParameters(fitParamSideA);
634 fitterGC->GetParameters(fitParamSideC);
635 fitterGA->GetCovarianceMatrix(covMatrixSideA);
636 fitterGC->GetCovarianceMatrix(covMatrixSideC);
638 delete arrFitFormulas;
644 AliTPCCalPad *AliTPCCalPad::CreateCalPadFit(const char* fitFormula, const TVectorD &fitParamSideA, const TVectorD &fitParamSideC)
649 TObjArray *arrFitFormulas=CreateFormulaArray(fitFormula);
650 Int_t ndim = arrFitFormulas->GetEntries();
651 //check if dimension of fit formula and fit parameters agree
652 if (ndim!=fitParamSideA.GetNrows()||ndim!=fitParamSideC.GetNrows()){
653 printf("AliTPCCalPad::CreateCalPadFit: Dimensions of fit formula and fit Parameters does not match!");
657 AliTPCCalPad *pad=new AliTPCCalPad("fitResultPad",Form("Fit result: %s",fitFormula));
658 //fill cal pad with fit results if requested
659 for (UInt_t isec = 0; isec<kNsec; ++isec){
660 AliTPCCalROC *roc=pad->GetCalROC(isec);
661 for (UInt_t irow = 0; irow < roc->GetNrows(); irow++) {
662 for (UInt_t ipad = 0; ipad < roc->GetNPads(irow); ipad++) {
663 const TVectorD *fitPar=0;
664 TVectorD fitResArray;
666 fitPar=&fitParamSideA;
668 fitPar=&fitParamSideC;
670 EvalFormulaArray(*arrFitFormulas,fitResArray, isec, irow, ipad);
671 for (Int_t idim=0;idim<ndim;++idim)
672 fitResArray(idim)*=(*fitPar)(idim);
673 roc->SetValue(irow,ipad,fitResArray.Sum());
677 delete arrFitFormulas;
681 void 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){
683 // Makes a GlobalFit over each side and return fit-parameters, covariance and chi2 for each side
684 // fitType == 0: fit plane function
685 // fitType == 1: fit parabolic function
686 // PadOutliers - pads with value !=0 are not used in fitting procedure
687 // chi2Threshold: Threshold for chi2 when EvalRobust is called
688 // robustFraction: Fraction of data that will be used in EvalRobust
690 TLinearFitter* fitterGA = 0;
691 TLinearFitter* fitterGC = 0;
694 fitterGA = new TLinearFitter (6,"x0++x1++x2++x3++x4++x5");
695 fitterGC = new TLinearFitter (6,"x0++x1++x2++x3++x4++x5");
698 fitterGA = new TLinearFitter(3,"x0++x1++x2");
699 fitterGC = new TLinearFitter(3,"x0++x1++x2");
701 fitterGA->StoreData(kTRUE);
702 fitterGC->StoreData(kTRUE);
703 fitterGA->ClearPoints();
704 fitterGC->ClearPoints();
709 Float_t localXY[3] = {0}; // pad's position, needed to get the pad's position
710 Float_t lx, ly; // pads position
712 AliTPCROC* tpcROCinstance = AliTPCROC::Instance(); // to calculate the pad's position
714 // loop over all sectors and pads and read data into fitterGA and fitterGC
717 fitParamSideA.ResizeTo(6);
718 fitParamSideC.ResizeTo(6);
719 covMatrixSideA.ResizeTo(6,6);
720 covMatrixSideC.ResizeTo(6,6);
721 for (UInt_t isec = 0; isec<72; isec++){
722 for (UInt_t irow = 0; irow < GetCalROC(isec)->GetNrows(); irow++) {
723 for (UInt_t ipad = 0; ipad < GetCalROC(isec)->GetNPads(irow); ipad++) {
725 tpcROCinstance->GetPositionLocal(isec, irow, ipad, localXY); // calculate position localXY by sector, pad and row number
734 if (!PadOutliers || PadOutliers->GetCalROC(isec)->GetValue(irow, ipad) != 1) {
735 // if given pad is no outlier, add it to TLinearFitter, decide to which of both
736 // sector 0 - 17: IROC, A
737 // sector 18 - 35: IROC, C
738 // sector 36 - 53: OROC, A
739 // sector 54 - 71: CROC, C
740 if (isec <= 17 || (isec >= 36 && isec <= 53)) { // Side A
742 fitterGA->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);
746 fitterGC->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);
755 fitParamSideA.ResizeTo(3);
756 fitParamSideC.ResizeTo(3);
757 covMatrixSideA.ResizeTo(3,3);
758 covMatrixSideC.ResizeTo(3,3);
760 for (UInt_t isec = 0; isec<72; isec++){
761 for (UInt_t irow = 0; irow < GetCalROC(isec)->GetNrows(); irow++) {
762 for (UInt_t ipad = 0; ipad < GetCalROC(isec)->GetNPads(irow); ipad++) {
764 tpcROCinstance->GetPositionLocal(isec, irow, ipad, localXY); // calculate position localXY by sector, pad and row number
770 if (!PadOutliers || PadOutliers->GetCalROC(isec)->GetValue(irow, ipad) != 1) {
771 // if given pad is no outlier, add it to TLinearFitter, decide to which of both
772 // sector 0 - 17: IROC, A
773 // sector 18 - 35: IROC, C
774 // sector 36 - 53: OROC, A
775 // sector 54 - 71: CROC, C
776 if (isec <= 17 || (isec >= 36 && isec <= 53)) {
779 fitterGA->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);
784 fitterGC->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);
794 fitterGA->GetParameters(fitParamSideA);
795 fitterGC->GetParameters(fitParamSideC);
796 fitterGA->GetCovarianceMatrix(covMatrixSideA);
797 fitterGC->GetCovarianceMatrix(covMatrixSideC);
799 chi2SideA = fitterGA->GetChisquare()/(npointsA-6.);
800 chi2SideC = fitterGC->GetChisquare()/(npointsC-6.);
803 chi2SideA = fitterGA->GetChisquare()/(npointsA-3.);
804 chi2SideC = fitterGC->GetChisquare()/(npointsC-3.);
806 if (robust && chi2SideA > chi2Threshold) {
807 // std::cout << "robust fitter called... " << std::endl;
808 fitterGA->EvalRobust(robustFraction);
809 fitterGA->GetParameters(fitParamSideA);
811 if (robust && chi2SideC > chi2Threshold) {
812 // std::cout << "robust fitter called... " << std::endl;
813 fitterGC->EvalRobust(robustFraction);
814 fitterGC->GetParameters(fitParamSideC);