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>
48 #include "AliTPCPreprocessorOnline.h"
49 #include "AliTPCCalibViewer.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));
162 Bool_t AliTPCCalPad::MedianFilter(Int_t deltaRow, Int_t deltaPad, AliTPCCalPad*outlierPad, Bool_t doEdge){
164 // replace constent with median in the neigborhood
167 for (Int_t isec = 0; isec < kNsec; isec++) {
168 AliTPCCalROC *outlierROC=(outlierPad==NULL)?NULL:outlierPad->GetCalROC(isec);
170 isOK&=fROC[isec]->MedianFilter(deltaRow,deltaPad,outlierROC,doEdge);
176 Bool_t AliTPCCalPad::LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type, AliTPCCalPad*outlierPad, Bool_t doEdge){
178 // replace constent with LTM statistic in neigborhood
181 for (Int_t isec = 0; isec < kNsec; isec++) {
182 AliTPCCalROC *outlierROC=(outlierPad==NULL)?NULL:outlierPad->GetCalROC(isec);
184 isOK&=fROC[isec]->LTMFilter(deltaRow, deltaPad,fraction,type,outlierROC,doEdge);
190 Bool_t AliTPCCalPad::Convolute(Double_t sigmaPad, Double_t sigmaRow, AliTPCCalPad*outlierPad, TF1 *fpad, TF1 *frow){
192 // replace constent with median in the neigborhood
195 for (Int_t isec = 0; isec < kNsec; isec++) {
196 AliTPCCalROC *outlierROC=(outlierPad==NULL)?NULL:outlierPad->GetCalROC(isec);
198 isOK&=fROC[isec]->Convolute(sigmaPad,sigmaRow,outlierROC,fpad,frow);
205 //_____________________________________________________________________________
206 void AliTPCCalPad::Add(Float_t c1)
209 // add constant c1 to all channels of all ROCs
212 for (Int_t isec = 0; isec < kNsec; isec++) {
219 //_____________________________________________________________________________
220 void AliTPCCalPad::Multiply(Float_t c1)
223 // multiply each channel of all ROCs with c1
225 for (Int_t isec = 0; isec < kNsec; isec++) {
227 fROC[isec]->Multiply(c1);
232 //_____________________________________________________________________________
233 void AliTPCCalPad::Add(const AliTPCCalPad * pad, Double_t c1)
236 // multiply AliTPCCalPad 'pad' by c1 and add each channel to the coresponing channel in all ROCs
239 for (Int_t isec = 0; isec < kNsec; isec++) {
240 if (fROC[isec] && pad->GetCalROC(isec)){
241 fROC[isec]->Add(pad->GetCalROC(isec),c1);
246 //_____________________________________________________________________________
247 void AliTPCCalPad::Multiply(const AliTPCCalPad * pad)
250 // multiply each channel of all ROCs with the coresponding channel of 'pad'
253 for (Int_t isec = 0; isec < kNsec; isec++) {
255 fROC[isec]->Multiply(pad->GetCalROC(isec));
260 //_____________________________________________________________________________
261 void AliTPCCalPad::Divide(const AliTPCCalPad * pad)
264 // divide each channel of all ROCs by the coresponding channel of 'pad'
267 for (Int_t isec = 0; isec < kNsec; isec++) {
269 fROC[isec]->Divide(pad->GetCalROC(isec));
274 //_____________________________________________________________________________
275 TGraph * AliTPCCalPad::MakeGraph(Int_t type, Float_t ratio){
282 for (Int_t i=0;i<72;i++) if (fROC[i]) npoints++;
283 TGraph * graph = new TGraph(npoints);
285 for (Int_t isec=0;isec<72;isec++){
286 if (!fROC[isec]) continue;
287 if (type==0) graph->SetPoint(npoints,isec,fROC[isec]->GetMean());
288 if (type==1) graph->SetPoint(npoints,isec,fROC[isec]->GetMedian());
289 if (type==2) graph->SetPoint(npoints,isec,fROC[isec]->GetLTM(0,ratio));
293 graph->GetXaxis()->SetTitle("Sector");
295 graph->GetYaxis()->SetTitle("Mean");
296 graph->SetMarkerStyle(22);
299 graph->GetYaxis()->SetTitle("Median");
300 graph->SetMarkerStyle(22);
303 graph->GetYaxis()->SetTitle(Form("Mean%f",ratio));
304 graph->SetMarkerStyle(24);
310 //_____________________________________________________________________________
311 Double_t AliTPCCalPad::GetMeanRMS(Double_t &rms)
314 // Calculates mean and RMS of all ROCs
316 Double_t sum = 0, sum2 = 0, n=0, val=0;
317 for (Int_t isec = 0; isec < kNsec; isec++) {
318 AliTPCCalROC *calRoc = fROC[isec];
320 for (UInt_t irow=0; irow<calRoc->GetNrows(); irow++){
321 for (UInt_t ipad=0; ipad<calRoc->GetNPads(irow); ipad++){
322 val = calRoc->GetValue(irow,ipad);
332 Double_t mean = sum*n1;
333 rms = TMath::Sqrt(TMath::Abs(sum2*n1-mean*mean));
338 //_____________________________________________________________________________
339 Double_t AliTPCCalPad::GetMean(AliTPCCalPad* outlierPad)
342 // return mean of the mean of all ROCs
346 for (Int_t isec = 0; isec < kNsec; isec++) {
347 AliTPCCalROC *calRoc = fROC[isec];
349 AliTPCCalROC* outlierROC = 0;
350 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
351 arr[n] = calRoc->GetMean(outlierROC);
355 return TMath::Mean(n,arr);
358 //_____________________________________________________________________________
359 Double_t AliTPCCalPad::GetRMS(AliTPCCalPad* outlierPad)
362 // return mean of the RMS of all ROCs
366 for (Int_t isec = 0; isec < kNsec; isec++) {
367 AliTPCCalROC *calRoc = fROC[isec];
369 AliTPCCalROC* outlierROC = 0;
370 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
371 arr[n] = calRoc->GetRMS(outlierROC);
375 return TMath::Mean(n,arr);
378 //_____________________________________________________________________________
379 Double_t AliTPCCalPad::GetMedian(AliTPCCalPad* outlierPad)
382 // return mean of the median of all ROCs
386 for (Int_t isec = 0; isec < kNsec; isec++) {
387 AliTPCCalROC *calRoc = fROC[isec];
389 AliTPCCalROC* outlierROC = 0;
390 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
391 arr[n] = calRoc->GetMedian(outlierROC);
395 return TMath::Mean(n,arr);
398 //_____________________________________________________________________________
399 Double_t AliTPCCalPad::GetLTM(Double_t *sigma, Double_t fraction, AliTPCCalPad* outlierPad)
402 // return mean of the LTM and sigma of all ROCs
404 Double_t arrm[kNsec];
405 Double_t arrs[kNsec];
409 for (Int_t isec = 0; isec < kNsec; isec++) {
410 AliTPCCalROC *calRoc = fROC[isec];
412 if ( sigma ) sTemp=arrs+n;
413 AliTPCCalROC* outlierROC = 0;
414 if (outlierPad) outlierROC = outlierPad->GetCalROC(isec);
415 arrm[n] = calRoc->GetLTM(sTemp,fraction, outlierROC);
419 if ( sigma ) *sigma = TMath::Mean(n,arrs);
420 return TMath::Mean(n,arrm);
423 //_____________________________________________________________________________
424 TH1F * AliTPCCalPad::MakeHisto1D(Float_t min, Float_t max,Int_t type, Int_t side){
427 // type -1 = user defined range
428 // 0 = nsigma cut nsigma=min
433 Float_t mean = GetMean();
434 Float_t sigma = GetRMS();
435 Float_t nsigma = TMath::Abs(min);
436 min = mean-nsigma*sigma;
437 max = mean+nsigma*sigma;
441 Float_t mean = GetMedian();
448 // LTM mean +- nsigma
451 Float_t mean = GetLTM(&sigma,max);
457 TString name=Form("%s Pad 1D",GetTitle());
458 TH1F * his = new TH1F(name.Data(),name.Data(),100, min,max);
459 for (Int_t isec = 0; isec < kNsec; isec++) {
460 if (side==1 && isec%36>18) continue;
461 if (side==-1 && isec%36<18) continue;
463 for (UInt_t irow=0; irow<fROC[isec]->GetNrows(); irow++){
464 UInt_t npads = (Int_t)fROC[isec]->GetNPads(irow);
465 for (UInt_t ipad=0; ipad<npads; ipad++){
466 his->Fill(fROC[isec]->GetValue(irow,ipad));
474 //_____________________________________________________________________________
475 TH2F *AliTPCCalPad::MakeHisto2D(Int_t side){
478 // side - specify the side A = 0 C = 1
479 // type - used types of determination of boundaries in z
481 Float_t kEpsilon = 0.000000000001;
482 TH2F * his = new TH2F(GetName(), GetName(), 250,-250,250,250,-250,250);
483 AliTPCROC * roc = AliTPCROC::Instance();
484 for (Int_t isec=0; isec<72; isec++){
485 if (side==0 && isec%36>=18) continue;
486 if (side>0 && isec%36<18) continue;
488 AliTPCCalROC * calRoc = fROC[isec];
489 for (UInt_t irow=0; irow<calRoc->GetNrows(); irow++)
490 for (UInt_t ipad=0; ipad<calRoc->GetNPads(irow); ipad++)
491 if (TMath::Abs(calRoc->GetValue(irow,ipad))>kEpsilon){
493 roc->GetPositionGlobal(isec,irow,ipad,xyz);
494 Int_t binx = 1+TMath::Nint((xyz[0]+250.)*0.5);
495 Int_t biny = 1+TMath::Nint((xyz[1]+250.)*0.5);
496 Float_t value = calRoc->GetValue(irow,ipad);
497 his->SetBinContent(binx,biny,value);
501 his->SetXTitle("x (cm)");
502 his->SetYTitle("y (cm)");
507 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 {
509 // Loops over all AliTPCCalROCs and performs a localFit in each ROC
510 // AliTPCCalPad with fit-data is returned
511 // rowRadius and padRadius specifies a window around a given pad in one sector.
512 // The data of this window are fitted with a parabolic function.
513 // This function is evaluated at the pad's position.
514 // At the edges the window is shifted, so that the specified pad is not anymore in the center of the window.
515 // rowRadius - radius - rows to be used for smoothing
516 // padradius - radius - pads to be used for smoothing
517 // ROCoutlier - map of outliers - pads not to be used for local smoothing
518 // robust - robust method of fitting - (much slower)
519 // chi2Threshold: Threshold for chi2 when EvalRobust is called
520 // robustFraction: Fraction of data that will be used in EvalRobust
523 AliTPCCalPad* pad = new AliTPCCalPad(padName, padName);
524 for (Int_t isec = 0; isec < 72; isec++){
525 if (printCurrentSector) std::cout << "LocalFit in sector " << isec << "\r" << std::flush;
527 pad->SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, PadOutliers->GetCalROC(isec), robust, chi2Threshold, robustFraction));
529 pad->SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, 0, robust, chi2Threshold, robustFraction));
535 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){
537 // Loops over all AliTPCCalROCs and performs a globalFit in each ROC
538 // AliTPCCalPad with fit-data is returned
539 // chi2Threshold: Threshold for chi2 when EvalRobust is called
540 // robustFraction: Fraction of data that will be used in EvalRobust
541 // chi2Threshold: Threshold for chi2 when EvalRobust is called
542 // robustFraction: Fraction of data that will be used in EvalRobust
543 // err: error of the data points
544 // if fitParArr and/or fitCovArr is given, write fitParameters and/or covariance Matrices into the array
546 AliTPCCalPad* pad = new AliTPCCalPad(padName, padName);
547 TVectorD fitParam(0);
548 TMatrixD covMatrix(0,0);
550 for (Int_t isec = 0; isec < 72; isec++){
552 GetCalROC(isec)->GlobalFit(PadOutliers->GetCalROC(isec), robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction, err);
554 GetCalROC(isec)->GlobalFit(0, robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction, err);
556 AliTPCCalROC *roc=AliTPCCalROC::CreateGlobalFitCalROC(fitParam, isec);
559 if ( fitParArr ) fitParArr->AddAtAndExpand(new TVectorD(fitParam), isec);
560 if ( fitCovArr ) fitCovArr->AddAtAndExpand(new TMatrixD(covMatrix), isec);
564 //_____________________________________________________________________________
565 TObjArray* AliTPCCalPad::CreateFormulaArray(const char *fitFormula)
568 // create an array of TFormulas for the each parameter of the fit function
571 // split fit string in single parameters
572 // find dimension of the fit:
573 TString fitString(fitFormula);
574 fitString.ReplaceAll("++","#");
575 fitString.ReplaceAll(" ","");
576 TObjArray *arrFitParams = fitString.Tokenize("#");
577 Int_t ndim = arrFitParams->GetEntries();
578 //create array of TFormulas to evaluate the parameters
579 TObjArray *arrFitFormulas = new TObjArray(ndim);
580 arrFitFormulas->SetOwner(kTRUE);
581 for (Int_t idim=0;idim<ndim;++idim){
582 TString s=((TObjString*)arrFitParams->At(idim))->GetString();
583 s.ReplaceAll("gx","[0]");
584 s.ReplaceAll("gy","[1]");
585 s.ReplaceAll("lx","[2]");
586 s.ReplaceAll("ly","[3]");
587 s.ReplaceAll("sector","[4]");
588 arrFitFormulas->AddAt(new TFormula(Form("param%02d",idim),s.Data()),idim);
592 return arrFitFormulas;
594 //_____________________________________________________________________________
595 void AliTPCCalPad::EvalFormulaArray(const TObjArray &arrFitFormulas, TVectorD &results,
596 const Int_t sec, const Int_t row, const Int_t pad)
599 // evaluate the fit formulas
601 Int_t ndim=arrFitFormulas.GetEntries();
602 results.ResizeTo(ndim);
604 AliTPCROC* tpcROCinstance = AliTPCROC::Instance(); // to calculate the pad's position
606 Float_t globalXYZ[3];
607 tpcROCinstance->GetPositionLocal(sec, row, pad, localXYZ);
608 tpcROCinstance->GetPositionGlobal(sec, row, pad, globalXYZ);
609 //calculate parameter values
610 for (Int_t idim=0;idim<ndim;++idim){
611 TFormula *f=(TFormula*)arrFitFormulas.At(idim);
612 f->SetParameters(globalXYZ[0],globalXYZ[1],localXYZ[0],localXYZ[1],sec);
613 results[idim]=f->Eval(0);
616 //_____________________________________________________________________________
617 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){
619 // Performs a fit on both sides.
620 // Valid information for the fitFormula are the variables
621 // - gx, gy, lx ,ly: meaning global x, global y, local x, local y value of the padName
622 // - sector: the sector number.
623 // eg. a formula might look 'gy' or '(sector<36) ++ gy' or 'gx ++ gy' or 'gx ++ gy ++ lx ++ lx^2' and so on
625 // PadOutliers - pads with value !=0 are not used in fitting procedure
626 // chi2Threshold: Threshold for chi2 when EvalRobust is called
627 // robustFraction: Fraction of data that will be used in EvalRobust
630 TObjArray* arrFitFormulas=CreateFormulaArray(fitFormula);
631 Int_t ndim = arrFitFormulas->GetEntries();
632 //resize output data arrays
633 fitParamSideA.ResizeTo(ndim+1);
634 fitParamSideC.ResizeTo(ndim+1);
635 covMatrixSideA.ResizeTo(ndim+1,ndim+1);
636 covMatrixSideC.ResizeTo(ndim+1,ndim+1);
637 // create linear fitter for A- and C- Side
638 TLinearFitter* fitterGA = new TLinearFitter(ndim+1,Form("hyp%d",ndim));
639 TLinearFitter* fitterGC = new TLinearFitter(ndim+1,Form("hyp%d",ndim));
640 fitterGA->StoreData(kTRUE);
641 fitterGC->StoreData(kTRUE);
643 TVectorD parValues(ndim);
645 AliTPCCalROC *rocErr=0x0;
647 for (UInt_t isec = 0; isec<kNsec; ++isec){
648 AliTPCCalROC *rocOut=PadOutliers->GetCalROC(isec);
649 AliTPCCalROC *rocData=GetCalROC(isec);
650 if (pointError) rocErr=pointError->GetCalROC(isec);
651 if (!rocData) continue;
652 for (UInt_t irow = 0; irow < GetCalROC(isec)->GetNrows(); irow++) {
653 for (UInt_t ipad = 0; ipad < GetCalROC(isec)->GetNPads(irow); ipad++) {
655 if (rocOut && rocOut->GetValue(irow,ipad)) continue;
656 //calculate parameter values
657 EvalFormulaArray(*arrFitFormulas,parValues,isec,irow,ipad);
659 Float_t value=rocData->GetValue(irow,ipad);
663 err=TMath::Nint(rocErr->GetValue(irow,ipad));
666 //add points to the fitters
668 fitterGA->AddPoint(parValues.GetMatrixArray(),value,err);
670 fitterGC->AddPoint(parValues.GetMatrixArray(),value,err);
676 fitterGA->EvalRobust(robustFraction);
677 fitterGC->EvalRobust(robustFraction);
682 chi2SideA=fitterGA->GetChisquare()/(fitterGA->GetNpoints()-(ndim+1));
683 chi2SideC=fitterGC->GetChisquare()/(fitterGC->GetNpoints()-(ndim+1));
684 fitterGA->GetParameters(fitParamSideA);
685 fitterGC->GetParameters(fitParamSideC);
686 fitterGA->GetCovarianceMatrix(covMatrixSideA);
687 fitterGC->GetCovarianceMatrix(covMatrixSideC);
689 delete arrFitFormulas;
695 AliTPCCalPad *AliTPCCalPad::CreateCalPadFit(const char* fitFormula, const TVectorD &fitParamSideA, const TVectorD &fitParamSideC)
700 TObjArray *arrFitFormulas=CreateFormulaArray(fitFormula);
701 Int_t ndim = arrFitFormulas->GetEntries();
702 //check if dimension of fit formula and fit parameters agree
703 if (ndim!=fitParamSideA.GetNrows()||ndim!=fitParamSideC.GetNrows()){
704 printf("AliTPCCalPad::CreateCalPadFit: Dimensions of fit formula and fit Parameters does not match!");
708 AliTPCCalPad *pad=new AliTPCCalPad("fitResultPad",Form("Fit result: %s",fitFormula));
709 //fill cal pad with fit results if requested
710 for (UInt_t isec = 0; isec<kNsec; ++isec){
711 AliTPCCalROC *roc=pad->GetCalROC(isec);
712 for (UInt_t irow = 0; irow < roc->GetNrows(); irow++) {
713 for (UInt_t ipad = 0; ipad < roc->GetNPads(irow); ipad++) {
714 const TVectorD *fitPar=0;
715 TVectorD fitResArray;
717 fitPar=&fitParamSideA;
719 fitPar=&fitParamSideC;
721 EvalFormulaArray(*arrFitFormulas,fitResArray, isec, irow, ipad);
722 for (Int_t idim=0;idim<ndim;++idim)
723 fitResArray(idim)*=(*fitPar)(idim);
724 roc->SetValue(irow,ipad,fitResArray.Sum());
728 delete arrFitFormulas;
734 TCanvas * AliTPCCalPad::MakeReportPadSector(TTree *chain, const char* varName, const char*varTitle, const char *axisTitle, Float_t min, Float_t max, const char *cutUser){
736 // Make a report - cal pads per sector
737 // mean valeus per sector and local X
741 TCanvas *canvas = new TCanvas(Form("Sector: %s",varTitle),Form("Sector: %s",varTitle),1500,1100);
744 chain->SetAlias("lX","lx.fElements");
747 TString strDraw=varName;
749 legend = new TLegend(0.5,0.50,0.9,0.9, Form("%s TPC A side", varTitle));
750 for (Int_t isec=-1; isec<18; isec+=1){
751 TCut cutSec=Form("sector%%36==%d",isec);
753 if (isec==-1) cutSec="sector%36<18";
754 chain->SetMarkerColor(1+(isec+2)%5);
755 chain->SetLineColor(1+(isec+2)%5);
756 chain->SetMarkerStyle(25+(isec+2)%4);
758 chain->Draw(strDraw.Data(),cutSec,"profgoff");
759 his=(TH1*)chain->GetHistogram()->Clone();
760 delete chain->GetHistogram();
761 his->SetMaximum(max);
762 his->SetMinimum(min);
763 his->GetXaxis()->SetTitle("R (cm)");
764 his->GetYaxis()->SetTitle(axisTitle);
765 his->SetTitle(Form("%s- sector %d",varTitle, isec));
766 his->SetName(Form("%s- sector %d",varTitle, isec));
767 if (isec==-1) his->SetTitle(Form("%s A side",varTitle));
768 if (isec==-1) his->Draw();
770 legend->AddEntry(his);
775 legend = new TLegend(0.5,0.50,0.9,0.9, Form("%s TPC C side", varTitle));
776 for (Int_t isec=-1; isec<18; isec+=1){
777 TCut cutSec=Form("(sector+18)%%36==%d",isec);
779 if (isec==-1) cutSec="sector%36>18";
780 chain->SetMarkerColor(1+(isec+2)%5);
781 chain->SetLineColor(1+(isec+2)%5);
782 chain->SetMarkerStyle(25+isec%4);
784 chain->Draw(strDraw.Data(),cutSec,"profgoff");
785 his=(TH1*)chain->GetHistogram()->Clone();
786 delete chain->GetHistogram();
787 his->SetMaximum(max);
788 his->SetMinimum(min);
789 his->GetXaxis()->SetTitle("R (cm)");
790 his->GetYaxis()->SetTitle(axisTitle);
791 his->SetTitle(Form("%s- sector %d",varTitle,isec));
792 his->SetName(Form("%s- sector %d",varTitle,isec));
793 if (isec==-1) his->SetTitle(Form("%s C side",varTitle));
794 if (isec==-1) his->Draw();
796 legend->AddEntry(his);
805 TCanvas * AliTPCCalPad::MakeReportPadSector2D(TTree *chain, const char* varName, const char*varTitle, const char *axisTitle, Float_t min, Float_t max, const char *cutUser){
807 // Make a report - cal pads per sector
809 // Input tree should be created using AliPreprocesorOnline before
812 TCanvas *canvas = new TCanvas(Form("%s2D",varTitle),Form("%s2D",varTitle),1500,1100);
815 TString strDraw=varName;
816 strDraw+=":gy.fElements:gx.fElements>>his(250,-250,250,250,-250,250)";
820 pad->SetMargin(0.15,0.15,0.15,0.15);
822 chain->Draw(strDraw.Data(),"sector%36<18"+cut,"profgoffcolz2");
823 his=(TH1*)chain->GetHistogram()->Clone();
824 delete chain->GetHistogram();
825 his->SetMaximum(max);
826 his->SetMinimum(min);
827 his->GetXaxis()->SetTitle("x (cm)");
828 his->GetYaxis()->SetTitle("y (cm)");
829 his->GetZaxis()->SetTitle(axisTitle);
830 his->SetTitle(Form("%s A side",varTitle));
831 his->SetName(Form("%s A side",varTitle));
835 pad->SetMargin(0.15,0.15,0.15,0.15);
837 chain->Draw(strDraw.Data(),"sector%36>=18"+cut,"profgoffcolz2");
838 his=(TH1*)chain->GetHistogram()->Clone();
839 delete chain->GetHistogram();
840 his->SetMaximum(max);
841 his->SetMinimum(min);
842 his->GetXaxis()->SetTitle("x (cm)");
843 his->GetYaxis()->SetTitle("y (cm)");
844 his->GetZaxis()->SetTitle(axisTitle);
845 his->SetTitle(Form("%s C side",varTitle));
846 his->SetName(Form("%s C side",varTitle));
853 void AliTPCCalPad::Draw(Option_t* option){
855 // Draw function - standard 2D view
858 TCanvas *canvas = new TCanvas(Form("%s2D",GetTitle()),Form("%s2D",GetTitle()),900,900);
864 pad->SetMargin(0.15,0.15,0.15,0.15);
866 his->GetXaxis()->SetTitle("x (cm)");
867 his->GetYaxis()->SetTitle("y (cm)");
868 his->GetZaxis()->SetTitle(GetTitle());
869 his->SetTitle(Form("%s A side",GetTitle()));
870 his->SetName(Form("%s A side",GetTitle()));
874 pad->SetMargin(0.15,0.15,0.15,0.15);
876 his->GetXaxis()->SetTitle("x (cm)");
877 his->GetYaxis()->SetTitle("y (cm)");
878 his->GetZaxis()->SetTitle(GetTitle());
879 his->SetTitle(Form("%s C side",GetTitle()));
880 his->SetName(Form("%s C side",GetTitle()));
884 pad->SetMargin(0.15,0.15,0.15,0.15);
885 his=MakeHisto1D(-8,8,0,1);
886 his->GetXaxis()->SetTitle(GetTitle());
887 his->SetTitle(Form("%s A side",GetTitle()));
888 his->SetName(Form("%s A side",GetTitle()));
892 pad->SetMargin(0.15,0.15,0.15,0.15);
893 his=MakeHisto1D(-8,8,0,-1);
894 his->GetXaxis()->SetTitle(GetTitle());
895 his->SetTitle(Form("%s C side",GetTitle()));
896 his->SetName(Form("%s C side",GetTitle()));
903 AliTPCCalPad * AliTPCCalPad::MakeCalPadFromHistoRPHI(TH2 * hisA, TH2* hisC){
905 // Make cal pad from r-phi histograms
907 AliTPCROC *proc= AliTPCROC::Instance();
908 AliTPCCalPad *calPad = new AliTPCCalPad("his","his");
909 Float_t globalPos[3];
910 for (Int_t isec=0; isec<72; isec++){
911 AliTPCCalROC* calRoc = calPad->GetCalROC(isec);
912 TH2 * his = ((isec%36<18) ? hisA:hisC);
913 for (UInt_t irow=0; irow<calRoc->GetNrows(); irow+=1){
915 if (isec>=36) jrow+=63;
916 for (UInt_t ipad=0;ipad<proc->GetNPads(isec,irow);ipad+=1){
917 proc->GetPositionGlobal(isec,irow,ipad, globalPos);
918 Double_t phi=TMath::ATan2(globalPos[1],globalPos[0]);
919 //if (phi<0) phi+=TMath::Pi()*2;
920 Int_t bin=his->FindBin(phi,jrow);
921 Float_t value= his->GetBinContent(bin);
922 calRoc->SetValue(irow,ipad,value);
929 AliTPCCalPad *AliTPCCalPad::MakePadFromTree(TTree * treePad, const char *query, const char* name, Bool_t doFast){
931 // make cal pad from the tree
934 ::Error("AliTPCCalPad::MakePadFromTree(TTree * treePad, const char *query, const char* name)","Input tree is missing");
937 if (treePad->GetEntries()!=kNsec) return 0;
938 AliTPCCalPad * calPad= new AliTPCCalPad(name,name);
939 if (name) calPad->SetName(name);
941 for (Int_t iSec=0; iSec<72; iSec++){
942 AliTPCCalROC* calROC = calPad->GetCalROC(iSec);
943 UInt_t nchannels = (UInt_t)treePad->Draw(query,"1","goff",1,iSec);
944 if (nchannels!=calROC->GetNchannels()) {
945 ::Error("AliTPCCalPad::MakePad",TString::Format("%s\t:Wrong query sector\t%d\t%d",treePad->GetName(),iSec,nchannels).Data());
948 for (UInt_t index=0; index<nchannels; index++) calROC->SetValue(index,treePad->GetV1()[index]);
951 UInt_t nchannelsTree = (UInt_t)treePad->Draw(query,"1","goff");
952 UInt_t nchannelsAll=0;
953 for (Int_t iSec=0; iSec<72; iSec++){
954 AliTPCCalROC* calROC = calPad->GetCalROC(iSec);
955 UInt_t nchannels=calROC->GetNchannels();
956 for (UInt_t index=0; index<nchannels; index++) {
957 if (nchannelsAll<=nchannelsTree)calROC->SetValue(index,treePad->GetV1()[nchannelsAll]);
961 if (nchannelsAll>nchannelsTree){
962 ::Error("AliTPCCalPad::MakePad",TString::Format("%s\t:Wrong query: cout mismatch\t%d\t%d",query, nchannelsAll,nchannelsTree).Data());
968 void AliTPCCalPad::AddFriend(TTree * treePad, const char *friendName, const char *fname){
972 TObjArray *fArray = new TObjArray(1);
973 fArray->AddLast(this);
974 this->SetName(friendName);
975 AliTPCCalibViewer::MakeTree(fname, fArray,0);
976 TFile * f = TFile::Open(fname);
977 TTree * tree = (TTree*)f->Get("calPads");
978 treePad->AddFriend(tree,friendName);
979 // tree->AddFriend(TString::Format("%s = calPads",friendName).Data(),fname);