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>
47 #include <TVirtualPad.h>
50 ClassImp(AliTPCCalPad)
52 //_____________________________________________________________________________
53 AliTPCCalPad::AliTPCCalPad():TNamed()
56 // AliTPCCalPad default constructor
59 for (Int_t isec = 0; isec < kNsec; isec++) {
65 //_____________________________________________________________________________
66 AliTPCCalPad::AliTPCCalPad(const Text_t *name, const Text_t *title)
70 // AliTPCCalPad constructor
72 for (Int_t isec = 0; isec < kNsec; isec++) {
73 fROC[isec] = new AliTPCCalROC(isec);
78 //_____________________________________________________________________________
79 AliTPCCalPad::AliTPCCalPad(const AliTPCCalPad &c):TNamed(c)
82 // AliTPCCalPad copy constructor
85 for (Int_t isec = 0; isec < kNsec; isec++) {
88 fROC[isec] = new AliTPCCalROC(*(c.fROC[isec]));
92 //_____________________________________________________________________________
93 AliTPCCalPad::AliTPCCalPad(TObjArray * array):TNamed()
96 // AliTPCCalPad default constructor
99 for (Int_t isec = 0; isec < kNsec; isec++) {
100 fROC[isec] = (AliTPCCalROC *)array->At(isec);
106 ///_____________________________________________________________________________
107 AliTPCCalPad::~AliTPCCalPad()
110 // AliTPCCalPad destructor
113 for (Int_t isec = 0; isec < kNsec; isec++) {
122 //_____________________________________________________________________________
123 AliTPCCalPad &AliTPCCalPad::operator=(const AliTPCCalPad &c)
126 // Assignment operator
129 if (this != &c) ((AliTPCCalPad &) c).Copy(*this);
134 //_____________________________________________________________________________
135 void AliTPCCalPad::Copy(TObject &c) const
141 for (Int_t isec = 0; isec < kNsec; isec++) {
143 fROC[isec]->Copy(*((AliTPCCalPad &) c).fROC[isec]);
150 void AliTPCCalPad::SetCalROC(AliTPCCalROC* roc, Int_t sector){
152 // Set AliTPCCalROC copies values from 'roc'
153 // if sector == -1 the sector specified in 'roc' is used
154 // else sector specified in 'roc' is ignored and specified sector is filled
156 if (sector == -1) sector = roc->GetSector();
157 if (!fROC[sector]) fROC[sector] = new AliTPCCalROC(sector);
158 for (UInt_t ichannel = 0; ichannel < roc->GetNchannels(); ichannel++)
159 fROC[sector]->SetValue(ichannel, roc->GetValue(ichannel));
164 //_____________________________________________________________________________
165 void AliTPCCalPad::Add(Float_t c1)
168 // add constant c1 to all channels of all ROCs
171 for (Int_t isec = 0; isec < kNsec; isec++) {
178 //_____________________________________________________________________________
179 void AliTPCCalPad::Multiply(Float_t c1)
182 // multiply each channel of all ROCs with c1
184 for (Int_t isec = 0; isec < kNsec; isec++) {
186 fROC[isec]->Multiply(c1);
191 //_____________________________________________________________________________
192 void AliTPCCalPad::Add(const AliTPCCalPad * pad, Double_t c1)
195 // multiply AliTPCCalPad 'pad' by c1 and add each channel to the coresponing channel in all ROCs
198 for (Int_t isec = 0; isec < kNsec; isec++) {
199 if (fROC[isec] && pad->GetCalROC(isec)){
200 fROC[isec]->Add(pad->GetCalROC(isec),c1);
205 //_____________________________________________________________________________
206 void AliTPCCalPad::Multiply(const AliTPCCalPad * pad)
209 // multiply each channel of all ROCs with the coresponding channel of 'pad'
212 for (Int_t isec = 0; isec < kNsec; isec++) {
214 fROC[isec]->Multiply(pad->GetCalROC(isec));
219 //_____________________________________________________________________________
220 void AliTPCCalPad::Divide(const AliTPCCalPad * pad)
223 // divide each channel of all ROCs by the coresponding channel of 'pad'
226 for (Int_t isec = 0; isec < kNsec; isec++) {
228 fROC[isec]->Divide(pad->GetCalROC(isec));
233 //_____________________________________________________________________________
234 TGraph * AliTPCCalPad::MakeGraph(Int_t type, Float_t ratio){
241 for (Int_t i=0;i<72;i++) if (fROC[i]) npoints++;
242 TGraph * graph = new TGraph(npoints);
244 for (Int_t isec=0;isec<72;isec++){
245 if (!fROC[isec]) continue;
246 if (type==0) graph->SetPoint(npoints,isec,fROC[isec]->GetMean());
247 if (type==1) graph->SetPoint(npoints,isec,fROC[isec]->GetMedian());
248 if (type==2) graph->SetPoint(npoints,isec,fROC[isec]->GetLTM(0,ratio));
252 graph->GetXaxis()->SetTitle("Sector");
254 graph->GetYaxis()->SetTitle("Mean");
255 graph->SetMarkerStyle(22);
258 graph->GetYaxis()->SetTitle("Median");
259 graph->SetMarkerStyle(22);
262 graph->GetYaxis()->SetTitle(Form("Mean%f",ratio));
263 graph->SetMarkerStyle(24);
269 //_____________________________________________________________________________
270 Double_t AliTPCCalPad::GetMeanRMS(Double_t &rms)
273 // Calculates mean and RMS of all ROCs
275 Double_t sum = 0, sum2 = 0, n=0, val=0;
276 for (Int_t isec = 0; isec < kNsec; isec++) {
277 AliTPCCalROC *calRoc = fROC[isec];
279 for (UInt_t irow=0; irow<calRoc->GetNrows(); irow++){
280 for (UInt_t ipad=0; ipad<calRoc->GetNPads(irow); ipad++){
281 val = calRoc->GetValue(irow,ipad);
291 Double_t mean = sum*n1;
292 rms = TMath::Sqrt(TMath::Abs(sum2*n1-mean*mean));
297 //_____________________________________________________________________________
298 Double_t AliTPCCalPad::GetMean(AliTPCCalPad* outlierPad)
301 // return mean of the mean of all ROCs
305 for (Int_t isec = 0; isec < kNsec; isec++) {
306 AliTPCCalROC *calRoc = fROC[isec];
308 AliTPCCalROC* outlierROC = 0;
309 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
310 arr[n] = calRoc->GetMean(outlierROC);
314 return TMath::Mean(n,arr);
317 //_____________________________________________________________________________
318 Double_t AliTPCCalPad::GetRMS(AliTPCCalPad* outlierPad)
321 // return mean of the RMS of all ROCs
325 for (Int_t isec = 0; isec < kNsec; isec++) {
326 AliTPCCalROC *calRoc = fROC[isec];
328 AliTPCCalROC* outlierROC = 0;
329 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
330 arr[n] = calRoc->GetRMS(outlierROC);
334 return TMath::Mean(n,arr);
337 //_____________________________________________________________________________
338 Double_t AliTPCCalPad::GetMedian(AliTPCCalPad* outlierPad)
341 // return mean of the median of all ROCs
345 for (Int_t isec = 0; isec < kNsec; isec++) {
346 AliTPCCalROC *calRoc = fROC[isec];
348 AliTPCCalROC* outlierROC = 0;
349 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
350 arr[n] = calRoc->GetMedian(outlierROC);
354 return TMath::Mean(n,arr);
357 //_____________________________________________________________________________
358 Double_t AliTPCCalPad::GetLTM(Double_t *sigma, Double_t fraction, AliTPCCalPad* outlierPad)
361 // return mean of the LTM and sigma of all ROCs
363 Double_t arrm[kNsec];
364 Double_t arrs[kNsec];
368 for (Int_t isec = 0; isec < kNsec; isec++) {
369 AliTPCCalROC *calRoc = fROC[isec];
371 if ( sigma ) sTemp=arrs+n;
372 AliTPCCalROC* outlierROC = 0;
373 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
374 arrm[n] = calRoc->GetLTM(sTemp,fraction, outlierROC);
378 if ( sigma ) *sigma = TMath::Mean(n,arrs);
379 return TMath::Mean(n,arrm);
382 //_____________________________________________________________________________
383 TH1F * AliTPCCalPad::MakeHisto1D(Float_t min, Float_t max,Int_t type, Int_t side){
386 // type -1 = user defined range
387 // 0 = nsigma cut nsigma=min
392 Float_t mean = GetMean();
393 Float_t sigma = GetRMS();
394 Float_t nsigma = TMath::Abs(min);
395 min = mean-nsigma*sigma;
396 max = mean+nsigma*sigma;
400 Float_t mean = GetMedian();
407 // LTM mean +- nsigma
410 Float_t mean = GetLTM(&sigma,max);
416 TString name=Form("%s Pad 1D",GetTitle());
417 TH1F * his = new TH1F(name.Data(),name.Data(),100, min,max);
418 for (Int_t isec = 0; isec < kNsec; isec++) {
419 if (side==1 && isec%36>18) continue;
420 if (side==-1 && isec%36<18) continue;
422 for (UInt_t irow=0; irow<fROC[isec]->GetNrows(); irow++){
423 UInt_t npads = (Int_t)fROC[isec]->GetNPads(irow);
424 for (UInt_t ipad=0; ipad<npads; ipad++){
425 his->Fill(fROC[isec]->GetValue(irow,ipad));
433 //_____________________________________________________________________________
434 TH2F *AliTPCCalPad::MakeHisto2D(Int_t side){
437 // side - specify the side A = 0 C = 1
438 // type - used types of determination of boundaries in z
440 Float_t kEpsilon = 0.000000000001;
441 TH2F * his = new TH2F(GetName(), GetName(), 250,-250,250,250,-250,250);
442 AliTPCROC * roc = AliTPCROC::Instance();
443 for (Int_t isec=0; isec<72; isec++){
444 if (side==0 && isec%36>=18) continue;
445 if (side>0 && isec%36<18) continue;
447 AliTPCCalROC * calRoc = fROC[isec];
448 for (UInt_t irow=0; irow<calRoc->GetNrows(); irow++)
449 for (UInt_t ipad=0; ipad<calRoc->GetNPads(irow); ipad++)
450 if (TMath::Abs(calRoc->GetValue(irow,ipad))>kEpsilon){
452 roc->GetPositionGlobal(isec,irow,ipad,xyz);
453 Int_t binx = 1+TMath::Nint((xyz[0]+250.)*0.5);
454 Int_t biny = 1+TMath::Nint((xyz[1]+250.)*0.5);
455 Float_t value = calRoc->GetValue(irow,ipad);
456 his->SetBinContent(binx,biny,value);
460 his->SetXTitle("x (cm)");
461 his->SetYTitle("y (cm)");
466 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 {
468 // Loops over all AliTPCCalROCs and performs a localFit in each ROC
469 // AliTPCCalPad with fit-data is returned
470 // rowRadius and padRadius specifies a window around a given pad in one sector.
471 // The data of this window are fitted with a parabolic function.
472 // This function is evaluated at the pad's position.
473 // At the edges the window is shifted, so that the specified pad is not anymore in the center of the window.
474 // rowRadius - radius - rows to be used for smoothing
475 // padradius - radius - pads to be used for smoothing
476 // ROCoutlier - map of outliers - pads not to be used for local smoothing
477 // robust - robust method of fitting - (much slower)
478 // chi2Threshold: Threshold for chi2 when EvalRobust is called
479 // robustFraction: Fraction of data that will be used in EvalRobust
482 AliTPCCalPad* pad = new AliTPCCalPad(padName, padName);
483 for (Int_t isec = 0; isec < 72; isec++){
484 if (printCurrentSector) std::cout << "LocalFit in sector " << isec << "\r" << std::flush;
486 pad->SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, PadOutliers->GetCalROC(isec), robust, chi2Threshold, robustFraction));
488 pad->SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, 0, robust, chi2Threshold, robustFraction));
494 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){
496 // Loops over all AliTPCCalROCs and performs a globalFit in each ROC
497 // AliTPCCalPad with fit-data is returned
498 // chi2Threshold: Threshold for chi2 when EvalRobust is called
499 // robustFraction: Fraction of data that will be used in EvalRobust
500 // chi2Threshold: Threshold for chi2 when EvalRobust is called
501 // robustFraction: Fraction of data that will be used in EvalRobust
502 // err: error of the data points
503 // if fitParArr and/or fitCovArr is given, write fitParameters and/or covariance Matrices into the array
505 AliTPCCalPad* pad = new AliTPCCalPad(padName, padName);
506 TVectorD fitParam(0);
507 TMatrixD covMatrix(0,0);
509 for (Int_t isec = 0; isec < 72; isec++){
511 GetCalROC(isec)->GlobalFit(PadOutliers->GetCalROC(isec), robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction, err);
513 GetCalROC(isec)->GlobalFit(0, robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction, err);
515 AliTPCCalROC *roc=AliTPCCalROC::CreateGlobalFitCalROC(fitParam, isec);
518 if ( fitParArr ) fitParArr->AddAtAndExpand(new TVectorD(fitParam), isec);
519 if ( fitCovArr ) fitCovArr->AddAtAndExpand(new TMatrixD(covMatrix), isec);
523 //_____________________________________________________________________________
524 TObjArray* AliTPCCalPad::CreateFormulaArray(const char *fitFormula)
527 // create an array of TFormulas for the each parameter of the fit function
530 // split fit string in single parameters
531 // find dimension of the fit:
532 TString fitString(fitFormula);
533 fitString.ReplaceAll("++","#");
534 fitString.ReplaceAll(" ","");
535 TObjArray *arrFitParams = fitString.Tokenize("#");
536 Int_t ndim = arrFitParams->GetEntries();
537 //create array of TFormulas to evaluate the parameters
538 TObjArray *arrFitFormulas = new TObjArray(ndim);
539 arrFitFormulas->SetOwner(kTRUE);
540 for (Int_t idim=0;idim<ndim;++idim){
541 TString s=((TObjString*)arrFitParams->At(idim))->GetString();
542 s.ReplaceAll("gx","[0]");
543 s.ReplaceAll("gy","[1]");
544 s.ReplaceAll("lx","[2]");
545 s.ReplaceAll("ly","[3]");
546 s.ReplaceAll("sector","[4]");
547 arrFitFormulas->AddAt(new TFormula(Form("param%02d",idim),s.Data()),idim);
551 return arrFitFormulas;
553 //_____________________________________________________________________________
554 void AliTPCCalPad::EvalFormulaArray(const TObjArray &arrFitFormulas, TVectorD &results,
555 const Int_t sec, const Int_t row, const Int_t pad)
558 // evaluate the fit formulas
560 Int_t ndim=arrFitFormulas.GetEntries();
561 results.ResizeTo(ndim);
563 AliTPCROC* tpcROCinstance = AliTPCROC::Instance(); // to calculate the pad's position
565 Float_t globalXYZ[3];
566 tpcROCinstance->GetPositionLocal(sec, row, pad, localXYZ);
567 tpcROCinstance->GetPositionGlobal(sec, row, pad, globalXYZ);
568 //calculate parameter values
569 for (Int_t idim=0;idim<ndim;++idim){
570 TFormula *f=(TFormula*)arrFitFormulas.At(idim);
571 f->SetParameters(globalXYZ[0],globalXYZ[1],localXYZ[0],localXYZ[1],sec);
572 results[idim]=f->Eval(0);
575 //_____________________________________________________________________________
576 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){
578 // Performs a fit on both sides.
579 // Valid information for the fitFormula are the variables
580 // - gx, gy, lx ,ly: meaning global x, global y, local x, local y value of the padName
581 // - sector: the sector number.
582 // eg. a formula might look 'gy' or '(sector<36) ++ gy' or 'gx ++ gy' or 'gx ++ gy ++ lx ++ lx^2' and so on
584 // PadOutliers - pads with value !=0 are not used in fitting procedure
585 // chi2Threshold: Threshold for chi2 when EvalRobust is called
586 // robustFraction: Fraction of data that will be used in EvalRobust
589 TObjArray* arrFitFormulas=CreateFormulaArray(fitFormula);
590 Int_t ndim = arrFitFormulas->GetEntries();
591 //resize output data arrays
592 fitParamSideA.ResizeTo(ndim+1);
593 fitParamSideC.ResizeTo(ndim+1);
594 covMatrixSideA.ResizeTo(ndim+1,ndim+1);
595 covMatrixSideC.ResizeTo(ndim+1,ndim+1);
596 // create linear fitter for A- and C- Side
597 TLinearFitter* fitterGA = new TLinearFitter(ndim+1,Form("hyp%d",ndim));
598 TLinearFitter* fitterGC = new TLinearFitter(ndim+1,Form("hyp%d",ndim));
599 fitterGA->StoreData(kTRUE);
600 fitterGC->StoreData(kTRUE);
602 TVectorD parValues(ndim);
604 AliTPCCalROC *rocErr=0x0;
606 for (UInt_t isec = 0; isec<kNsec; ++isec){
607 AliTPCCalROC *rocOut=PadOutliers->GetCalROC(isec);
608 AliTPCCalROC *rocData=GetCalROC(isec);
609 if (pointError) rocErr=pointError->GetCalROC(isec);
610 if (!rocData) continue;
611 for (UInt_t irow = 0; irow < GetCalROC(isec)->GetNrows(); irow++) {
612 for (UInt_t ipad = 0; ipad < GetCalROC(isec)->GetNPads(irow); ipad++) {
614 if (rocOut && rocOut->GetValue(irow,ipad)) continue;
615 //calculate parameter values
616 EvalFormulaArray(*arrFitFormulas,parValues,isec,irow,ipad);
618 Float_t value=rocData->GetValue(irow,ipad);
622 err=TMath::Nint(rocErr->GetValue(irow,ipad));
625 //add points to the fitters
627 fitterGA->AddPoint(parValues.GetMatrixArray(),value,err);
629 fitterGC->AddPoint(parValues.GetMatrixArray(),value,err);
635 fitterGA->EvalRobust(robustFraction);
636 fitterGC->EvalRobust(robustFraction);
641 chi2SideA=fitterGA->GetChisquare()/(fitterGA->GetNpoints()-(ndim+1));
642 chi2SideC=fitterGC->GetChisquare()/(fitterGC->GetNpoints()-(ndim+1));
643 fitterGA->GetParameters(fitParamSideA);
644 fitterGC->GetParameters(fitParamSideC);
645 fitterGA->GetCovarianceMatrix(covMatrixSideA);
646 fitterGC->GetCovarianceMatrix(covMatrixSideC);
648 delete arrFitFormulas;
654 AliTPCCalPad *AliTPCCalPad::CreateCalPadFit(const char* fitFormula, const TVectorD &fitParamSideA, const TVectorD &fitParamSideC)
659 TObjArray *arrFitFormulas=CreateFormulaArray(fitFormula);
660 Int_t ndim = arrFitFormulas->GetEntries();
661 //check if dimension of fit formula and fit parameters agree
662 if (ndim!=fitParamSideA.GetNrows()||ndim!=fitParamSideC.GetNrows()){
663 printf("AliTPCCalPad::CreateCalPadFit: Dimensions of fit formula and fit Parameters does not match!");
667 AliTPCCalPad *pad=new AliTPCCalPad("fitResultPad",Form("Fit result: %s",fitFormula));
668 //fill cal pad with fit results if requested
669 for (UInt_t isec = 0; isec<kNsec; ++isec){
670 AliTPCCalROC *roc=pad->GetCalROC(isec);
671 for (UInt_t irow = 0; irow < roc->GetNrows(); irow++) {
672 for (UInt_t ipad = 0; ipad < roc->GetNPads(irow); ipad++) {
673 const TVectorD *fitPar=0;
674 TVectorD fitResArray;
676 fitPar=&fitParamSideA;
678 fitPar=&fitParamSideC;
680 EvalFormulaArray(*arrFitFormulas,fitResArray, isec, irow, ipad);
681 for (Int_t idim=0;idim<ndim;++idim)
682 fitResArray(idim)*=(*fitPar)(idim);
683 roc->SetValue(irow,ipad,fitResArray.Sum());
687 delete arrFitFormulas;
693 TCanvas * AliTPCCalPad::MakeReportPadSector(TTree *chain, const char* varName, const char*varTitle, const char *axisTitle, Float_t min, Float_t max, const char *cutUser){
695 // Make a report - cal pads per sector
696 // mean valeus per sector and local X
700 TCanvas *canvas = new TCanvas(Form("Sector: %s",varTitle),Form("Sector: %s",varTitle),1500,1100);
703 chain->SetAlias("lX","lx.fElements");
706 TString strDraw=varName;
708 legend = new TLegend(0.5,0.50,0.9,0.9, Form("%s TPC A side", varTitle));
709 for (Int_t isec=-1; isec<18; isec+=1){
710 TCut cutSec=Form("sector%%36==%d",isec);
712 if (isec==-1) cutSec="sector%36<18";
713 chain->SetMarkerColor(1+(isec+2)%5);
714 chain->SetLineColor(1+(isec+2)%5);
715 chain->SetMarkerStyle(25+(isec+2)%4);
717 chain->Draw(strDraw.Data(),cutSec,"profgoff");
718 his=(TH1*)chain->GetHistogram()->Clone();
719 delete chain->GetHistogram();
720 his->SetMaximum(max);
721 his->SetMinimum(min);
722 his->GetXaxis()->SetTitle("R (cm)");
723 his->GetYaxis()->SetTitle(axisTitle);
724 his->SetTitle(Form("%s- sector %d",varTitle, isec));
725 his->SetName(Form("%s- sector %d",varTitle, isec));
726 if (isec==-1) his->SetTitle(Form("%s A side",varTitle));
727 if (isec==-1) his->Draw();
729 legend->AddEntry(his);
734 legend = new TLegend(0.5,0.50,0.9,0.9, Form("%s TPC C side", varTitle));
735 for (Int_t isec=-1; isec<18; isec+=1){
736 TCut cutSec=Form("(sector+18)%%36==%d",isec);
738 if (isec==-1) cutSec="sector%36>18";
739 chain->SetMarkerColor(1+(isec+2)%5);
740 chain->SetLineColor(1+(isec+2)%5);
741 chain->SetMarkerStyle(25+isec%4);
743 chain->Draw(strDraw.Data(),cutSec,"profgoff");
744 his=(TH1*)chain->GetHistogram()->Clone();
745 delete chain->GetHistogram();
746 his->SetMaximum(max);
747 his->SetMinimum(min);
748 his->GetXaxis()->SetTitle("R (cm)");
749 his->GetYaxis()->SetTitle(axisTitle);
750 his->SetTitle(Form("%s- sector %d",varTitle,isec));
751 his->SetName(Form("%s- sector %d",varTitle,isec));
752 if (isec==-1) his->SetTitle(Form("%s C side",varTitle));
753 if (isec==-1) his->Draw();
755 legend->AddEntry(his);
764 TCanvas * AliTPCCalPad::MakeReportPadSector2D(TTree *chain, const char* varName, const char*varTitle, const char *axisTitle, Float_t min, Float_t max, const char *cutUser){
766 // Make a report - cal pads per sector
768 // Input tree should be created using AliPreprocesorOnline before
771 TCanvas *canvas = new TCanvas(Form("%s2D",varTitle),Form("%s2D",varTitle),1500,1100);
774 TString strDraw=varName;
775 strDraw+=":gy.fElements:gx.fElements>>his(250,-250,250,250,-250,250)";
779 pad->SetMargin(0.15,0.15,0.15,0.15);
781 chain->Draw(strDraw.Data(),"sector%36<18"+cut,"profgoffcolz2");
782 his=(TH1*)chain->GetHistogram()->Clone();
783 delete chain->GetHistogram();
784 his->SetMaximum(max);
785 his->SetMinimum(min);
786 his->GetXaxis()->SetTitle("x (cm)");
787 his->GetYaxis()->SetTitle("y (cm)");
788 his->GetZaxis()->SetTitle(axisTitle);
789 his->SetTitle(Form("%s A side",varTitle));
790 his->SetName(Form("%s A side",varTitle));
794 pad->SetMargin(0.15,0.15,0.15,0.15);
796 chain->Draw(strDraw.Data(),"sector%36>=18"+cut,"profgoffcolz2");
797 his=(TH1*)chain->GetHistogram()->Clone();
798 delete chain->GetHistogram();
799 his->SetMaximum(max);
800 his->SetMinimum(min);
801 his->GetXaxis()->SetTitle("x (cm)");
802 his->GetYaxis()->SetTitle("y (cm)");
803 his->GetZaxis()->SetTitle(axisTitle);
804 his->SetTitle(Form("%s C side",varTitle));
805 his->SetName(Form("%s C side",varTitle));
812 void AliTPCCalPad::Draw(Option_t* option){
814 // Draw function - standard 2D view
817 TCanvas *canvas = new TCanvas(Form("%s2D",GetTitle()),Form("%s2D",GetTitle()),900,900);
823 pad->SetMargin(0.15,0.15,0.15,0.15);
825 his->GetXaxis()->SetTitle("x (cm)");
826 his->GetYaxis()->SetTitle("y (cm)");
827 his->GetZaxis()->SetTitle(GetTitle());
828 his->SetTitle(Form("%s A side",GetTitle()));
829 his->SetName(Form("%s A side",GetTitle()));
833 pad->SetMargin(0.15,0.15,0.15,0.15);
835 his->GetXaxis()->SetTitle("x (cm)");
836 his->GetYaxis()->SetTitle("y (cm)");
837 his->GetZaxis()->SetTitle(GetTitle());
838 his->SetTitle(Form("%s C side",GetTitle()));
839 his->SetName(Form("%s C side",GetTitle()));
843 pad->SetMargin(0.15,0.15,0.15,0.15);
844 his=MakeHisto1D(-8,8,0,1);
845 his->GetXaxis()->SetTitle(GetTitle());
846 his->SetTitle(Form("%s A side",GetTitle()));
847 his->SetName(Form("%s A side",GetTitle()));
851 pad->SetMargin(0.15,0.15,0.15,0.15);
852 his=MakeHisto1D(-8,8,0,-1);
853 his->GetXaxis()->SetTitle(GetTitle());
854 his->SetTitle(Form("%s C side",GetTitle()));
855 his->SetName(Form("%s C side",GetTitle()));
862 AliTPCCalPad * AliTPCCalPad::MakeCalPadFromHistoRPHI(TH2 * hisA, TH2* hisC){
864 // Make cal pad from r-phi histograms
866 AliTPCROC *proc= AliTPCROC::Instance();
867 AliTPCCalPad *calPad = new AliTPCCalPad("his","his");
868 Float_t globalPos[3];
869 for (Int_t isec=0; isec<72; isec++){
870 AliTPCCalROC* calRoc = calPad->GetCalROC(isec);
871 TH2 * his = ((isec%36<18) ? hisA:hisC);
872 for (UInt_t irow=0; irow<calRoc->GetNrows(); irow+=1){
874 if (isec>=36) jrow+=63;
875 for (UInt_t ipad=0;ipad<proc->GetNPads(isec,irow);ipad+=1){
876 proc->GetPositionGlobal(isec,irow,ipad, globalPos);
877 Double_t phi=TMath::ATan2(globalPos[1],globalPos[0]);
878 //if (phi<0) phi+=TMath::Pi()*2;
879 Int_t bin=his->FindBin(phi,jrow);
880 Float_t value= his->GetBinContent(bin);
881 calRoc->SetValue(irow,ipad,value);