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"
36 ClassImp(AliTPCCalPad)
38 //_____________________________________________________________________________
39 AliTPCCalPad::AliTPCCalPad():TNamed()
42 // AliTPCCalPad default constructor
45 for (Int_t isec = 0; isec < kNsec; isec++) {
51 //_____________________________________________________________________________
52 AliTPCCalPad::AliTPCCalPad(const Text_t *name, const Text_t *title)
56 // AliTPCCalPad constructor
58 for (Int_t isec = 0; isec < kNsec; isec++) {
59 fROC[isec] = new AliTPCCalROC(isec);
64 //_____________________________________________________________________________
65 AliTPCCalPad::AliTPCCalPad(const AliTPCCalPad &c):TNamed(c)
68 // AliTPCCalPad copy constructor
71 for (Int_t isec = 0; isec < kNsec; isec++) {
74 fROC[isec] = new AliTPCCalROC(*(c.fROC[isec]));
78 //_____________________________________________________________________________
79 AliTPCCalPad::AliTPCCalPad(TObjArray * array):TNamed()
82 // AliTPCCalPad default constructor
85 for (Int_t isec = 0; isec < kNsec; isec++) {
86 fROC[isec] = (AliTPCCalROC *)array->At(isec);
92 ///_____________________________________________________________________________
93 AliTPCCalPad::~AliTPCCalPad()
96 // AliTPCCalPad destructor
99 for (Int_t isec = 0; isec < kNsec; isec++) {
108 //_____________________________________________________________________________
109 AliTPCCalPad &AliTPCCalPad::operator=(const AliTPCCalPad &c)
112 // Assignment operator
115 if (this != &c) ((AliTPCCalPad &) c).Copy(*this);
120 //_____________________________________________________________________________
121 void AliTPCCalPad::Copy(TObject &c) const
127 for (Int_t isec = 0; isec < kNsec; isec++) {
129 fROC[isec]->Copy(*((AliTPCCalPad &) c).fROC[isec]);
136 void AliTPCCalPad::SetCalROC(AliTPCCalROC* roc, Int_t sector){
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
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));
149 //_____________________________________________________________________________
150 void AliTPCCalPad::Add(Float_t c1)
153 // add constant c1 to all channels of all ROCs
156 for (Int_t isec = 0; isec < kNsec; isec++) {
163 //_____________________________________________________________________________
164 void AliTPCCalPad::Multiply(Float_t c1)
167 // multiply each channel of all ROCs with c1
169 for (Int_t isec = 0; isec < kNsec; isec++) {
171 fROC[isec]->Multiply(c1);
176 //_____________________________________________________________________________
177 void AliTPCCalPad::Add(const AliTPCCalPad * pad, Double_t c1)
180 // multiply AliTPCCalPad 'pad' by c1 and add each channel to the coresponing channel in all ROCs
183 for (Int_t isec = 0; isec < kNsec; isec++) {
185 fROC[isec]->Add(pad->GetCalROC(isec),c1);
190 //_____________________________________________________________________________
191 void AliTPCCalPad::Multiply(const AliTPCCalPad * pad)
194 // multiply each channel of all ROCs with the coresponding channel of 'pad'
197 for (Int_t isec = 0; isec < kNsec; isec++) {
199 fROC[isec]->Multiply(pad->GetCalROC(isec));
204 //_____________________________________________________________________________
205 void AliTPCCalPad::Divide(const AliTPCCalPad * pad)
208 // divide each channel of all ROCs by the coresponding channel of 'pad'
211 for (Int_t isec = 0; isec < kNsec; isec++) {
213 fROC[isec]->Divide(pad->GetCalROC(isec));
218 //_____________________________________________________________________________
219 TGraph * AliTPCCalPad::MakeGraph(Int_t type, Float_t ratio){
226 for (Int_t i=0;i<72;i++) if (fROC[i]) npoints++;
227 TGraph * graph = new TGraph(npoints);
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));
237 graph->GetXaxis()->SetTitle("Sector");
239 graph->GetYaxis()->SetTitle("Mean");
240 graph->SetMarkerStyle(22);
243 graph->GetYaxis()->SetTitle("Median");
244 graph->SetMarkerStyle(22);
247 graph->GetYaxis()->SetTitle(Form("Mean%f",ratio));
248 graph->SetMarkerStyle(24);
254 //_____________________________________________________________________________
255 Double_t AliTPCCalPad::GetMeanRMS(Double_t &rms)
258 // Calculates mean and RMS of all ROCs
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];
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);
276 Double_t mean = sum*n1;
277 rms = TMath::Sqrt(TMath::Abs(sum2*n1-mean*mean));
282 //_____________________________________________________________________________
283 Double_t AliTPCCalPad::GetMean(AliTPCCalPad* outlierPad)
286 // return mean of the mean of all ROCs
290 for (Int_t isec = 0; isec < kNsec; isec++) {
291 AliTPCCalROC *calRoc = fROC[isec];
293 AliTPCCalROC* outlierROC = 0;
294 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
295 arr[n] = calRoc->GetMean(outlierROC);
299 return TMath::Mean(n,arr);
302 //_____________________________________________________________________________
303 Double_t AliTPCCalPad::GetRMS(AliTPCCalPad* outlierPad)
306 // return mean of the RMS of all ROCs
310 for (Int_t isec = 0; isec < kNsec; isec++) {
311 AliTPCCalROC *calRoc = fROC[isec];
313 AliTPCCalROC* outlierROC = 0;
314 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
315 arr[n] = calRoc->GetRMS(outlierROC);
319 return TMath::Mean(n,arr);
322 //_____________________________________________________________________________
323 Double_t AliTPCCalPad::GetMedian(AliTPCCalPad* outlierPad)
326 // return mean of the median of all ROCs
330 for (Int_t isec = 0; isec < kNsec; isec++) {
331 AliTPCCalROC *calRoc = fROC[isec];
333 AliTPCCalROC* outlierROC = 0;
334 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
335 arr[n] = calRoc->GetMedian(outlierROC);
339 return TMath::Mean(n,arr);
342 //_____________________________________________________________________________
343 Double_t AliTPCCalPad::GetLTM(Double_t *sigma, Double_t fraction, AliTPCCalPad* outlierPad)
346 // return mean of the LTM and sigma of all ROCs
348 Double_t arrm[kNsec];
349 Double_t arrs[kNsec];
353 for (Int_t isec = 0; isec < kNsec; isec++) {
354 AliTPCCalROC *calRoc = fROC[isec];
356 if ( sigma ) sTemp=arrs+n;
357 AliTPCCalROC* outlierROC = 0;
358 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
359 arrm[n] = calRoc->GetLTM(sTemp,fraction, outlierROC);
363 if ( sigma ) *sigma = TMath::Mean(n,arrs);
364 return TMath::Mean(n,arrm);
367 //_____________________________________________________________________________
368 TH1F * AliTPCCalPad::MakeHisto1D(Float_t min, Float_t max,Int_t type){
371 // type -1 = user defined range
372 // 0 = nsigma cut nsigma=min
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;
385 Float_t mean = GetMedian();
392 // LTM mean +- nsigma
395 Float_t mean = GetLTM(&sigma,max);
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++) {
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));
417 //_____________________________________________________________________________
418 TH2F *AliTPCCalPad::MakeHisto2D(Int_t side){
421 // side - specify the side A = 0 C = 1
422 // type - used types of determination of boundaries in z
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;
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){
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);
444 his->SetXTitle("x (cm)");
445 his->SetYTitle("y (cm)");
450 AliTPCCalPad* AliTPCCalPad::LocalFit(const char* padName, Int_t rowRadius, Int_t padRadius, AliTPCCalPad* PadOutliers, Bool_t robust, Double_t chi2Threshold, Double_t robustFraction){
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
466 AliTPCCalPad* pad = new AliTPCCalPad(padName, padName);
467 for (Int_t isec = 0; isec < 72; isec++){
469 SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, PadOutliers->GetCalROC(isec), robust, chi2Threshold, robustFraction));
471 SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, 0, robust, chi2Threshold, robustFraction));
477 AliTPCCalPad* AliTPCCalPad::GlobalFit(const char* padName, AliTPCCalPad* PadOutliers, Bool_t robust, Int_t fitType, Double_t chi2Threshold, Double_t robustFraction){
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
486 AliTPCCalPad* pad = new AliTPCCalPad(padName, padName);
487 TVectorD fitParam(0);
488 TMatrixD covMatrix(0,0);
490 for (Int_t isec = 0; isec < 72; isec++){
492 GetCalROC(isec)->GlobalFit(PadOutliers->GetCalROC(isec), robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction);
494 GetCalROC(isec)->GlobalFit(0, robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction);
495 pad->SetCalROC(AliTPCCalROC::CreateGlobalFitCalROC(fitParam, isec));
501 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){
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
510 TLinearFitter* fitterGA = 0;
511 TLinearFitter* fitterGC = 0;
514 fitterGA = new TLinearFitter (6,"x0++x1++x2++x3++x4++x5");
515 fitterGC = new TLinearFitter (6,"x0++x1++x2++x3++x4++x5");
518 fitterGA = new TLinearFitter(3,"x0++x1++x2");
519 fitterGC = new TLinearFitter(3,"x0++x1++x2");
521 fitterGA->StoreData(kTRUE);
522 fitterGC->StoreData(kTRUE);
523 fitterGA->ClearPoints();
524 fitterGC->ClearPoints();
529 Float_t localXY[3] = {0}; // pad's position, needed to get the pad's position
530 Float_t lx, ly; // pads position
532 AliTPCROC* tpcROCinstance = AliTPCROC::Instance(); // to calculate the pad's position
534 // loop over all sectors and pads and read data into fitterGA and fitterGC
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++) {
545 tpcROCinstance->GetPositionLocal(isec, irow, ipad, localXY); // calculate position localXY by sector, pad and row number
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
562 fitterGA->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);
566 fitterGC->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);
575 fitParamSideA.ResizeTo(3);
576 fitParamSideC.ResizeTo(3);
577 covMatrixSideA.ResizeTo(3,3);
578 covMatrixSideC.ResizeTo(3,3);
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++) {
584 tpcROCinstance->GetPositionLocal(isec, irow, ipad, localXY); // calculate position localXY by sector, pad and row number
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)) {
599 fitterGA->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);
604 fitterGC->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);
614 fitterGA->GetParameters(fitParamSideA);
615 fitterGC->GetParameters(fitParamSideC);
616 fitterGA->GetCovarianceMatrix(covMatrixSideA);
617 fitterGC->GetCovarianceMatrix(covMatrixSideC);
619 chi2SideA = fitterGA->GetChisquare()/(npointsA-6.);
620 chi2SideC = fitterGC->GetChisquare()/(npointsC-6.);
623 chi2SideA = fitterGA->GetChisquare()/(npointsA-3.);
624 chi2SideC = fitterGC->GetChisquare()/(npointsC-3.);
626 if (robust && chi2SideA > chi2Threshold) {
627 // std::cout << "robust fitter called... " << std::endl;
628 fitterGA->EvalRobust(robustFraction);
629 fitterGA->GetParameters(fitParamSideA);
631 if (robust && chi2SideC > chi2Threshold) {
632 // std::cout << "robust fitter called... " << std::endl;
633 fitterGC->EvalRobust(robustFraction);
634 fitterGC->GetParameters(fitParamSideC);
641 void AliTPCCalPad::MakeTree(const char * fileName, TObjArray * array, const char * mapFileName, AliTPCCalPad* outlierPad, Float_t ltmFraction) {
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 -
649 AliTPCROC* tpcROCinstance = AliTPCROC::Instance();
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;
659 TFile mapFile(mapFileName, "read");
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];
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);
677 for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
678 mapIROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(0));
679 mapOROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(36));
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);
687 } // if (mapFileName)
689 TTreeSRedirector cstream(fileName);
690 Int_t arrayEntries = array->GetEntries();
692 TString* names = new TString[arrayEntries];
693 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
694 names[ivalue].Append(((AliTPCCalPad*)array->At(ivalue))->GetName());
696 for (UInt_t isector = 0; isector < tpcROCinstance->GetNSectors(); isector++) {
698 // get statistic for given sector
700 TVectorF median(arrayEntries);
701 TVectorF mean(arrayEntries);
702 TVectorF rms(arrayEntries);
703 TVectorF ltm(arrayEntries);
704 TVectorF ltmrms(arrayEntries);
705 TVectorF medianWithOut(arrayEntries);
706 TVectorF meanWithOut(arrayEntries);
707 TVectorF rmsWithOut(arrayEntries);
708 TVectorF ltmWithOut(arrayEntries);
709 TVectorF ltmrmsWithOut(arrayEntries);
711 TVectorF *vectorArray = new TVectorF[arrayEntries];
712 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
713 vectorArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
715 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
716 AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
717 AliTPCCalROC* calROC = calPad->GetCalROC(isector);
718 AliTPCCalROC* outlierROC = 0;
719 if (outlierPad) outlierROC = outlierPad->GetCalROC(isector);
721 median[ivalue] = calROC->GetMedian();
722 mean[ivalue] = calROC->GetMean();
723 rms[ivalue] = calROC->GetRMS();
724 Double_t ltmrmsValue = 0;
725 ltm[ivalue] = calROC->GetLTM(<mrmsValue, ltmFraction);
726 ltmrms[ivalue] = ltmrmsValue;
728 medianWithOut[ivalue] = calROC->GetMedian(outlierROC);
729 meanWithOut[ivalue] = calROC->GetMean(outlierROC);
730 rmsWithOut[ivalue] = calROC->GetRMS(outlierROC);
732 ltmWithOut[ivalue] = calROC->GetLTM(<mrmsValue, ltmFraction, outlierROC);
733 ltmrmsWithOut[ivalue] = ltmrmsValue;
742 medianWithOut[ivalue] = 0.;
743 meanWithOut[ivalue] = 0.;
744 rmsWithOut[ivalue] = 0.;
745 ltmWithOut[ivalue] = 0.;
746 ltmrmsWithOut[ivalue] = 0.;
751 // fill vectors of variable per pad
753 TVectorF *posArray = new TVectorF[8];
754 for (Int_t ivalue = 0; ivalue < 8; ivalue++)
755 posArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
757 Float_t posG[3] = {0};
758 Float_t posL[3] = {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];
770 posArray[6][ichannel] = (Int_t)(ipad - (Double_t)(tpcROCinstance->GetNPads(isector, irow))/2);
771 posArray[7][ichannel] = ichannel;
773 // loop over array containing AliTPCCalPads
774 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
775 AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
776 AliTPCCalROC* calROC = calPad->GetCalROC(isector);
778 (vectorArray[ivalue])[ichannel] = calROC->GetValue(irow, ipad);
780 (vectorArray[ivalue])[ichannel] = 0;
786 cstream << "calPads" <<
787 "sector=" << isector;
789 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
790 cstream << "calPads" <<
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] <<
794 (Char_t*)((names[ivalue] + "_LTM=").Data()) << ltm[ivalue] <<
795 (Char_t*)((names[ivalue] + "_RMS_LTM=").Data()) << ltmrms[ivalue];
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];
806 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
807 cstream << "calPads" <<
808 (Char_t*)((names[ivalue] + ".=").Data()) << &vectorArray[ivalue];
812 for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
814 cstream << "calPads" <<
815 (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapIROCArray[ivalue];
817 cstream << "calPads" <<
818 (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapOROCArray[ivalue];
822 cstream << "calPads" <<
823 "row.=" << &posArray[0] <<
824 "pad.=" << &posArray[1] <<
825 "lx.=" << &posArray[2] <<
826 "ly.=" << &posArray[3] <<
827 "gx.=" << &posArray[4] <<
828 "gy.=" << &posArray[5] <<
829 "rpad.=" << &posArray[6] <<
830 "channel.=" << &posArray[7];
832 cstream << "calPads" <<
836 delete[] vectorArray;
842 delete[] mapIROCArray;
843 delete[] mapOROCArray;