]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPCCalPad.cxx
Changes for #89427: Porting modification for material budget issues to the Release...
[u/mrichter/AliRoot.git] / TPC / AliTPCCalPad.cxx
index 10f33d93a8037793d166e5e073e04824eae099b1..e91e4e8d36c192138e63fbaac2ec556f4aae0d27 100644 (file)
 #include "TTreeStream.h"
 #include "TFile.h"
 #include "TKey.h"
+#include <TFormula.h>
+#include <TString.h>
+#include <TObjString.h>
+#include <iostream>
+#include <AliLog.h>
+
+//graphic includes
+#include <TTree.h>
+#include <TH1.h>
+#include <TCanvas.h>
+#include <TLegend.h>
+#include <TCut.h>
+#include <TVirtualPad.h>
+
 
 ClassImp(AliTPCCalPad)
 
@@ -140,6 +154,7 @@ void AliTPCCalPad::SetCalROC(AliTPCCalROC* roc, Int_t sector){
    // else sector specified in 'roc' is ignored and specified sector is filled
    //
    if (sector == -1) sector = roc->GetSector();
+   if (!fROC[sector]) fROC[sector] = new AliTPCCalROC(sector);
    for (UInt_t ichannel = 0; ichannel < roc->GetNchannels(); ichannel++) 
       fROC[sector]->SetValue(ichannel, roc->GetValue(ichannel));
 }
@@ -181,7 +196,7 @@ void AliTPCCalPad::Add(const AliTPCCalPad * pad, Double_t c1)
   //  - pad by pad -
   //
     for (Int_t isec = 0; isec < kNsec; isec++) {
-       if (fROC[isec]){
+       if (fROC[isec] && pad->GetCalROC(isec)){
            fROC[isec]->Add(pad->GetCalROC(isec),c1);
        }
     }
@@ -365,7 +380,7 @@ Double_t AliTPCCalPad::GetLTM(Double_t *sigma, Double_t fraction, AliTPCCalPad*
 }
 
 //_____________________________________________________________________________
-TH1F * AliTPCCalPad::MakeHisto1D(Float_t min, Float_t max,Int_t type){
+TH1F * AliTPCCalPad::MakeHisto1D(Float_t min, Float_t max,Int_t type, Int_t side){
   //
   // make 1D histo
   // type -1 = user defined range
@@ -398,10 +413,11 @@ TH1F * AliTPCCalPad::MakeHisto1D(Float_t min, Float_t max,Int_t type){
       max = mean+sigma;
     }
   }
-  char  name[1000];
-  sprintf(name,"%s Pad 1D",GetTitle());
-  TH1F * his = new TH1F(name,name,100, min,max);
+  TString name=Form("%s Pad 1D",GetTitle());
+  TH1F * his = new TH1F(name.Data(),name.Data(),100, min,max);
     for (Int_t isec = 0; isec < kNsec; isec++) {
+      if (side==1 && isec%36>18) continue;
+      if (side==-1 && isec%36<18) continue;
        if (fROC[isec]){
            for (UInt_t irow=0; irow<fROC[isec]->GetNrows(); irow++){
                UInt_t npads = (Int_t)fROC[isec]->GetNPads(irow);
@@ -447,7 +463,7 @@ TH2F *AliTPCCalPad::MakeHisto2D(Int_t side){
 }
 
 
-AliTPCCalPad* AliTPCCalPad::LocalFit(const char* padName, Int_t rowRadius, Int_t padRadius, AliTPCCalPad* PadOutliers, Bool_t robust, Double_t chi2Threshold, Double_t robustFraction){
+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 {
    //
    // Loops over all AliTPCCalROCs and performs a localFit in each ROC
    // AliTPCCalPad with fit-data is returned
@@ -465,16 +481,17 @@ AliTPCCalPad* AliTPCCalPad::LocalFit(const char* padName, Int_t rowRadius, Int_t
    //
    AliTPCCalPad* pad = new AliTPCCalPad(padName, padName);
    for (Int_t isec = 0; isec < 72; isec++){
+      if (printCurrentSector) std::cout << "LocalFit in sector " << isec << "\r" << std::flush;
       if (PadOutliers)
-         SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, PadOutliers->GetCalROC(isec), robust, chi2Threshold, robustFraction));
+         pad->SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, PadOutliers->GetCalROC(isec), robust, chi2Threshold, robustFraction));
       else 
-         SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, 0, robust, chi2Threshold, robustFraction));
+         pad->SetCalROC(GetCalROC(isec)->LocalFit(rowRadius, padRadius, 0, robust, chi2Threshold, robustFraction));
    }
    return pad;
 }
 
 
-AliTPCCalPad* AliTPCCalPad::GlobalFit(const char* padName, AliTPCCalPad* PadOutliers, Bool_t robust, Int_t fitType, Double_t chi2Threshold, Double_t robustFraction){
+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){
    //
    // Loops over all AliTPCCalROCs and performs a globalFit in each ROC
    // AliTPCCalPad with fit-data is returned
@@ -482,6 +499,8 @@ AliTPCCalPad* AliTPCCalPad::GlobalFit(const char* padName, AliTPCCalPad* PadOutl
    // robustFraction: Fraction of data that will be used in EvalRobust
    // chi2Threshold: Threshold for chi2 when EvalRobust is called
    // robustFraction: Fraction of data that will be used in EvalRobust
+   // err: error of the data points
+   // if fitParArr and/or fitCovArr is given, write fitParameters and/or covariance Matrices into the array
    //
    AliTPCCalPad* pad = new AliTPCCalPad(padName, padName);
    TVectorD fitParam(0);
@@ -489,360 +508,379 @@ AliTPCCalPad* AliTPCCalPad::GlobalFit(const char* padName, AliTPCCalPad* PadOutl
    Float_t chi2 = 0;
    for (Int_t isec = 0; isec < 72; isec++){
       if (PadOutliers)
-         GetCalROC(isec)->GlobalFit(PadOutliers->GetCalROC(isec), robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction);
+         GetCalROC(isec)->GlobalFit(PadOutliers->GetCalROC(isec), robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction, err);
       else 
-         GetCalROC(isec)->GlobalFit(0, robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction);
-      pad->SetCalROC(AliTPCCalROC::CreateGlobalFitCalROC(fitParam, isec));
+         GetCalROC(isec)->GlobalFit(0, robust, fitParam, covMatrix, chi2, fitType, chi2Threshold, robustFraction, err);
+
+      AliTPCCalROC *roc=AliTPCCalROC::CreateGlobalFitCalROC(fitParam, isec);
+      pad->SetCalROC(roc);
+      delete roc;
+      if ( fitParArr ) fitParArr->AddAtAndExpand(new TVectorD(fitParam), isec);
+      if ( fitCovArr ) fitCovArr->AddAtAndExpand(new TMatrixD(covMatrix), isec);
    }
    return pad;
 }
-
-
-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){
+//_____________________________________________________________________________
+TObjArray* AliTPCCalPad::CreateFormulaArray(const char *fitFormula)
+{
   //
-  // Makes a  GlobalFit over each side and return fit-parameters, covariance and chi2 for each side
-  // fitType == 0: fit plane function
-  // fitType == 1: fit parabolic function
-  // PadOutliers - pads with value !=0 are not used in fitting procedure
-  // chi2Threshold: Threshold for chi2 when EvalRobust is called
-  // robustFraction: Fraction of data that will be used in EvalRobust
+  // create an array of TFormulas for the each parameter of the fit function
   //
-  TLinearFitter* fitterGA = 0;
-  TLinearFitter* fitterGC = 0;
-  
-  if (fitType  == 1) {
-    fitterGA = new TLinearFitter (6,"x0++x1++x2++x3++x4++x5");
-    fitterGC = new TLinearFitter (6,"x0++x1++x2++x3++x4++x5");
-  }
-  else {
-    fitterGA = new TLinearFitter(3,"x0++x1++x2");
-    fitterGC = new TLinearFitter(3,"x0++x1++x2");
+
+  // split fit string in single parameters
+  // find dimension of the fit:
+  TString fitString(fitFormula);
+  fitString.ReplaceAll("++","#");
+  fitString.ReplaceAll(" ","");
+  TObjArray *arrFitParams = fitString.Tokenize("#");
+  Int_t ndim = arrFitParams->GetEntries();
+  //create array of TFormulas to evaluate the parameters
+  TObjArray *arrFitFormulas = new TObjArray(ndim);
+  arrFitFormulas->SetOwner(kTRUE);
+  for (Int_t idim=0;idim<ndim;++idim){
+    TString s=((TObjString*)arrFitParams->At(idim))->GetString();
+    s.ReplaceAll("gx","[0]");
+    s.ReplaceAll("gy","[1]");
+    s.ReplaceAll("lx","[2]");
+    s.ReplaceAll("ly","[3]");
+    s.ReplaceAll("sector","[4]");
+    arrFitFormulas->AddAt(new TFormula(Form("param%02d",idim),s.Data()),idim);
   }
-  fitterGA->StoreData(kTRUE);   
-  fitterGC->StoreData(kTRUE);   
-  fitterGA->ClearPoints();
-  fitterGC->ClearPoints();
-  Double_t xx[6];  
-  Int_t    npointsA=0;
-  Int_t    npointsC=0;
+  delete arrFitParams;
   
-  Float_t localXY[3] = {0}; // pad's position, needed to get the pad's position
-  Float_t lx, ly;  // pads position
+  return arrFitFormulas;
+}
+//_____________________________________________________________________________
+void AliTPCCalPad::EvalFormulaArray(const TObjArray &arrFitFormulas, TVectorD &results,
+                                    const Int_t sec, const Int_t row, const Int_t pad)
+{
+  //
+  // evaluate the fit formulas
+  //
+  Int_t ndim=arrFitFormulas.GetEntries();
+  results.ResizeTo(ndim);
   
   AliTPCROC* tpcROCinstance = AliTPCROC::Instance();  // to calculate the pad's position
+  Float_t localXYZ[3];
+  Float_t globalXYZ[3];
+  tpcROCinstance->GetPositionLocal(sec, row, pad, localXYZ);
+  tpcROCinstance->GetPositionGlobal(sec, row, pad, globalXYZ);
+  //calculate parameter values
+  for (Int_t idim=0;idim<ndim;++idim){
+    TFormula *f=(TFormula*)arrFitFormulas.At(idim);
+    f->SetParameters(globalXYZ[0],globalXYZ[1],localXYZ[0],localXYZ[1],sec);
+    results[idim]=f->Eval(0);
+  }
+}
+//_____________________________________________________________________________
+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){
+  //
+  // Performs a fit on both sides.
+  // Valid information for the fitFormula are the variables
+  // - gx, gy, lx ,ly: meaning global x, global y, local x, local y value of the padName
+  // - sector:         the sector number.
+  //  eg. a formula might look 'gy' or '(sector<36) ++ gy' or 'gx ++ gy' or 'gx ++ gy ++ lx ++ lx^2' and so on
+  //
+  // PadOutliers - pads with value !=0 are not used in fitting procedure
+  // chi2Threshold: Threshold for chi2 when EvalRobust is called
+  // robustFraction: Fraction of data that will be used in EvalRobust
+  //
+
+  TObjArray* arrFitFormulas=CreateFormulaArray(fitFormula);
+  Int_t ndim = arrFitFormulas->GetEntries();
+  //resize output data arrays
+  fitParamSideA.ResizeTo(ndim+1);
+  fitParamSideC.ResizeTo(ndim+1);
+  covMatrixSideA.ResizeTo(ndim+1,ndim+1);
+  covMatrixSideC.ResizeTo(ndim+1,ndim+1);
+  // create linear fitter for A- and C- Side
+  TLinearFitter* fitterGA = new TLinearFitter(ndim+1,Form("hyp%d",ndim));
+  TLinearFitter* fitterGC = new TLinearFitter(ndim+1,Form("hyp%d",ndim));
+  fitterGA->StoreData(kTRUE);
+  fitterGC->StoreData(kTRUE);
+  //parameter values
+  TVectorD parValues(ndim);
+
+  AliTPCCalROC *rocErr=0x0;
   
-  // loop over all sectors and pads and read data into fitterGA and fitterGC 
-  if (fitType == 1) {  
-  // parabolic fit
-    fitParamSideA.ResizeTo(6);
-    fitParamSideC.ResizeTo(6);
-    covMatrixSideA.ResizeTo(6,6);
-    covMatrixSideC.ResizeTo(6,6);
-    for (UInt_t isec = 0; isec<72; isec++){
-      for (UInt_t irow = 0; irow < GetCalROC(isec)->GetNrows(); irow++) {
-         for (UInt_t ipad = 0; ipad < GetCalROC(isec)->GetNPads(irow); ipad++) {
-            // fill fitterG
-            tpcROCinstance->GetPositionLocal(isec, irow, ipad, localXY);   // calculate position localXY by sector, pad and row number
-            lx = localXY[0];
-            ly = localXY[1];
-            xx[0] = 1;
-            xx[1] = lx;
-            xx[2] = ly;
-            xx[3] = lx*lx;
-            xx[4] = ly*ly;
-            xx[5] = lx*ly;
-            if (!PadOutliers || PadOutliers->GetCalROC(isec)->GetValue(irow, ipad) != 1) {
-            // if given pad is no outlier, add it to TLinearFitter, decide to which of both
-//                sector  0 - 17: IROC, A
-//                sector 18 - 35: IROC, C
-//                sector 36 - 53: OROC, A
-//                sector 54 - 71: CROC, C
-               if (isec <= 17 || (isec >= 36 && isec <= 53)) { // Side A
-                  npointsA++;
-                  fitterGA->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);  
-               }
-               else { // side C
-                  npointsC++;
-                  fitterGC->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);  
-               }
-            }
-         }
+  for (UInt_t isec = 0; isec<kNsec; ++isec){
+    AliTPCCalROC *rocOut=PadOutliers->GetCalROC(isec);
+    AliTPCCalROC *rocData=GetCalROC(isec);
+    if (pointError) rocErr=pointError->GetCalROC(isec);
+    if (!rocData) continue;
+    for (UInt_t irow = 0; irow < GetCalROC(isec)->GetNrows(); irow++) {
+      for (UInt_t ipad = 0; ipad < GetCalROC(isec)->GetNPads(irow); ipad++) {
+        //check for outliers
+        if (rocOut && rocOut->GetValue(irow,ipad)) continue;
+        //calculate parameter values
+        EvalFormulaArray(*arrFitFormulas,parValues,isec,irow,ipad);
+        //get value
+        Float_t value=rocData->GetValue(irow,ipad);
+        //point error
+        Int_t err=1;
+        if (rocErr) {
+          err=TMath::Nint(rocErr->GetValue(irow,ipad));
+          if (err==0) err=1;
+        }
+        //add points to the fitters
+        if (isec/18%2==0){
+          fitterGA->AddPoint(parValues.GetMatrixArray(),value,err);
+        }else{
+          fitterGC->AddPoint(parValues.GetMatrixArray(),value,err);
+        }
       }
     }
   }
-  else {   
-  // linear fit
-    fitParamSideA.ResizeTo(3);
-    fitParamSideC.ResizeTo(3);
-    covMatrixSideA.ResizeTo(3,3);
-    covMatrixSideC.ResizeTo(3,3);
-    
-    for (UInt_t isec = 0; isec<72; isec++){
-      for (UInt_t irow = 0; irow < GetCalROC(isec)->GetNrows(); irow++) {
-         for (UInt_t ipad = 0; ipad < GetCalROC(isec)->GetNPads(irow); ipad++) {
-            // fill fitterG
-            tpcROCinstance->GetPositionLocal(isec, irow, ipad, localXY);   // calculate position localXY by sector, pad and row number
-            lx = localXY[0];
-            ly = localXY[1];
-            xx[0] = 1;
-            xx[1] = lx;
-            xx[2] = ly;
-            if (!PadOutliers || PadOutliers->GetCalROC(isec)->GetValue(irow, ipad) != 1) {
-            // if given pad is no outlier, add it to TLinearFitter, decide to which of both
-//                sector  0 - 17: IROC, A
-//                sector 18 - 35: IROC, C
-//                sector 36 - 53: OROC, A
-//                sector 54 - 71: CROC, C
-               if (isec <= 17 || (isec >= 36 && isec <= 53)) { 
-               // Side A
-                  npointsA++;
-                  fitterGA->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);  
-               }
-               else { 
-               // side C
-                  npointsC++;
-                  fitterGC->AddPoint(xx, GetCalROC(isec)->GetValue(irow, ipad), 1);  
-               }
-            }
-         }
-      }
-    }
-  }    
-  
-  fitterGA->Eval();
-  fitterGC->Eval();
+  if (robust){
+    fitterGA->EvalRobust(robustFraction);
+    fitterGC->EvalRobust(robustFraction);
+  } else {
+    fitterGA->Eval();
+    fitterGC->Eval();
+  }
+  chi2SideA=fitterGA->GetChisquare()/(fitterGA->GetNpoints()-(ndim+1));
+  chi2SideC=fitterGC->GetChisquare()/(fitterGC->GetNpoints()-(ndim+1));
   fitterGA->GetParameters(fitParamSideA);
   fitterGC->GetParameters(fitParamSideC);
   fitterGA->GetCovarianceMatrix(covMatrixSideA);
   fitterGC->GetCovarianceMatrix(covMatrixSideC);
-  if (fitType == 1){
-    chi2SideA = fitterGA->GetChisquare()/(npointsA-6.);
-    chi2SideC = fitterGC->GetChisquare()/(npointsC-6.);
+  
+  delete arrFitFormulas;
+  delete fitterGA;
+  delete fitterGC;
+  
+}
+//
+AliTPCCalPad *AliTPCCalPad::CreateCalPadFit(const char* fitFormula, const TVectorD &fitParamSideA, const TVectorD &fitParamSideC)
+{
+  //
+  //
+  //
+  TObjArray *arrFitFormulas=CreateFormulaArray(fitFormula);
+  Int_t ndim = arrFitFormulas->GetEntries();
+  //check if dimension of fit formula and fit parameters agree
+  if (ndim!=fitParamSideA.GetNrows()||ndim!=fitParamSideC.GetNrows()){
+    printf("AliTPCCalPad::CreateCalPadFit: Dimensions of fit formula and fit Parameters does not match!");
+    return 0;
   }
-  else {
-   chi2SideA = fitterGA->GetChisquare()/(npointsA-3.);
-   chi2SideC = fitterGC->GetChisquare()/(npointsC-3.);
+  //create cal pad
+  AliTPCCalPad *pad=new AliTPCCalPad("fitResultPad",Form("Fit result: %s",fitFormula));
+  //fill cal pad with fit results if requested
+  for (UInt_t isec = 0; isec<kNsec; ++isec){
+    AliTPCCalROC *roc=pad->GetCalROC(isec);
+    for (UInt_t irow = 0; irow < roc->GetNrows(); irow++) {
+      for (UInt_t ipad = 0; ipad < roc->GetNPads(irow); ipad++) {
+        const TVectorD *fitPar=0;
+        TVectorD fitResArray;
+        if (isec/18%2==0){
+          fitPar=&fitParamSideA;
+        }else{
+          fitPar=&fitParamSideC;
+        }
+        EvalFormulaArray(*arrFitFormulas,fitResArray, isec, irow, ipad);
+        for (Int_t idim=0;idim<ndim;++idim)
+          fitResArray(idim)*=(*fitPar)(idim);
+        roc->SetValue(irow,ipad,fitResArray.Sum());
+      }
+    }
   }
-  if (robust && chi2SideA > chi2Threshold) {
-    //    std::cout << "robust fitter called... " << std::endl;
-    fitterGA->EvalRobust(robustFraction);
-    fitterGA->GetParameters(fitParamSideA);
+  delete arrFitFormulas;
+  return pad;
+}
+
+
+
+TCanvas * AliTPCCalPad::MakeReportPadSector(TTree *chain, const char* varName, const char*varTitle, const char *axisTitle, Float_t min, Float_t max, const char *cutUser){
+  //
+  // Make a report - cal pads per sector
+  // mean valeus per sector and local X
+  //
+  TH1* his=0; 
+  TLegend *legend = 0;
+  TCanvas *canvas = new TCanvas(Form("Sector: %s",varTitle),Form("Sector: %s",varTitle),1500,1100);
+
+  canvas->Divide(2);
+  chain->SetAlias("lX","lx.fElements"); 
+  //
+  canvas->cd(1);
+  TString strDraw=varName;
+  strDraw+=":lX";
+  legend = new TLegend(0.5,0.50,0.9,0.9, Form("%s TPC A side", varTitle));
+  for (Int_t isec=-1; isec<18; isec+=1){
+    TCut cutSec=Form("sector%%36==%d",isec);
+    cutSec+=cutUser;
+    if (isec==-1) cutSec="sector%36<18";
+    chain->SetMarkerColor(1+(isec+2)%5);
+    chain->SetLineColor(1+(isec+2)%5);
+    chain->SetMarkerStyle(25+(isec+2)%4);
+    //
+    chain->Draw(strDraw.Data(),cutSec,"profgoff");
+    his=(TH1*)chain->GetHistogram()->Clone();
+    delete chain->GetHistogram();
+    his->SetMaximum(max);
+    his->SetMinimum(min);
+    his->GetXaxis()->SetTitle("R (cm)");
+    his->GetYaxis()->SetTitle(axisTitle);
+    his->SetTitle(Form("%s- sector %d",varTitle, isec));
+    his->SetName(Form("%s- sector %d",varTitle, isec));
+    if (isec==-1) his->SetTitle(Form("%s A side",varTitle));
+    if (isec==-1) his->Draw();
+    his->Draw("same");
+    legend->AddEntry(his);
   }
-  if (robust && chi2SideC > chi2Threshold) {
-    //    std::cout << "robust fitter called... " << std::endl;
-    fitterGC->EvalRobust(robustFraction);
-    fitterGC->GetParameters(fitParamSideC);
+  legend->Draw();
+  canvas->cd(2);
+  //
+  legend = new TLegend(0.5,0.50,0.9,0.9, Form("%s TPC C side", varTitle));
+  for (Int_t isec=-1; isec<18; isec+=1){
+    TCut cutSec=Form("(sector+18)%%36==%d",isec);
+    cutSec+=cutUser;
+    if (isec==-1) cutSec="sector%36>18";
+    chain->SetMarkerColor(1+(isec+2)%5);
+    chain->SetLineColor(1+(isec+2)%5);
+    chain->SetMarkerStyle(25+isec%4);
+    //
+    chain->Draw(strDraw.Data(),cutSec,"profgoff");
+    his=(TH1*)chain->GetHistogram()->Clone();
+    delete chain->GetHistogram();
+    his->SetMaximum(max);
+    his->SetMinimum(min);
+    his->GetXaxis()->SetTitle("R (cm)");
+    his->GetYaxis()->SetTitle(axisTitle);
+    his->SetTitle(Form("%s- sector %d",varTitle,isec));
+    his->SetName(Form("%s- sector %d",varTitle,isec));
+    if (isec==-1) his->SetTitle(Form("%s C side",varTitle));
+    if (isec==-1) his->Draw();
+    his->Draw("same");
+    legend->AddEntry(his);
   }
-  delete fitterGA;
-  delete fitterGC;
+  legend->Draw();
+  //
+  //
+  return canvas;
 }
 
 
-void AliTPCCalPad::MakeTree(const char * fileName, TObjArray * array, const char * mapFileName, AliTPCCalPad* outlierPad, Float_t ltmFraction) {
-  //
-  // Write a tree with all available information
-  // im mapFileName is speciefied, the Map information are also written to the tree
-  // pads specified in outlierPad are not used for calculating statistics
-  //  - this function will be moved to AliTPCCalibViewer - 
-  //  - DO NOT USE THIS FUNCTION ANY MORE - 
-  //
-   AliTPCROC* tpcROCinstance = AliTPCROC::Instance();
-
-   TObjArray* mapIROCs = 0;
-   TObjArray* mapOROCs = 0;
-   TVectorF *mapIROCArray = 0;
-   TVectorF *mapOROCArray = 0;
-   Int_t mapEntries = 0;
-   TString* mapNames = 0;
-   
-   if (mapFileName) {
-      TFile mapFile(mapFileName, "read");
-      
-      TList* listOfROCs = mapFile.GetListOfKeys();
-      mapEntries = listOfROCs->GetEntries()/2;
-      mapIROCs = new TObjArray(mapEntries*2);
-      mapOROCs = new TObjArray(mapEntries*2);
-      mapIROCArray = new TVectorF[mapEntries];
-      mapOROCArray = new TVectorF[mapEntries];
-      
-      mapNames = new TString[mapEntries];
-      for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
-         TString ROCname(((TKey*)(listOfROCs->At(ivalue*2)))->GetName());
-         ROCname.Remove(ROCname.Length()-4, 4);
-         mapIROCs->AddAt((AliTPCCalROC*)mapFile.Get((ROCname + "IROC").Data()), ivalue);
-         mapOROCs->AddAt((AliTPCCalROC*)mapFile.Get((ROCname + "OROC").Data()), ivalue);
-         mapNames[ivalue].Append(ROCname);
-      }
-      
-      for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
-         mapIROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(0));
-         mapOROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(36));
-      
-         for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(0); ichannel++)
-            (mapIROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapIROCs->At(ivalue)))->GetValue(ichannel);
-         for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(36); ichannel++)
-            (mapOROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapOROCs->At(ivalue)))->GetValue(ichannel);
-      }
-
-   } //  if (mapFileName)
-  
-   TTreeSRedirector cstream(fileName);
-   Int_t arrayEntries = array->GetEntries();
-   
-   TString* names = new TString[arrayEntries];
-   for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
-      names[ivalue].Append(((AliTPCCalPad*)array->At(ivalue))->GetName());
-
-   for (UInt_t isector = 0; isector < tpcROCinstance->GetNSectors(); isector++) {
-      //
-      // get statistic for given sector
-      //
-      TVectorF median(arrayEntries);
-      TVectorF mean(arrayEntries);
-      TVectorF rms(arrayEntries);
-      TVectorF ltm(arrayEntries);
-      TVectorF ltmrms(arrayEntries);
-      TVectorF medianWithOut(arrayEntries);
-      TVectorF meanWithOut(arrayEntries);
-      TVectorF rmsWithOut(arrayEntries);
-      TVectorF ltmWithOut(arrayEntries);
-      TVectorF ltmrmsWithOut(arrayEntries);
-      
-      TVectorF *vectorArray = new TVectorF[arrayEntries];
-      for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
-         vectorArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
-      
-      for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
-         AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
-         AliTPCCalROC* calROC = calPad->GetCalROC(isector);
-         AliTPCCalROC* outlierROC = 0;
-         if (outlierPad) outlierROC = outlierPad->GetCalROC(isector);
-         if (calROC) {
-            median[ivalue] = calROC->GetMedian();
-            mean[ivalue] = calROC->GetMean();
-            rms[ivalue] = calROC->GetRMS();
-            Double_t ltmrmsValue = 0;
-            ltm[ivalue] = calROC->GetLTM(&ltmrmsValue, ltmFraction);
-            ltmrms[ivalue] = ltmrmsValue;
-            if (outlierROC) {
-               medianWithOut[ivalue] = calROC->GetMedian(outlierROC);
-               meanWithOut[ivalue] = calROC->GetMean(outlierROC);
-               rmsWithOut[ivalue] = calROC->GetRMS(outlierROC);
-               ltmrmsValue = 0;
-               ltmWithOut[ivalue] = calROC->GetLTM(&ltmrmsValue, ltmFraction, outlierROC);
-               ltmrmsWithOut[ivalue] = ltmrmsValue;
-            }
-         }
-         else {
-            median[ivalue] = 0.;
-            mean[ivalue] = 0.;
-            rms[ivalue] = 0.;
-            ltm[ivalue] = 0.;
-            ltmrms[ivalue] = 0.;
-            medianWithOut[ivalue] = 0.;
-            meanWithOut[ivalue] = 0.;
-            rmsWithOut[ivalue] = 0.;
-            ltmWithOut[ivalue] = 0.;
-            ltmrmsWithOut[ivalue] = 0.;
-         }
-      }
-      
-      //
-      // fill vectors of variable per pad
-      //
-      TVectorF *posArray = new TVectorF[8];
-      for (Int_t ivalue = 0; ivalue < 8; ivalue++)
-         posArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
-
-      Float_t posG[3] = {0};
-      Float_t posL[3] = {0};
-      Int_t ichannel = 0;
-      for (UInt_t irow = 0; irow < tpcROCinstance->GetNRows(isector); irow++) {
-         for (UInt_t ipad = 0; ipad < tpcROCinstance->GetNPads(isector, irow); ipad++) {
-            tpcROCinstance->GetPositionLocal(isector, irow, ipad, posL);
-            tpcROCinstance->GetPositionGlobal(isector, irow, ipad, posG);
-            posArray[0][ichannel] = irow;
-            posArray[1][ichannel] = ipad;
-            posArray[2][ichannel] = posL[0];
-            posArray[3][ichannel] = posL[1];
-            posArray[4][ichannel] = posG[0];
-            posArray[5][ichannel] = posG[1];
-            posArray[6][ichannel] = (Int_t)(ipad - (Double_t)(tpcROCinstance->GetNPads(isector, irow))/2);
-            posArray[7][ichannel] = ichannel;
-            
-            // loop over array containing AliTPCCalPads
-            for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
-               AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
-               AliTPCCalROC* calROC = calPad->GetCalROC(isector);
-               if (calROC)
-                  (vectorArray[ivalue])[ichannel] = calROC->GetValue(irow, ipad);
-               else
-                  (vectorArray[ivalue])[ichannel] = 0;
-            }
-            ichannel++;
-         }
-      }
-      
-      cstream << "calPads" <<
-         "sector=" << isector;
-      
-      for  (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
-         cstream << "calPads" <<
-            (Char_t*)((names[ivalue] + "_Median=").Data()) << median[ivalue] <<
-            (Char_t*)((names[ivalue] + "_Mean=").Data()) << mean[ivalue] <<
-            (Char_t*)((names[ivalue] + "_RMS=").Data()) << rms[ivalue] <<
-            (Char_t*)((names[ivalue] + "_LTM=").Data()) << ltm[ivalue] <<
-            (Char_t*)((names[ivalue] + "_RMS_LTM=").Data()) << ltmrms[ivalue];
-         if (outlierPad) {
-            cstream << "calPads" <<
-               (Char_t*)((names[ivalue] + "_Median_OutlierCutted=").Data()) << medianWithOut[ivalue] <<
-               (Char_t*)((names[ivalue] + "_Mean_OutlierCutted=").Data()) << meanWithOut[ivalue] <<
-               (Char_t*)((names[ivalue] + "_RMS_OutlierCutted=").Data()) << rmsWithOut[ivalue] <<
-               (Char_t*)((names[ivalue] + "_LTM_OutlierCutted=").Data()) << ltmWithOut[ivalue] <<
-               (Char_t*)((names[ivalue] + "_RMS_LTM_OutlierCutted=").Data()) << ltmrmsWithOut[ivalue];
-         }
-      }
+TCanvas * AliTPCCalPad::MakeReportPadSector2D(TTree *chain, const char* varName, const char*varTitle, const char *axisTitle, Float_t min, Float_t max, const char *cutUser){
+  //
+  // Make a report - cal pads per sector
+  // 2D view
+  // Input tree should be created using AliPreprocesorOnline before
+  // 
+  TH1* his=0; 
+  TCanvas *canvas = new TCanvas(Form("%s2D",varTitle),Form("%s2D",varTitle),1500,1100);
+  canvas->Divide(2);
+  //
+  TString strDraw=varName;
+  strDraw+=":gy.fElements:gx.fElements>>his(250,-250,250,250,-250,250)";
+  //
+  TVirtualPad * pad=0;
+  pad=canvas->cd(1);
+  pad->SetMargin(0.15,0.15,0.15,0.15);
+  TCut cut=cutUser;
+  chain->Draw(strDraw.Data(),"sector%36<18"+cut,"profgoffcolz2");
+  his=(TH1*)chain->GetHistogram()->Clone();
+  delete chain->GetHistogram();
+  his->SetMaximum(max);
+  his->SetMinimum(min);
+  his->GetXaxis()->SetTitle("x (cm)");
+  his->GetYaxis()->SetTitle("y (cm)");
+  his->GetZaxis()->SetTitle(axisTitle);
+  his->SetTitle(Form("%s A side",varTitle));
+  his->SetName(Form("%s A side",varTitle));
+  his->Draw("colz2");
+  //
+  pad=canvas->cd(2);
+  pad->SetMargin(0.15,0.15,0.15,0.15);
+
+  chain->Draw(strDraw.Data(),"sector%36>=18"+cut,"profgoffcolz2");
+  his=(TH1*)chain->GetHistogram()->Clone();
+  delete chain->GetHistogram();
+  his->SetMaximum(max);
+  his->SetMinimum(min);
+  his->GetXaxis()->SetTitle("x (cm)");
+  his->GetYaxis()->SetTitle("y (cm)");
+  his->GetZaxis()->SetTitle(axisTitle);
+  his->SetTitle(Form("%s C side",varTitle));
+  his->SetName(Form("%s C side",varTitle));
+  his->Draw("colz2");
+  //
+  //
+  return canvas;
+}
 
-      for  (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
-         cstream << "calPads" <<
-            (Char_t*)((names[ivalue] + ".=").Data()) << &vectorArray[ivalue];
-      }
+void  AliTPCCalPad::Draw(Option_t* option){
+  // 
+  // Draw function - standard 2D view
+  //
+  TH1* his=0; 
+  TCanvas *canvas = new TCanvas(Form("%s2D",GetTitle()),Form("%s2D",GetTitle()),900,900);
+  canvas->Divide(2,2);
+  //
+  //
+  TVirtualPad * pad=0;
+  pad=canvas->cd(1);
+  pad->SetMargin(0.15,0.15,0.15,0.15);
+  his=MakeHisto2D(0);
+  his->GetXaxis()->SetTitle("x (cm)");
+  his->GetYaxis()->SetTitle("y (cm)");
+  his->GetZaxis()->SetTitle(GetTitle());
+  his->SetTitle(Form("%s A side",GetTitle()));
+  his->SetName(Form("%s A side",GetTitle()));
+  his->Draw(option);
+  //
+  pad=canvas->cd(2);
+  pad->SetMargin(0.15,0.15,0.15,0.15);
+  his=MakeHisto2D(1);
+  his->GetXaxis()->SetTitle("x (cm)");
+  his->GetYaxis()->SetTitle("y (cm)");
+  his->GetZaxis()->SetTitle(GetTitle());
+  his->SetTitle(Form("%s C side",GetTitle()));
+  his->SetName(Form("%s C side",GetTitle()));
+  his->Draw(option);
+  //
+  pad=canvas->cd(3);
+  pad->SetMargin(0.15,0.15,0.15,0.15);
+  his=MakeHisto1D(-8,8,0,1);
+  his->GetXaxis()->SetTitle(GetTitle());
+  his->SetTitle(Form("%s A side",GetTitle()));
+  his->SetName(Form("%s A side",GetTitle()));
+  his->Draw("err");
+  //
+  pad=canvas->cd(4);
+  pad->SetMargin(0.15,0.15,0.15,0.15);
+  his=MakeHisto1D(-8,8,0,-1);
+  his->GetXaxis()->SetTitle(GetTitle());
+  his->SetTitle(Form("%s C side",GetTitle()));
+  his->SetName(Form("%s C side",GetTitle()));
+  his->Draw("err");
 
-      if (mapFileName) {
-         for  (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
-            if (isector < 36)
-               cstream << "calPads" <<
-                  (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapIROCArray[ivalue];
-            else
-               cstream << "calPads" <<
-                  (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapOROCArray[ivalue];
-         }
-      }
 
-      cstream << "calPads" <<
-         "row.=" << &posArray[0] <<
-         "pad.=" << &posArray[1] <<
-         "lx.=" << &posArray[2] <<
-         "ly.=" << &posArray[3] <<
-         "gx.=" << &posArray[4] <<
-         "gy.=" << &posArray[5] <<
-         "rpad.=" << &posArray[6] <<
-         "channel.=" << &posArray[7];
-         
-      cstream << "calPads" <<
-         "\n";
-
-      delete[] posArray;
-      delete[] vectorArray;
-   }
-   delete[] names;
-   if (mapFileName) {
-      delete mapIROCs;
-      delete mapOROCs;
-      delete[] mapIROCArray;
-      delete[] mapOROCArray;
-      delete[] mapNames;
-   }
 }
 
 
+AliTPCCalPad * AliTPCCalPad::MakeCalPadFromHistoRPHI(TH2 * hisA, TH2* hisC){
+  //
+  // Make cal pad from r-phi histograms
+  //
+  AliTPCROC *proc= AliTPCROC::Instance();
+  AliTPCCalPad *calPad = new AliTPCCalPad("his","his");
+  Float_t globalPos[3];
+  for (Int_t isec=0; isec<72; isec++){
+    AliTPCCalROC* calRoc  = calPad->GetCalROC(isec);
+    TH2 * his = ((isec%36<18) ? hisA:hisC);
+    for (UInt_t irow=0; irow<calRoc->GetNrows(); irow+=1){
+      Int_t jrow=irow;
+      if (isec>=36) jrow+=63;
+      for (UInt_t ipad=0;ipad<proc->GetNPads(isec,irow);ipad+=1){
+        proc->GetPositionGlobal(isec,irow,ipad, globalPos);
+        Double_t phi=TMath::ATan2(globalPos[1],globalPos[0]);
+        //if (phi<0) phi+=TMath::Pi()*2;
+        Int_t bin=his->FindBin(phi,jrow);
+        Float_t value= his->GetBinContent(bin);
+       calRoc->SetValue(irow,ipad,value);
+      }
+    }
+  }
+  return calPad;
+}