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 **************************************************************************/
17 ///////////////////////////////////////////////////////////////////////////////
19 // Calibration base class for a single ROC //
20 // Contains one float value per pad //
21 // mapping of the pads taken form AliTPCROC //
23 ///////////////////////////////////////////////////////////////////////////////
36 #include "AliTPCCalROC.h"
37 #include "AliMathBase.h"
39 #include "TRandom3.h" // only needed by the AliTPCCalROCTest() method
41 ClassImp(AliTPCCalROC)
44 //_____________________________________________________________________________
45 AliTPCCalROC::AliTPCCalROC()
54 // Default constructor
59 //_____________________________________________________________________________
60 AliTPCCalROC::AliTPCCalROC(UInt_t sector)
69 // Constructor that initializes a given sector
72 fNChannels = AliTPCROC::Instance()->GetNChannels(fSector);
73 fNRows = AliTPCROC::Instance()->GetNRows(fSector);
74 fIndexes = AliTPCROC::Instance()->GetRowIndexes(fSector);
75 fData = new Float_t[fNChannels];
76 for (UInt_t idata = 0; idata< fNChannels; idata++) fData[idata] = 0.;
79 //_____________________________________________________________________________
80 AliTPCCalROC::AliTPCCalROC(const AliTPCCalROC &c)
89 // AliTPCCalROC copy constructor
92 fNChannels = AliTPCROC::Instance()->GetNChannels(fSector);
93 fNRows = AliTPCROC::Instance()->GetNRows(fSector);
94 fIndexes = AliTPCROC::Instance()->GetRowIndexes(fSector);
96 fData = new Float_t[fNChannels];
97 for (UInt_t idata = 0; idata< fNChannels; idata++) fData[idata] = c.fData[idata];
99 //____________________________________________________________________________
100 AliTPCCalROC & AliTPCCalROC::operator =(const AliTPCCalROC & param)
103 // assignment operator - dummy
110 //_____________________________________________________________________________
111 AliTPCCalROC::~AliTPCCalROC()
114 // AliTPCCalROC destructor
124 void AliTPCCalROC::Streamer(TBuffer &R__b)
127 // Stream an object of class AliTPCCalROC.
129 if (R__b.IsReading()) {
130 AliTPCCalROC::Class()->ReadBuffer(R__b, this);
131 fIndexes = AliTPCROC::Instance()->GetRowIndexes(fSector);
133 AliTPCCalROC::Class()->WriteBuffer(R__b,this);
140 void AliTPCCalROC::Add(Float_t c1){
142 // add c1 to each channel of the ROC
144 for (UInt_t idata = 0; idata< fNChannels; idata++) fData[idata]+=c1;
148 void AliTPCCalROC::Multiply(Float_t c1){
150 // multiply each channel of the ROC with c1
152 for (UInt_t idata = 0; idata< fNChannels; idata++) fData[idata]*=c1;
156 void AliTPCCalROC::Add(const AliTPCCalROC * roc, Double_t c1){
158 // multiply AliTPCCalROC roc by c1 and add each channel to the coresponing channel in the ROC
161 for (UInt_t idata = 0; idata< fNChannels; idata++){
162 fData[idata]+=roc->fData[idata]*c1;
167 void AliTPCCalROC::Multiply(const AliTPCCalROC* roc) {
169 // multiply each channel of the ROC with the coresponding channel of 'roc'
172 for (UInt_t idata = 0; idata< fNChannels; idata++){
173 fData[idata]*=roc->fData[idata];
178 void AliTPCCalROC::Divide(const AliTPCCalROC* roc) {
180 // divide each channel of the ROC by the coresponding channel of 'roc'
183 Float_t kEpsilon=0.00000000000000001;
184 for (UInt_t idata = 0; idata< fNChannels; idata++){
185 if (TMath::Abs(roc->fData[idata])>kEpsilon)
186 fData[idata]/=roc->fData[idata];
190 Double_t AliTPCCalROC::GetMean(AliTPCCalROC* outlierROC) {
192 // returns the mean value of the ROC
193 // pads with value != 0 in outlierROC are not used for the calculation
195 if (!outlierROC) return TMath::Mean(fNChannels, fData);
196 Double_t *ddata = new Double_t[fNChannels];
198 for (UInt_t i=0;i<fNChannels;i++) {
199 if (!(outlierROC->GetValue(i))) {
200 ddata[NPoints]= fData[NPoints];
204 Double_t mean = TMath::Mean(NPoints, ddata);
209 Double_t AliTPCCalROC::GetMedian(AliTPCCalROC* outlierROC) {
211 // returns the median value of the ROC
212 // pads with value != 0 in outlierROC are not used for the calculation
214 if (!outlierROC) return TMath::Median(fNChannels, fData);
215 Double_t *ddata = new Double_t[fNChannels];
217 for (UInt_t i=0;i<fNChannels;i++) {
218 if (!(outlierROC->GetValue(i))) {
219 ddata[NPoints]= fData[NPoints];
223 Double_t mean = TMath::Median(NPoints, ddata);
228 Double_t AliTPCCalROC::GetRMS(AliTPCCalROC* outlierROC) {
230 // returns the RMS value of the ROC
231 // pads with value != 0 in outlierROC are not used for the calculation
233 if (!outlierROC) return TMath::RMS(fNChannels, fData);
234 Double_t *ddata = new Double_t[fNChannels];
236 for (UInt_t i=0;i<fNChannels;i++) {
237 if (!(outlierROC->GetValue(i))) {
238 ddata[NPoints]= fData[NPoints];
242 Double_t mean = TMath::RMS(NPoints, ddata);
247 Double_t AliTPCCalROC::GetLTM(Double_t *sigma, Double_t fraction, AliTPCCalROC* outlierROC){
249 // returns the LTM and sigma
250 // pads with value != 0 in outlierROC are not used for the calculation
252 Double_t *ddata = new Double_t[fNChannels];
253 Double_t mean=0, lsigma=0;
255 for (UInt_t i=0;i<fNChannels;i++) {
256 if (!outlierROC || !(outlierROC->GetValue(i))) {
257 ddata[NPoints]= fData[NPoints];
261 Int_t hh = TMath::Min(TMath::Nint(fraction *NPoints), Int_t(NPoints));
262 AliMathBase::EvaluateUni(NPoints,ddata, mean, lsigma, hh);
263 if (sigma) *sigma=lsigma;
268 TH1F * AliTPCCalROC::MakeHisto1D(Float_t min, Float_t max,Int_t type){
271 // type -1 = user defined range
272 // 0 = nsigma cut nsigma=min
273 // 1 = delta cut around median delta=min
278 Float_t mean = GetMean();
279 Float_t sigma = GetRMS();
280 Float_t nsigma = TMath::Abs(min);
281 min = mean-nsigma*sigma;
282 max = mean+nsigma*sigma;
286 Float_t mean = GetMedian();
293 // LTM mean +- nsigma
296 Float_t mean = GetLTM(&sigma,max);
303 sprintf(name,"%s ROC 1D%d",GetTitle(),fSector);
304 TH1F * his = new TH1F(name,name,100, min,max);
305 for (UInt_t irow=0; irow<fNRows; irow++){
306 UInt_t npads = (Int_t)GetNPads(irow);
307 for (UInt_t ipad=0; ipad<npads; ipad++){
308 his->Fill(GetValue(irow,ipad));
316 TH2F * AliTPCCalROC::MakeHisto2D(Float_t min, Float_t max,Int_t type){
319 // type -1 = user defined range
320 // 0 = nsigma cut nsigma=min
321 // 1 = delta cut around median delta=min
326 Float_t mean = GetMean();
327 Float_t sigma = GetRMS();
328 Float_t nsigma = TMath::Abs(min);
329 min = mean-nsigma*sigma;
330 max = mean+nsigma*sigma;
334 Float_t mean = GetMedian();
341 Float_t mean = GetLTM(&sigma,max);
349 for (UInt_t irow=0; irow<fNRows; irow++){
350 if (GetNPads(irow)>maxPad) maxPad = GetNPads(irow);
353 sprintf(name,"%s ROC%d",GetTitle(),fSector);
354 TH2F * his = new TH2F(name,name,fNRows+10,-5, fNRows+5, maxPad+10, -(Int_t(maxPad/2))-5, maxPad/2+5);
355 for (UInt_t irow=0; irow<fNRows; irow++){
356 UInt_t npads = (Int_t)GetNPads(irow);
357 for (UInt_t ipad=0; ipad<npads; ipad++){
358 his->Fill(irow+0.5,Int_t(ipad)-Int_t(npads/2)+0.5,GetValue(irow,ipad));
361 his->SetMaximum(max);
362 his->SetMinimum(min);
366 TH2F * AliTPCCalROC::MakeHistoOutliers(Float_t delta, Float_t fraction, Int_t type){
368 // Make Histogram with outliers
369 // mode = 0 - sigma cut used
370 // mode = 1 - absolute cut used
371 // fraction - fraction of values used to define sigma
372 // delta - in mode 0 - nsigma cut
373 // mode 1 - delta cut
376 Float_t mean = GetLTM(&sigma,fraction);
377 if (type==0) delta*=sigma;
379 for (UInt_t irow=0; irow<fNRows; irow++){
380 if (GetNPads(irow)>maxPad) maxPad = GetNPads(irow);
384 sprintf(name,"%s ROC Outliers%d",GetTitle(),fSector);
385 TH2F * his = new TH2F(name,name,fNRows+10,-5, fNRows+5, maxPad+10, -(Int_t(maxPad/2))-5, maxPad/2+5);
386 for (UInt_t irow=0; irow<fNRows; irow++){
387 UInt_t npads = (Int_t)GetNPads(irow);
388 for (UInt_t ipad=0; ipad<npads; ipad++){
389 if (TMath::Abs(GetValue(irow,ipad)-mean)>delta)
390 his->Fill(irow+0.5,Int_t(ipad)-Int_t(npads/2)+0.5,1);
398 void AliTPCCalROC::Draw(Option_t* opt){
400 // create histogram with values and draw it
405 if (option.Contains("1D")){
418 void AliTPCCalROC::Test() {
420 // example function to show functionality and test AliTPCCalROC
423 Float_t kEpsilon=0.00001;
425 // create CalROC with defined values
426 AliTPCCalROC roc0(0);
427 for (UInt_t irow = 0; irow <roc0.GetNrows(); irow++){
428 for (UInt_t ipad = 0; ipad <roc0.GetNPads(irow); ipad++){
429 Float_t value = irow+ipad/1000.;
430 roc0.SetValue(irow,ipad,value);
434 // copy CalROC, readout values and compare to original
435 AliTPCCalROC roc1(roc0);
436 for (UInt_t irow = 0; irow <roc1.GetNrows(); irow++){
437 for (UInt_t ipad = 0; ipad <roc1.GetNPads(irow); ipad++){
438 Float_t value = irow+ipad/1000.;
439 if (roc1.GetValue(irow,ipad)!=value){
440 printf("Read/Write error\trow=%d\tpad=%d\n",irow,ipad);
445 // write original CalROC to file
446 TFile f("calcTest.root","recreate");
448 AliTPCCalROC * roc2 = (AliTPCCalROC*)f.Get("Roc0");
451 // read CalROC from file and compare to original CalROC
452 for (UInt_t irow = 0; irow <roc0.GetNrows(); irow++){
453 if (roc0.GetNPads(irow)!=roc2->GetNPads(irow))
454 printf("NPads - Read/Write error\trow=%d\n",irow);
455 for (UInt_t ipad = 0; ipad <roc1.GetNPads(irow); ipad++){
456 Float_t value = irow+ipad/1000.;
457 if (roc2->GetValue(irow,ipad)!=value){
458 printf("Read/Write error\trow=%d\tpad=%d\n",irow,ipad);
468 AliTPCCalROC roc3(roc0);
470 for (UInt_t irow = 0; irow <roc3.GetNrows(); irow++){
471 for (UInt_t ipad = 0; ipad <roc3.GetNPads(irow); ipad++){
472 Float_t value = irow+ipad/1000. + 1.5;
473 if (TMath::Abs(roc3.GetValue(irow,ipad)-value) > kEpsilon){
474 printf("Add constant - error\trow=%d\tpad=%d\n",irow,ipad);
479 // Add another CalROC
480 AliTPCCalROC roc4(roc0);
481 roc4.Add(&roc0, -1.5);
482 for (UInt_t irow = 0; irow <roc4.GetNrows(); irow++){
483 for (UInt_t ipad = 0; ipad <roc4.GetNPads(irow); ipad++){
484 Float_t value = irow+ipad/1000. - 1.5 * (irow+ipad/1000.);
485 if (TMath::Abs(roc4.GetValue(irow,ipad)-value) > kEpsilon){
486 printf("Add CalROC - error\trow=%d\tpad=%d\n",irow,ipad);
491 // Multiply with constant
492 AliTPCCalROC roc5(roc0);
494 for (UInt_t irow = 0; irow <roc5.GetNrows(); irow++){
495 for (UInt_t ipad = 0; ipad <roc5.GetNPads(irow); ipad++){
496 Float_t value = (irow+ipad/1000.) * (-1.4);
497 if (TMath::Abs(roc5.GetValue(irow,ipad)-value) > kEpsilon){
498 printf("Multiply with constant - error\trow=%d\tpad=%d\n",irow,ipad);
503 // Multiply another CalROC
504 AliTPCCalROC roc6(roc0);
505 roc6.Multiply(&roc0);
506 for (UInt_t irow = 0; irow <roc6.GetNrows(); irow++){
507 for (UInt_t ipad = 0; ipad <roc6.GetNPads(irow); ipad++){
508 Float_t value = (irow+ipad/1000.) * (irow+ipad/1000.);
509 if (TMath::Abs(roc6.GetValue(irow,ipad)-value) > kEpsilon*100){
510 printf("Multiply with CalROC - error\trow=%d\tpad=%d\n",irow,ipad);
517 AliTPCCalROC roc7(roc0);
519 for (UInt_t irow = 0; irow <roc7.GetNrows(); irow++){
520 for (UInt_t ipad = 0; ipad <roc7.GetNPads(irow); ipad++){
522 if (irow+ipad == 0) value = roc0.GetValue(irow,ipad);
523 if (TMath::Abs(roc7.GetValue(irow,ipad)-value) > kEpsilon){
524 printf("Multiply with CalROC - error\trow=%d\tpad=%d\n",irow,ipad);
533 // create CalROC with defined values
535 AliTPCCalROC sroc0(0);
536 for (UInt_t ichannel = 0; ichannel < sroc0.GetNchannels(); ichannel++){
537 Float_t value = rnd.Gaus(10., 2.);
538 sroc0.SetValue(ichannel,value);
541 printf("Mean (should be close to 10): %f\n", sroc0.GetMean());
542 printf("RMS (should be close to 2): %f\n", sroc0.GetRMS());
543 printf("Median (should be close to 10): %f\n", sroc0.GetMedian());
544 printf("LTM (should be close to 10): %f\n", sroc0.GetLTM());
546 //AliTPCCalROC* sroc1 = sroc0.LocalFit(4, 8);
550 // std::cout << TMath::Abs(roc5.GetValue(irow,ipad)-value) << std::endl;
554 AliTPCCalROC * AliTPCCalROC::LocalFit(Int_t rowRadius, Int_t padRadius, AliTPCCalROC* ROCoutliers, Bool_t robust, Double_t chi2Threshold, Double_t robustFraction) {
556 // MakeLocalFit - smoothing
557 // returns a AliTPCCalROC with smoothed data
558 // rowRadius and padRadius specifies a window around a given pad.
559 // The data of this window are fitted with a parabolic function.
560 // This function is evaluated at the pad's position.
561 // At the edges the window is shifted, so that the specified pad is not anymore in the center of the window.
562 // rowRadius - radius - rows to be used for smoothing
563 // padradius - radius - pads to be used for smoothing
564 // ROCoutlier - map of outliers - pads not to be used for local smoothing
565 // robust - robust method of fitting - (much slower)
566 // chi2Threshold: Threshold for chi2 when EvalRobust is called
567 // robustFraction: Fraction of data that will be used in EvalRobust
569 AliTPCCalROC * ROCfitted = new AliTPCCalROC(fSector);
570 TLinearFitter fitterQ(6,"hyp5");
571 // TLinearFitter fitterQ(6,"x0++x1++x2++x3++x4++x5");
572 fitterQ.StoreData(kTRUE);
573 for (UInt_t row=0; row < GetNrows(); row++) {
574 //std::cout << "Entering row " << row << " of " << GetNrows() << " @ sector "<< fSector << " for local fitting... "<< std::endl;
575 for (UInt_t pad=0; pad < GetNPads(row); pad++)
576 ROCfitted->SetValue(row, pad, GetNeighbourhoodValue(&fitterQ, row, pad, rowRadius, padRadius, ROCoutliers, robust, chi2Threshold, robustFraction));
582 Double_t AliTPCCalROC::GetNeighbourhoodValue(TLinearFitter* fitterQ, Int_t row, Int_t pad, Int_t rRadius, Int_t pRadius, AliTPCCalROC* ROCoutliers, Bool_t robust, Double_t chi2Threshold, Double_t robustFraction) {
584 // AliTPCCalROC::GetNeighbourhoodValue - smoothing - PRIVATE
585 // in this function the fit for LocalFit is done
588 fitterQ->ClearPoints();
589 TVectorD fitParam(6);
593 Float_t lPad[3] = {0};
594 Float_t localXY[3] = {0};
596 AliTPCROC* tpcROCinstance = AliTPCROC::Instance();
597 tpcROCinstance->GetPositionLocal(fSector, row, pad, lPad); // calculate position lPad by pad and row number
599 TArrayI *neighbourhoodRows = 0;
600 TArrayI *neighbourhoodPads = 0;
602 //std::cerr << "Trying to get neighbourhood for row " << row << ", pad " << pad << std::endl;
603 GetNeighbourhood(neighbourhoodRows, neighbourhoodPads, row, pad, rRadius, pRadius);
604 //std::cerr << "Got neighbourhood for row " << row << ", pad " << pad << std::endl;
607 for (Int_t i=0; i < (2*rRadius+1)*(2*pRadius+1); i++) {
608 r = neighbourhoodRows->At(i);
609 p = neighbourhoodPads->At(i);
610 if (r == -1 || p == -1) continue; // window is bigger than ROC
611 tpcROCinstance->GetPositionLocal(fSector, r, p, localXY); // calculate position localXY by pad and row number
612 dlx = lPad[0] - localXY[0];
613 dly = lPad[1] - localXY[1];
620 if (!ROCoutliers || ROCoutliers->GetValue(r,p) != 1) {
621 fitterQ->AddPoint(&xx[1], GetValue(r, p), 1);
626 delete neighbourhoodRows;
627 delete neighbourhoodPads;
629 if (npoints < 0.5 * ((2*rRadius+1)*(2*pRadius+1)) ) {
630 // std::cerr << "Too few data points for fitting @ row " << row << ", pad " << pad << " in sector " << fSector << std::endl;
631 return 0.; // for diagnostic
634 fitterQ->GetParameters(fitParam);
636 if (robust) chi2Q = fitterQ->GetChisquare()/(npoints-6.);
637 //if (robust) chi2Q = fitterQ->GetChisquare()/(npoints-6.);
638 if (robust && chi2Q > chi2Threshold) {
639 //std::cout << "robust fitter called... " << std::endl;
640 fitterQ->EvalRobust(robustFraction);
641 fitterQ->GetParameters(fitParam);
643 Double_t value = fitParam[0];
645 //if (value < 0) std::cerr << "negative fit-value " << value << " in sector "<< this->fSector << " @ row: " << row << " and pad: " << pad << ", with fitter Chi2 = " << chi2Q << std::endl;
652 void AliTPCCalROC::GetNeighbourhood(TArrayI* &rowArray, TArrayI* &padArray, Int_t row, Int_t pad, Int_t rRadius, Int_t pRadius) {
654 // AliTPCCalROC::GetNeighbourhood - PRIVATE
655 // in this function the window for LocalFit is determined
657 rowArray = new TArrayI((2*rRadius+1)*(2*pRadius+1));
658 padArray = new TArrayI((2*rRadius+1)*(2*pRadius+1));
660 Int_t rmin = row - rRadius;
661 UInt_t rmax = row + rRadius;
663 // if window goes out of ROC
668 if (rmax >= GetNrows()) {
669 rmin = rmin - (rmax - GetNrows()+1);
670 rmax = GetNrows() - 1;
671 if (rmin < 0 ) rmin = 0; // if the window is bigger than the ROC
677 for (UInt_t r = rmin; r <= rmax; r++) {
678 pmin = pad - pRadius;
679 pmax = pad + pRadius;
684 if (pmax >= GetNPads(r)) {
685 pmin = pmin - (pmax - GetNPads(r)+1);
686 pmax = GetNPads(r) - 1;
687 if (pmin < 0 ) pmin = 0; // if the window is bigger than the ROC
689 for (Int_t p = pmin; p <= pmax; p++) {
695 for (Int_t j = i; j < rowArray->GetSize(); j++){ // unused padArray-entries, in the case that the window is bigger than the ROC
696 //std::cout << "trying to write -1" << std::endl;
699 //std::cout << "writing -1" << std::endl;
705 void AliTPCCalROC::GlobalFit(const AliTPCCalROC* ROCoutliers, Bool_t robust, TVectorD &fitParam, TMatrixD &covMatrix, Float_t & chi2, Int_t fitType, Double_t chi2Threshold, Double_t robustFraction){
707 // Makes a GlobalFit for the given secotr and return fit-parameters, covariance and chi2
708 // The origin of the fit function is the center of the ROC!
709 // fitType == 0: fit plane function
710 // fitType == 1: fit parabolic function
711 // ROCoutliers - pads with value !=0 are not used in fitting procedure
712 // chi2Threshold: Threshold for chi2 when EvalRobust is called
713 // robustFraction: Fraction of data that will be used in EvalRobust
715 TLinearFitter* fitterG = 0;
719 fitterG = new TLinearFitter (6,"x0++x1++x2++x3++x4++x5");
721 fitterG = new TLinearFitter(3,"x0++x1++x2");
722 fitterG->StoreData(kTRUE);
723 fitterG->ClearPoints();
727 Float_t centerPad[3] = {0};
728 Float_t localXY[3] = {0};
730 AliTPCROC* tpcROCinstance = AliTPCROC::Instance();
731 tpcROCinstance->GetPositionLocal(fSector, GetNrows()/2, GetNPads(GetNrows()/2)/2, centerPad); // calculate center of ROC
733 // loop over all channels and read data into fitterG
734 if (fitType == 1) { // parabolic fit
735 fitParam.ResizeTo(6);
736 covMatrix.ResizeTo(6,6);
737 for (UInt_t irow = 0; irow < GetNrows(); irow++) {
738 for (UInt_t ipad = 0; ipad < GetNPads(irow); ipad++) {
740 tpcROCinstance->GetPositionLocal(fSector, irow, ipad, localXY); // calculate position localXY by pad and row number
741 dlx = centerPad[0] - localXY[0];
742 dly = centerPad[1] - localXY[1];
749 if (!ROCoutliers || ROCoutliers->GetValue(irow, ipad) != 1) {
751 fitterG->AddPoint(xx, GetValue(irow, ipad), 1);
757 fitParam.ResizeTo(3);
758 covMatrix.ResizeTo(3,3);
759 for (UInt_t irow = 0; irow < GetNrows(); irow++) {
760 for (UInt_t ipad = 0; ipad < GetNPads(irow); ipad++) {
762 tpcROCinstance->GetPositionLocal(fSector, irow, ipad, localXY); // calculate position localXY by pad and row number
763 dlx = centerPad[0] - localXY[0];
764 dly = centerPad[1] - localXY[1];
768 if (!ROCoutliers || ROCoutliers->GetValue(irow, ipad) != 1) {
770 fitterG->AddPoint(xx, GetValue(irow, ipad), 1);
776 fitterG->GetParameters(fitParam);
777 fitterG->GetCovarianceMatrix(covMatrix);
779 chi2 = fitterG->GetChisquare()/(npoints-6.);
780 else chi2 = fitterG->GetChisquare()/(npoints-3.);
781 if (robust && chi2 > chi2Threshold) {
782 // std::cout << "robust fitter called... " << std::endl;
783 fitterG->EvalRobust(robustFraction);
784 fitterG->GetParameters(fitParam);
790 AliTPCCalROC* AliTPCCalROC::CreateGlobalFitCalROC(TVectorD &fitParam, Int_t sector){
792 // Create ROC with global fit parameters
793 // The origin of the fit function is the center of the ROC!
794 // loop over all channels, write fit values into new ROC and return it
797 Float_t centerPad[3] = {0};
798 Float_t localXY[3] = {0};
799 AliTPCCalROC * ROCfitted = new AliTPCCalROC(sector);
800 AliTPCROC* tpcROCinstance = AliTPCROC::Instance();
801 tpcROCinstance->GetPositionLocal(sector, ROCfitted->GetNrows()/2, ROCfitted->GetNPads(ROCfitted->GetNrows()/2)/2, centerPad); // calculate center of ROC
803 if (fitParam.GetNoElements() == 6) fitType = 1;
806 if (fitType == 1) { // parabolic fit
807 for (UInt_t irow = 0; irow < ROCfitted->GetNrows(); irow++) {
808 for (UInt_t ipad = 0; ipad < ROCfitted->GetNPads(irow); ipad++) {
809 tpcROCinstance->GetPositionLocal(sector, irow, ipad, localXY); // calculate position localXY by pad and row number
810 dlx = centerPad[0] - localXY[0];
811 dly = centerPad[1] - localXY[1];
812 value = fitParam[0] + fitParam[1]*dlx + fitParam[2]*dly + fitParam[3]*dlx*dlx + fitParam[4]*dly*dly + fitParam[5]*dlx*dly;
813 ROCfitted->SetValue(irow, ipad, value);
818 for (UInt_t irow = 0; irow < ROCfitted->GetNrows(); irow++) {
819 for (UInt_t ipad = 0; ipad < ROCfitted->GetNPads(irow); ipad++) {
820 tpcROCinstance->GetPositionLocal(sector, irow, ipad, localXY); // calculate position localXY by pad and row number
821 dlx = centerPad[0] - localXY[0];
822 dly = centerPad[1] - localXY[1];
823 value = fitParam[0] + fitParam[1]*dlx + fitParam[2]*dly;
824 ROCfitted->SetValue(irow, ipad, value);