AliTPCcalibDB.cxx AliTPCcalibDB.h - spline fit creation
authormarian <marian@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 12 Oct 2009 05:16:04 +0000 (05:16 +0000)
committermarian <marian@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 12 Oct 2009 05:16:04 +0000 (05:16 +0000)
AliTPCcalibDButil.cxx AliTPCcalibDButil.h - adding more robust oultyer removal algorithm. + alternative CreatePadTime0
CalibMacros/CalibEnv.C                    - cleaning of the macro
                                          - adding pad gain
                                          - restructuring of the loops

TPC/AliTPCcalibDB.cxx
TPC/AliTPCcalibDB.h
TPC/AliTPCcalibDButil.cxx
TPC/AliTPCcalibDButil.h
TPC/CalibMacros/CalibEnv.C

index 65c0f26..ddf0e7e 100644 (file)
@@ -874,6 +874,21 @@ AliSplineFit* AliTPCcalibDB::GetVdriftSplineFit(const char* name, Int_t run){
   return dynamic_cast<AliSplineFit*>(arr->FindObject(name));
 }
 
+AliSplineFit* AliTPCcalibDB::CreateVdriftSplineFit(const char* graphName, Int_t run){
+  //
+  // create spline fit from the drift time graph in TimeDrift
+  //
+  TObjArray *arr=GetTimeVdriftSplineRun(run);
+  if (!arr) return 0;
+  TGraph *graph=dynamic_cast<TGraph*>(arr->FindObject(graphName));
+  if (!graph) return 0;
+  AliSplineFit *fit = new AliSplineFit();
+  fit->SetGraph(graph);
+  fit->SetMinPoints(graph->GetN()+1);
+  fit->InitKnots(graph,2,0,0.001);
+  fit->SplineFit(0);
+  return fit;
+}
 
 AliGRPObject *AliTPCcalibDB::GetGRP(Int_t run){
   //
index c7a53ea..7a92454 100644 (file)
@@ -128,6 +128,7 @@ class AliTPCcalibDB : public TObject
   static Float_t GetGain(Int_t sector, Int_t row, Int_t pad);
   //
   AliSplineFit* GetVdriftSplineFit(const char* name, Int_t run);
+  AliSplineFit* CreateVdriftSplineFit(const char* graphName, Int_t run);
   //
   static void     CreateObjectList(const Char_t *filename, TObjArray *calibObjects);
   static void MakeTree(const char * fileName, TObjArray * array, const char * mapFileName = 0, AliTPCCalPad* outlierPad = 0, Float_t ltmFraction = 0.9);
index 9e09e24..e85f195 100644 (file)
@@ -42,6 +42,7 @@
 #include "AliTPCCalibRaw.h"
 
 #include "AliTPCcalibDButil.h"
+#include "AliTPCPreprocessorOnline.h"
 
 ClassImp(AliTPCcalibDButil)
 AliTPCcalibDButil::AliTPCcalibDButil() :
@@ -124,7 +125,7 @@ void AliTPCcalibDButil::UpdateFromCalibDB()
 }
 //_____________________________________________________________________________________
 void AliTPCcalibDButil::ProcessCEdata(const char* fitFormula, TVectorD &fitResultsA, TVectorD &fitResultsC,
-                                      Int_t &noutliersCE, AliTPCCalPad *outCE)
+                                      Int_t &noutliersCE, Double_t & chi2A, Double_t &chi2C, AliTPCCalPad *outCE)
 {
   //
   // Process the CE data for this run
@@ -199,8 +200,10 @@ void AliTPCcalibDButil::ProcessCEdata(const char* fitFormula, TVectorD &fitResul
   }
   //perform fit
   TMatrixD dummy;
-  Float_t chi2A,chi2C;
-  fCETmean->GlobalSidesFit(out,fitFormula,fitResultsA,fitResultsC,dummy,dummy,chi2A,chi2C);
+  Float_t chi2Af,chi2Cf;
+  fCETmean->GlobalSidesFit(out,fitFormula,fitResultsA,fitResultsC,dummy,dummy,chi2Af,chi2Cf);
+  chi2A=chi2Af;
+  chi2C=chi2Cf;
   if (!outCE) delete out;
 }
 //_____________________________________________________________________________________
@@ -796,7 +799,7 @@ void AliTPCcalibDButil::PulserOutlierMap(AliTPCCalPad *pulOut, const AliTPCCalPa
   }
 }
 //_____________________________________________________________________________________
-AliTPCCalPad* AliTPCcalibDButil::CreatePadTime0(Int_t model)
+AliTPCCalPad* AliTPCcalibDButil::CreatePadTime0(Int_t model, Double_t &gyA, Double_t &gyC, Double_t &chi2A, Double_t &chi2C )
 {
   //
   // Create pad time0 object from pulser and/or CE data, depending on the selected model
@@ -804,7 +807,10 @@ AliTPCCalPad* AliTPCcalibDButil::CreatePadTime0(Int_t model)
   // Model 1: normalise IROCs/OROCs of each readout side to its mean, only Pulser
   // Model 2: use CE data and a combination CE fit + pulser in the outlier regions.
   //
-  
+  // In case model 2 is invoked - gy arival time gradient is also returned
+  //
+  gyA=0;
+  gyC=0;
   AliTPCCalPad *padTime0=new AliTPCCalPad("PadTime0",Form("PadTime0-Model_%d",model));
   // decide between different models
   if (model==0||model==1){
@@ -838,16 +844,23 @@ AliTPCCalPad* AliTPCcalibDButil::CreatePadTime0(Int_t model)
         }
       }
     }
-  } else if (model==2){
+  } else if (model==2){  
+    Double_t pgya,pgyc,pchi2a,pchi2c;
+    AliTPCCalPad * padPulser = CreatePadTime0(1,pgya,pgyc,pchi2a,pchi2c);
+    fCETmean->Add(padPulser,-1.);
     TVectorD vA,vC;
     AliTPCCalPad outCE("outCE","outCE");
     Int_t nOut;
-    ProcessCEdata("(sector<36)++gy++gx++(lx-134)++(sector<36)*(lx-134)",vA,vC,nOut,&outCE);
-    AliTPCCalPad *padFit=AliTPCCalPad::CreateCalPadFit("1++0++gy++0++(lx-134)++0",vA,vC);
+    ProcessCEdata("(sector<36)++gy++gx++(lx-134)++(sector<36)*(lx-134)++(ly/lx)^2",vA,vC,nOut,chi2A, chi2C,&outCE);
+    AliTPCCalPad *padFit=AliTPCCalPad::CreateCalPadFit("1++0++gy++0++(lx-134)++0++0",vA,vC);
 //     AliTPCCalPad *padFit=AliTPCCalPad::CreateCalPadFit("1++(sector<36)++gy++gx++(lx-134)++(sector<36)*(lx-134)",vA,vC);
-    if (!padFit) return 0;
+    if (!padFit) { delete padPulser; return 0;}
+    gyA=vA[2];
+    gyC=vC[2];
+    fCETmean->Add(padPulser,1.);
     padTime0->Add(fCETmean);
-    padTime0->Add(padFit,-1);
+    padTime0->Add(padFit,-1);  
+    delete padPulser;
     TVectorD vFitROC;
     TMatrixD mFitROC;
     Float_t chi2;
@@ -876,7 +889,8 @@ AliTPCCalPad* AliTPCcalibDButil::CreatePadTime0(Int_t model)
     }
     delete padFit;
   }
-
+  Double_t median = padTime0->GetMedian();
+  padTime0->Add(-median);  // normalize to median
   return padTime0;
 }
 //_____________________________________________________________________________________
@@ -933,3 +947,323 @@ void AliTPCcalibDButil::SetRefFile(const char* filename)
   currDir->cd();
 }
 
+
+
+
+AliTPCCalPad *AliTPCcalibDButil::CreateCEOutlyerMap( Int_t & noutliersCE, AliTPCCalPad *ceOut, Float_t minSignal, Float_t cutTrmsMin,  Float_t cutTrmsMax, Float_t cutMaxDistT){
+  //
+  // Author:  marian.ivanov@cern.ch
+  //
+  // Create outlier map for CE study
+  // Parameters:
+  //  Return value - outlyer map
+  //  noutlyersCE  - number of outlyers
+  //  minSignal    - minimal total Q signal
+  //  cutRMSMin    - minimal width of the signal in respect to the median 
+  //  cutRMSMax    - maximal width of the signal in respect to the median 
+  //  cutMaxDistT  - maximal deviation from time median per chamber
+  //
+  // Outlyers criteria:
+  // 0. Exclude masked pads
+  // 1. Exclude first two rows in IROC and last two rows in OROC
+  // 2. Exclude edge pads
+  // 3. Exclude channels with too large variations
+  // 4. Exclude pads with too small signal
+  // 5. Exclude signal with outlyers RMS
+  // 6. Exclude channels to far from the chamber median        
+  noutliersCE=0;
+  //create outlier map
+  AliTPCCalPad *out=ceOut;
+  if (!out)     out= new AliTPCCalPad("outCE","outCE");
+  AliTPCCalROC *rocMasked=0x0; 
+  if (!fCETmean) return 0;
+  if (!fCETrms) return 0;
+  if (!fCEQmean) return 0;
+  //
+  //loop over all channels
+  //
+  Double_t rmsMedian         = fCETrms->GetMedian();
+  for (UInt_t iroc=0;iroc<fCETmean->kNsec;++iroc){
+    AliTPCCalROC *rocData=fCETmean->GetCalROC(iroc);
+    if (fALTROMasked) rocMasked= fALTROMasked->GetCalROC(iroc);
+    AliTPCCalROC *rocOut       = out->GetCalROC(iroc);
+    AliTPCCalROC *rocCEQ       = fCEQmean->GetCalROC(iroc);
+    AliTPCCalROC *rocCETrms    = fCETrms->GetCalROC(iroc);
+    Double_t trocMedian        = rocData->GetMedian();
+    //
+    if (!rocData) {
+      noutliersCE+=AliTPCROC::Instance()->GetNChannels(iroc);
+      rocOut->Add(1.);
+      continue;
+    }
+    //
+    //select outliers
+    UInt_t nrows=rocData->GetNrows();
+    for (UInt_t irow=0;irow<nrows;++irow){
+      UInt_t npads=rocData->GetNPads(irow);
+      for (UInt_t ipad=0;ipad<npads;++ipad){
+        rocOut->SetValue(irow,ipad,0);
+        Float_t valTmean=rocData->GetValue(irow,ipad);
+        Float_t valQmean=rocCEQ->GetValue(irow,ipad);
+        Float_t valTrms =rocCETrms->GetValue(irow,ipad);
+        //0. exclude masked pads
+        if (rocMasked && rocMasked->GetValue(irow,ipad)) {
+          rocOut->SetValue(irow,ipad,1);
+          continue;
+        }
+        //1. exclude first two rows in IROC and last two rows in OROC
+        if (iroc<36){
+          if (irow<2) rocOut->SetValue(irow,ipad,1);
+        } else {
+          if (irow>nrows-3) rocOut->SetValue(irow,ipad,1);
+        }
+        //2. exclude edge pads
+        if (ipad==0||ipad==npads-1) rocOut->SetValue(irow,ipad,1);
+        //exclude values that are exactly 0
+        if (valTmean==0) {
+          rocOut->SetValue(irow,ipad,1);
+          ++noutliersCE;
+        }
+        //3.  exclude channels with too large variations
+        if (TMath::Abs(valTmean)>fCETmaxLimitAbs) {
+          rocOut->SetValue(irow,ipad,1);
+          ++noutliersCE;
+        }
+       //
+        //4.  exclude channels with too small signal
+        if (valQmean<minSignal) {
+          rocOut->SetValue(irow,ipad,1);
+          ++noutliersCE;
+        }
+        //
+       //5. exclude channels with too small rms
+       if (valTrms<cutTrmsMin*rmsMedian || valTrms>cutTrmsMax*rmsMedian){
+         rocOut->SetValue(irow,ipad,1);
+          ++noutliersCE;
+       }
+        //
+       //6. exclude channels to far from the chamber median    
+       if (TMath::Abs(valTmean-trocMedian)>cutMaxDistT){
+         rocOut->SetValue(irow,ipad,1);
+          ++noutliersCE;
+       }
+      }
+    }
+  }
+  //
+  return out;
+}
+
+
+AliTPCCalPad *AliTPCcalibDButil::CreatePulserOutlyerMap(Int_t &noutliersPulser, AliTPCCalPad *pulserOut,Float_t cutTime, Float_t cutnRMSQ, Float_t cutnRMSrms){
+  //
+  // Author: marian.ivanov@cern.ch
+  //
+  // Create outlier map for Pulser
+  // Parameters:
+  //  Return value     - outlyer map
+  //  noutlyersPulser  - number of outlyers
+  //  cutTime          - absolute cut - distance to the median of chamber
+  //  cutnRMSQ         - nsigma cut from median  q distribution per chamber
+  //  cutnRMSrms       - nsigma cut from median  rms distribution 
+  // Outlyers criteria:
+  // 0. Exclude masked pads
+  // 1. Exclude time outlyers (default 3 time bins)
+  // 2. Exclude q outlyers    (default 5 sigma)
+  // 3. Exclude rms outlyers  (default 5 sigma)
+  noutliersPulser=0;
+  AliTPCCalPad *out=pulserOut;
+  if (!out)     out= new AliTPCCalPad("outPulser","outPulser");
+  AliTPCCalROC *rocMasked=0x0; 
+  if (!fPulserTmean) return 0;
+  if (!fPulserTrms) return 0;
+  if (!fPulserQmean) return 0;
+  //
+  //loop over all channels
+  //
+  for (UInt_t iroc=0;iroc<fCETmean->kNsec;++iroc){
+    if (fALTROMasked)   rocMasked= fALTROMasked->GetCalROC(iroc);
+    AliTPCCalROC *rocData       = fPulserTmean->GetCalROC(iroc);
+    AliTPCCalROC *rocOut        = out->GetCalROC(iroc);
+    AliTPCCalROC *rocPulserQ    = fPulserQmean->GetCalROC(iroc);
+    AliTPCCalROC *rocPulserTrms = fPulserTrms->GetCalROC(iroc);
+    //
+    Double_t rocMedianT         = rocData->GetMedian();
+    Double_t rocMedianQ         = rocPulserQ->GetMedian();
+    Double_t rocRMSQ            = rocPulserQ->GetRMS();
+    Double_t rocMedianTrms      = rocPulserTrms->GetMedian();
+    Double_t rocRMSTrms         = rocPulserTrms->GetRMS();
+    for (UInt_t ichannel=0;ichannel<rocData->GetNchannels();++ichannel){
+      rocOut->SetValue(ichannel,0);
+      Float_t valTmean=rocData->GetValue(ichannel);
+      Float_t valQmean=rocPulserQ->GetValue(ichannel);
+      Float_t valTrms =rocPulserTrms->GetValue(ichannel);
+      Int_t isOut=0;
+      if (TMath::Abs(valTmean-rocMedianT)>cutTime) isOut=1;
+      if (TMath::Abs(valQmean-rocMedianQ)>cutnRMSQ*rocRMSQ) isOut=1;
+      if (TMath::Abs(valTrms-rocMedianTrms)>cutnRMSrms*rocRMSTrms) isOut=1;
+      rocOut->SetValue(ichannel,isOut);
+      if (isOut) noutliersPulser++;
+    }
+  }
+  return out;
+}
+
+
+AliTPCCalPad *AliTPCcalibDButil::CreatePadTime0CE(TVectorD &fitResultsA, TVectorD&fitResultsC, Int_t &nOut, Double_t &chi2A, Double_t &chi2C, const char *dumpfile){
+  //
+  // Author : Marian Ivanov
+  // Create pad time0 correction map using information from the CE and from pulser
+  //
+  //
+  // Return PadTime0 to be used for time0 relative alignment
+  // if dump file specified intermediat results are dumped to the fiel and can be visualized 
+  // using $ALICE_ROOT/TPC/script/gui application
+  //
+  // fitResultsA - fitParameters A side
+  // fitResultsC - fitParameters C side
+  // chi2A       - chi2/ndf for A side (assuming error 1 time bin)
+  // chi2C       - chi2/ndf for C side (assuming error 1 time bin)
+  //
+  //
+  // Algorithm:
+  // 1. Find outlier map for CE
+  // 2. Find outlier map for Pulser
+  // 3. Replace outlier by median at given sector  (median without outliers)
+  // 4. Substract from the CE data pulser
+  // 5. Fit the CE with formula
+  //    5.1) (IROC-OROC) offset
+  //    5.2) gx
+  //    5.3) gy
+  //    5.4) (lx-xmid)
+  //    5.5) (IROC-OROC)*(lx-xmid)
+  //    5.6) (ly/lx)^2
+  // 6. Substract gy fit dependence from the CE data
+  // 7. Add pulser back to CE data  
+  // 8. Replace outliers by fit value - median of diff per given chamber -GY fit
+  // 9. return CE data
+  //
+  // Time0 <= padCE = padCEin  -padCEfitGy  - if not outlier
+  // Time0 <= padCE = padFitAll-padCEfitGy  - if outlier 
+
+  // fit formula
+  const char *formulaIn="(-1.+2.*(sector<36))*0.5++gx++gy++(lx-134.)++(-1.+2.*(sector<36))*0.5*(lx-134)++((ly/lx)^2/(0.1763)^2)";
+  // output for fit formula
+  const char *formulaAll="1++(-1.+2.*(sector<36))*0.5++gx++gy++(lx-134.)++(-1.+2.*(sector<36))*0.5*(lx-134)++((ly/lx)^2/(0.1763)^2)";
+  // gy part of formula
+  const char *formulaOut="0++0*(-1.+2.*(sector<36))*0.5++0*gx++gy++0*(lx-134.)++0*(-1.+2.*(sector<36))*0.5*(lx-134)++0*((ly/lx)^2/(0.1763)^2)";
+  //
+  //
+  if (!fCETmean) return 0;
+  Double_t pgya,pgyc,pchi2a,pchi2c;
+  AliTPCCalPad * padPulserOut = CreatePulserOutlyerMap(nOut);
+  AliTPCCalPad * padCEOut     = CreateCEOutlyerMap(nOut);
+
+  AliTPCCalPad * padPulser    = CreatePadTime0(1,pgya,pgyc,pchi2a,pchi2c);
+  AliTPCCalPad * padCE        = new AliTPCCalPad(*fCETmean);
+  AliTPCCalPad * padCEIn      = new AliTPCCalPad(*fCETmean);
+  AliTPCCalPad * padOut       = new AliTPCCalPad("padOut","padOut");   
+  padPulser->SetName("padPulser");
+  padPulserOut->SetName("padPulserOut");
+  padCE->SetName("padCE");
+  padCEIn->SetName("padCEIn");
+  padCEOut->SetName("padCEOut");
+  padOut->SetName("padOut");
+
+  //
+  // make combined outlyers map
+  // and replace outlyers in maps with median for chamber
+  //
+  for (UInt_t iroc=0;iroc<fCETmean->kNsec;++iroc){  
+    AliTPCCalROC * rocOut       = padOut->GetCalROC(iroc);
+    AliTPCCalROC * rocPulser    = padPulser->GetCalROC(iroc);
+    AliTPCCalROC * rocPulserOut = padPulserOut->GetCalROC(iroc);
+    AliTPCCalROC * rocCEOut     = padCEOut->GetCalROC(iroc);
+    AliTPCCalROC * rocCE        = padCE->GetCalROC(iroc);
+    Double_t ceMedian           = rocCE->GetMedian(rocCEOut);
+    Double_t pulserMedian       = rocPulser->GetMedian(rocCEOut);
+    for (UInt_t ichannel=0;ichannel<rocOut->GetNchannels();++ichannel){
+      if (rocPulserOut->GetValue(ichannel)>0) {
+       rocPulser->SetValue(ichannel,pulserMedian);  
+       rocOut->SetValue(ichannel,1);
+      }
+      if (rocCEOut->GetValue(ichannel)>0) {
+       rocCE->SetValue(ichannel,ceMedian);
+       rocOut->SetValue(ichannel,1);
+      }
+    }
+  }
+  //
+  // remove pulser time 0
+  //
+  padCE->Add(padPulser,-1);
+  //
+  // Make fits
+  //
+  TMatrixD dummy;
+  Float_t chi2Af,chi2Cf;  
+  padCE->GlobalSidesFit(padOut,formulaIn,fitResultsA,fitResultsC,dummy,dummy,chi2Af,chi2Cf);
+  chi2A=chi2Af;
+  chi2C=chi2Cf;
+  //
+  AliTPCCalPad *padCEFitGY=AliTPCCalPad::CreateCalPadFit(formulaOut,fitResultsA,fitResultsC);
+  padCEFitGY->SetName("padCEFitGy");
+  //
+  AliTPCCalPad *padCEFit  =AliTPCCalPad::CreateCalPadFit(formulaAll,fitResultsA,fitResultsC);
+  padCEFit->SetName("padCEFit");
+  //
+  AliTPCCalPad* padCEDiff  = new AliTPCCalPad(*padCE);
+  padCEDiff->SetName("padCEDiff");
+  padCEDiff->Add(padCEFit,-1.);
+  //
+  // 
+  padCE->Add(padCEFitGY,-1.);
+
+  padCE->Add(padPulser,1.);  
+  Double_t padmedian = padCE->GetMedian();
+  padCE->Add(-padmedian);  // normalize to median
+  //
+  // Replace outliers by fit value - median of diff per given chamber -GY fit
+  //
+  for (UInt_t iroc=0;iroc<fCETmean->kNsec;++iroc){  
+    AliTPCCalROC * rocOut       = padOut->GetCalROC(iroc);
+    AliTPCCalROC * rocCE        = padCE->GetCalROC(iroc);
+    AliTPCCalROC * rocCEFit     = padCEFit->GetCalROC(iroc);
+    AliTPCCalROC * rocCEFitGY   = padCEFitGY->GetCalROC(iroc);
+    AliTPCCalROC * rocCEDiff    = padCEDiff->GetCalROC(iroc);
+    //
+    Double_t diffMedian         = rocCEDiff->GetMedian(rocOut);
+    for (UInt_t ichannel=0;ichannel<rocOut->GetNchannels();++ichannel){
+      if (rocOut->GetValue(ichannel)==0) continue;
+      Float_t value=rocCEFit->GetValue(ichannel)-rocCEFitGY->GetValue(ichannel)-diffMedian-padmedian;
+      rocCE->SetValue(ichannel,value);
+    }    
+  }
+  //
+  //
+  if (dumpfile){
+    //dump to the file - result can be visualized
+    AliTPCPreprocessorOnline preprocesor;
+    preprocesor.AddComponent(new AliTPCCalPad(*padCE));
+    preprocesor.AddComponent(new AliTPCCalPad(*padCEIn));
+    preprocesor.AddComponent(new AliTPCCalPad(*padCEFit));
+    preprocesor.AddComponent(new AliTPCCalPad(*padOut));
+    //
+    preprocesor.AddComponent(new AliTPCCalPad(*padCEFitGY));
+    preprocesor.AddComponent(new AliTPCCalPad(*padCEDiff));
+    //
+    preprocesor.AddComponent(new AliTPCCalPad(*padCEOut));
+    preprocesor.AddComponent(new AliTPCCalPad(*padPulser));
+    preprocesor.AddComponent(new AliTPCCalPad(*padPulserOut));
+    preprocesor.DumpToFile(dumpfile);
+  } 
+  delete padPulser;
+  delete padPulserOut;
+  delete padCEIn;
+  delete padCEOut;
+  delete padOut;
+  delete padCEDiff;
+  delete padCEFitGY;
+  return padCE;
+}
+
index d0f49a5..3bcfd63 100644 (file)
@@ -28,7 +28,7 @@ public:
   void UpdateFromCalibDB();
   //data processing functions
   void ProcessCEdata(const char* fitFormula, TVectorD &fitResultsA, TVectorD &fitResultsC,
-                     Int_t &noutliersCE, AliTPCCalPad *outCE=0);
+                     Int_t &noutliersCE, Double_t & chi2A, Double_t &chi2C, AliTPCCalPad *outCE=0);
   void ProcessCEgraphs(TVectorD &vecTEntries, TVectorD &vecTMean, TVectorD &vecTRMS, TVectorD &vecTMedian,
                        TVectorD &vecQEntries, TVectorD &vecQMean, TVectorD &vecQRMS, TVectorD &vecQMedian,
                        Float_t &driftTimeA, Float_t &driftTimeC );
@@ -68,8 +68,16 @@ public:
                 {fRefALTROMasked=masked;}
   
   //creation of derived pad by pad calibration data
-  AliTPCCalPad *CreatePadTime0(Int_t model=0);
+  AliTPCCalPad *CreatePadTime0(Int_t model, Double_t &gyA, Double_t &gyC, Double_t &chi2A, Double_t &chi2C);
   //
+  // create outlyer maps
+  //
+  AliTPCCalPad *CreateCEOutlyerMap(Int_t &noutliersCE, AliTPCCalPad *ceOut=0, Float_t minSignal=10, Float_t cutTrmsMin=0.9, Float_t cutTrmsMax=1.2, Float_t cutMaxDistT=0.7);
+  AliTPCCalPad *CreatePulserOutlyerMap(Int_t &noutliersPulser, AliTPCCalPad *pulserOut=0, Float_t cutTime=3, Float_t cutnRMSQ=5, Float_t cutnRMSrms=5);
+  //
+  AliTPCCalPad *CreatePadTime0CE(TVectorD &fitResultsA, TVectorD&fitResultsC, Int_t &nOut, Double_t &chi2A, Double_t &chi2C, const char *dumpfile=0);
+  //
+
   void UpdatePulserOutlierMap();
   void UpdateRefPulserOutlierMap();
   void PulserOutlierMap(AliTPCCalPad *pulOut, const AliTPCCalPad *pulT, const AliTPCCalPad *pulQ);
index 8633053..6a91ee4 100644 (file)
@@ -1,11 +1,9 @@
 /*
 //Make a tree dump of TPC calibration:
 
-
 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC");
 .L $ALICE_ROOT/TPC/CalibMacros/CalibEnv.C+
 
-
 CalibEnv("run.list");
 TFile f("dcsTime.root")
 */
@@ -44,10 +42,17 @@ TFile f("dcsTime.root")
 #include "TMatrixD.h"
 #include "AliTPCCalibRaw.h"
 #include "AliSplineFit.h"
+#include "TGraphErrors.h"
 
-
+//
+AliTPCcalibDB     *calibDB=0;
+AliTPCcalibDButil *dbutil =0;
 TTree * dcsTree=0;
-TString refFile="dummy.root";
+TString refFile="dummy.root"; 
+TTreeSRedirector *pcstream =0;
+//
+//
+void ProcessRun(Int_t irun, Int_t startTime, Int_t endTime);
 void GetProductionInfo(Int_t run, Int_t &nalien, Int_t &nRawAlien, Int_t &nlocal, Int_t &nRawLocal);
   
 
@@ -57,9 +62,9 @@ void CalibEnv(const char * runList, Int_t first=1, Int_t last=-1){
   // first   - first run to process
   // last    - last  to process
   // 
-
   refFile=gSystem->Getenv("REF_DATA_FILE");
-  AliTPCcalibDB * calibDB = AliTPCcalibDB::Instance();
+  calibDB = AliTPCcalibDB::Instance();
+  dbutil= new AliTPCcalibDButil;
   //
   // make list of runs
   //
@@ -78,264 +83,319 @@ void CalibEnv(const char * runList, Int_t first=1, Int_t last=-1){
     nruns++;
   }
   TMath::Sort(nruns, runArray.fArray, indexes,kFALSE);
-
-  TTreeSRedirector *pcstream = new TTreeSRedirector("dcsTime.root");
-  AliTPCcalibDButil dbutil;
-  dbutil.SetRefFile(refFile.Data());
+  pcstream = new TTreeSRedirector("dcsTime.root");
+  dbutil->SetRefFile(refFile.Data());
   Int_t startTime = 0;
-  Int_t endTime = 0;
-  Int_t startTimeGRP = 0;
-  Int_t stopTimeGRP  = 0;
-
-  AliSplineFit *fitVdrift=0x0;
+  Int_t endTime   = 0;
   for (Int_t run=0; run<nruns; run++){
     Int_t irun=runArray[indexes[run]];
     printf("Processing run %d ...\n",irun);
     AliTPCcalibDB::Instance()->SetRun(irun);
-    dbutil.UpdateFromCalibDB();
+    dbutil->UpdateFromCalibDB();
     //
     AliDCSSensorArray *arrHV=calibDB->GetVoltageSensors(irun);
     if (!arrHV) continue;
     for  (Int_t isenHV=0; isenHV<arrHV->NumSensors(); ++isenHV){
       AliDCSSensor *senHV=arrHV->GetSensorNum(isenHV);
-      if (!senHV) continue;
+      if (!senHV) {
+       printf("Not interesting OCDB info\n");
+       continue;
+      }
       startTime=senHV->GetStartTime();
       endTime  =senHV->GetEndTime();
       if (startTime>0&&endTime>0) break;
     }
-    if (calibDB->GetGRP(irun)){
-      startTimeGRP = AliTPCcalibDB::GetGRP(irun)->GetTimeStart();
-      stopTimeGRP  = AliTPCcalibDB::GetGRP(irun)->GetTimeEnd();
+    ProcessRun(irun, startTime,endTime);
+  }
+  delete pcstream;  
+}
+
+
+void ProcessRun(Int_t irun, Int_t startTime, Int_t endTime){
+  //
+  //
+  // 
+  AliSplineFit *fitVdrift=0x0;
+  Int_t startTimeGRP=0, stopTimeGRP=0;
+  if (calibDB->GetGRP(irun)){
+    startTimeGRP = AliTPCcalibDB::GetGRP(irun)->GetTimeStart();
+    stopTimeGRP  = AliTPCcalibDB::GetGRP(irun)->GetTimeEnd();
+  }
+  
+  AliTPCSensorTempArray * tempArray = AliTPCcalibDB::Instance()->GetTemperatureSensor(irun);
+  AliTPCTempMap * tempMap = new AliTPCTempMap(tempArray);
+  AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(irun);
+  //
+  Int_t dtime = TMath::Max((endTime-startTime)/20,10*60);
+  //
+  //Goofie statistical data
+  //
+  TVectorD vecEntries, vecMean, vecMedian,vecRMS;
+  dbutil->ProcessGoofie(vecEntries ,vecMedian, vecMean, vecRMS);
+  //
+  //CE data processing - see ProcessCEdata function for description of the results
+  //
+  TVectorD fitResultsA, fitResultsC;
+  Int_t nmaskedCE;
+  Double_t chi2ACE=0,chi2CCE=0;
+  //     dbutil->ProcessCEdata("(sector<36)++gx++gy++lx++lx**2",fitResultsA,fitResultsC,nmaskedCE);
+  dbutil->ProcessCEdata("(sector<36)++gx++gy++(lx-134)++(sector<36)*(lx-134)++(ly/lx)^2",fitResultsA,fitResultsC,nmaskedCE,chi2ACE,chi2CCE);
+
+  TVectorD fitCEResultsA(7), fitCEResultsC(7);
+  Int_t    noutCE;
+  Double_t chi2CEA=0,chi2CEC=0;
+  AliTPCCalPad *time0 = dbutil->CreatePadTime0CE(fitCEResultsA, fitCEResultsC, noutCE, chi2CEA, chi2CEC);
+  delete time0;
+  //  
+  //
+  TVectorD vecTEntries, vecTMean, vecTRMS, vecTMedian, vecQEntries, vecQMean, vecQRMS, vecQMedian;
+  Float_t driftTimeA, driftTimeC;
+  dbutil->ProcessCEgraphs(vecTEntries, vecTMean, vecTRMS, vecTMedian,
+                         vecQEntries, vecQMean, vecQRMS, vecQMedian,
+                         driftTimeA, driftTimeC );
+  //
+  //
+  //
+  //drift velocity using tracks
+  //
+  //     fitVdrift=calibDB->GetVdriftSplineFit("ALISPLINEFIT_MEAN_VDRIFT_COSMICS_ALL",irun);
+  fitVdrift=calibDB->CreateVdriftSplineFit("TGRAPHERRORS_MEAN_VDRIFT_COSMICS_ALL",irun);
+  //noise data Processing - see ProcessNoiseData function for description of the results
+  TVectorD vNoiseMean, vNoiseMeanSenRegions, vNoiseRMS, vNoiseRMSSenRegions;
+  Int_t nonMaskedZero=0;
+  dbutil->ProcessNoiseData(vNoiseMean, vNoiseMeanSenRegions, vNoiseRMS, vNoiseRMSSenRegions, nonMaskedZero);
+  //
+  // comparisons
+  //
+  TVectorF pedestalDeviations;
+  TVectorF noiseDeviations;
+  TVectorF pulserQdeviations;
+  Float_t varQMean;
+  Int_t npadsOutOneTB;
+  Int_t npadsOffAdd;
+  dbutil->ProcessPedestalVariations(pedestalDeviations);
+  dbutil->ProcessNoiseVariations(noiseDeviations);
+  dbutil->ProcessPulserVariations(pulserQdeviations,varQMean,npadsOutOneTB,npadsOffAdd);
+  //
+  //L3 data 
+  //
+  Float_t bz=AliTPCcalibDB::GetBz(irun);
+  Char_t  l3pol=AliTPCcalibDB::GetL3Polarity(irun);
+  //
+  //calibration Pulser data processing
+  //
+  Int_t nOffChannels=0;
+  TVectorD vTimePulser;
+  nOffChannels=dbutil->GetNPulserOutliers();
+  dbutil->ProcessPulser(vTimePulser);
+  //
+  //ALTRO data
+  //
+  Int_t nMasked=0;
+  dbutil->ProcessALTROConfig(nMasked);
+  //
+  //Calib RAW data
+  //
+  Int_t nFailL1=-1;
+  if (calibDB->GetCalibRaw()) nFailL1=calibDB->GetCalibRaw()->GetNFailL1Phase();
+  //
+  //production information
+  //
+  Int_t nalien=0,nRawAlien=0,nlocal=0,nRawLocal=0;
+  //     GetProductionInfo(irun, nalien, nRawAlien, nlocal,nRawLocal);
+  //run type
+  TObjString runType(AliTPCcalibDB::GetRunType(irun).Data());
+  //
+  //
+  //
+  
+  for (Int_t itime=startTime; itime<endTime; itime+=dtime){
+    //
+    TTimeStamp tstamp(itime);
+    Float_t valuePressure  = calibDB->GetPressure(tstamp,irun,0);
+    Float_t valuePressure2 = calibDB->GetPressure(tstamp,irun,1);
+    Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative((UInt_t)itime,irun,0);
+    Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative((UInt_t)itime,irun,1);
+    //temperature fits
+    TLinearFitter * fitter = 0;
+    TVectorD vecTemp[10];
+    for (Int_t itype=0; itype<5; itype++)
+      for (Int_t iside=0; iside<2; iside++){
+       fitter= tempMap->GetLinearFitter(itype,iside,tstamp);     
+       if (!fitter) continue;
+       fitter->Eval();
+       fitter->GetParameters(vecTemp[itype+iside*5]);
+       delete fitter;
+      }
+    //
+    //measured skirt temperatures
+    //
+    TVectorD vecSkirtTempA(18);
+    TVectorD vecSkirtTempC(18);
+    Int_t nsenTemp=tempArray->NumSensors();
+    for (Int_t isenTemp=0;isenTemp<nsenTemp;++isenTemp){
+      AliTPCSensorTemp *senTemp=(AliTPCSensorTemp*)tempArray->GetSensorNum(isenTemp);
+      if (senTemp->GetType()!=3) continue;
+      if (TMath::Sqrt(senTemp->GetX()*senTemp->GetX()+senTemp->GetY()*senTemp->GetY())<100) continue; //only skirt, outer FC vessel
+      Double_t val=senTemp->GetValue(tstamp);
+      if (senTemp->GetSide()==0)
+       vecSkirtTempA[senTemp->GetSector()]=val;
+      else
+       vecSkirtTempC[senTemp->GetSector()]=val;
+    }
+    //
+    //goofie data
+    //
+    TVectorD vecGoofie;
+    if (goofieArray){
+      vecGoofie.ResizeTo(goofieArray->NumSensors());
+      for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
+       AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
+       if (gsensor){
+         vecGoofie[isensor] = gsensor->GetValue(tstamp);
+       }
+      }
+    } else {
+      vecGoofie.ResizeTo(19);
     }
-    //    AliDCSSensor * sensorPressure = AliTPCcalibDB::Instance()->GetPressureSensor(irun);
-//     if (!sensorPressure) continue;
-//     Int_t startTime = sensorPressure->GetStartTime();
-//     Int_t endTime = sensorPressure->GetEndTime();
-//     Int_t startTimeGRP = AliTPCcalibDB::GetGRP(irun)->GetTimeStart();
-//     Int_t stopTimeGRP  = AliTPCcalibDB::GetGRP(irun)->GetTimeEnd();
-    AliTPCSensorTempArray * tempArray = AliTPCcalibDB::Instance()->GetTemperatureSensor(irun);
-    AliTPCTempMap * tempMap = new AliTPCTempMap(tempArray);
-    AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(irun);
     //
-    Int_t dtime = TMath::Max((endTime-startTime)/20,10*60);
-    //Goofie statistical data
-    TVectorD vecEntries, vecMean, vecMedian,vecRMS;
-    dbutil.ProcessGoofie(vecEntries ,vecMedian, vecMean, vecRMS);
-    //CE data processing - see ProcessCEdata function for description of the results
-    TVectorD fitResultsA, fitResultsC;
-    Int_t nmaskedCE;
-    Double_t chi2ACE=0,chi2CCE=0;
-//     dbutil.ProcessCEdata("(sector<36)++gx++gy++lx++lx**2",fitResultsA,fitResultsC,nmaskedCE);
-    dbutil.ProcessCEdata("(sector<36)++gy++gx++(lx-134)++(sector<36)*(lx-134)++(ly/lx)^2",fitResultsA,fitResultsC,nmaskedCE,chi2ACE,chi2CCE);
-    TVectorD vecTEntries, vecTMean, vecTRMS, vecTMedian, vecQEntries, vecQMean, vecQRMS, vecQMedian;
-    Float_t driftTimeA, driftTimeC;
-    dbutil.ProcessCEgraphs(vecTEntries, vecTMean, vecTRMS, vecTMedian,
-                           vecQEntries, vecQMean, vecQRMS, vecQMedian,
-                           driftTimeA, driftTimeC );
-    //drift velocity using tracks
-//     fitVdrift=calibDB->GetVdriftSplineFit("ALISPLINEFIT_MEAN_VDRIFT_COSMICS_ALL",irun);
-    fitVdrift=calibDB->CreateVdriftSplineFit("TGRAPHERRORS_MEAN_VDRIFT_COSMICS_ALL",irun);
-    //noise data Processing - see ProcessNoiseData function for description of the results
-    TVectorD vNoiseMean, vNoiseMeanSenRegions, vNoiseRMS, vNoiseRMSSenRegions;
-    Int_t nonMaskedZero=0;
-    dbutil.ProcessNoiseData(vNoiseMean, vNoiseMeanSenRegions, vNoiseRMS, vNoiseRMSSenRegions, nonMaskedZero);
-    // comparisons
-    TVectorF pedestalDeviations;
-    TVectorF noiseDeviations;
-    TVectorF pulserQdeviations;
-    Float_t varQMean;
-    Int_t npadsOutOneTB;
-    Int_t npadsOffAdd;
-    dbutil.ProcessPedestalVariations(pedestalDeviations);
-    dbutil.ProcessNoiseVariations(noiseDeviations);
-    dbutil.ProcessPulserVariations(pulserQdeviations,varQMean,npadsOutOneTB,npadsOffAdd);
+    TVectorD voltagesIROC(36);
+    TVectorD voltagesOROC(36);
+    for(Int_t j=1; j<36; j++) voltagesIROC[j-1] = AliTPCcalibDB::Instance()->GetChamberHighVoltage(irun, j,itime);
+    for(Int_t j=36; j<72; j++) voltagesOROC[j-36] = AliTPCcalibDB::Instance()->GetChamberHighVoltage(irun, j,itime);
+    Double_t voltIROC = TMath::Median(36, voltagesIROC.GetMatrixArray());
+    Double_t voltOROC = TMath::Median(36, voltagesOROC.GetMatrixArray());
+    //
+    Float_t  coverIA=AliTPCcalibDB::GetCoverVoltage(irun,0,itime);
+    Float_t  coverIC=AliTPCcalibDB::GetCoverVoltage(irun,18,itime);
+    Float_t  coverOA=AliTPCcalibDB::GetCoverVoltage(irun,36,itime);
+    Float_t  coverOC=AliTPCcalibDB::GetCoverVoltage(irun,54,itime);
+    Float_t  skirtA=AliTPCcalibDB::GetSkirtVoltage(irun,0,itime);
+    Float_t  skirtC=AliTPCcalibDB::GetSkirtVoltage(irun,18,itime);
+    Float_t  ggOffA=AliTPCcalibDB::GetGGoffsetVoltage(irun,0,itime);
+    Float_t  ggOffC=AliTPCcalibDB::GetGGoffsetVoltage(irun,18,itime);
+    //drift velocity
+    Float_t dvCorr=-5;
+    if (fitVdrift) dvCorr=fitVdrift->Eval(itime);
+    //
+    // Gain - Alexander Kalveit
+    //
+    Float_t  gainCosmic = 0;
+    Float_t  gainMIP = 0;
+    TObjArray * gainSplines = AliTPCcalibDB::Instance()->GetTimeGainSplinesRun(irun);
+    if (gainSplines) {
+      TGraphErrors * graphMIP = (TGraphErrors *) gainSplines->FindObject("TGRAPHERRORS_MEAN_GAIN_BEAM_ALL");
+      TGraphErrors * graphCosmic = (TGraphErrors *) gainSplines->FindObject("TGRAPHERRORS_MEAN_GAIN_COSMIC_ALL");
+      if (graphMIP) gainMIP = graphMIP->Eval(itime);
+      if (graphCosmic) gainCosmic = graphCosmic->Eval(itime);
+    }
     
-    //L3 data 
-    Float_t bz=AliTPCcalibDB::GetBz(irun);
-    Char_t  l3pol=AliTPCcalibDB::GetL3Polarity(irun);
-    //calibration Pulser data processing
-    Int_t nOffChannels=0;
-    TVectorD vTimePulser;
-    nOffChannels=dbutil.GetNPulserOutliers();
-    dbutil.ProcessPulser(vTimePulser);
-    //ALTRO data
-    Int_t nMasked=0;
-    dbutil.ProcessALTROConfig(nMasked);
-    //Calib RAW data
-    Int_t nFailL1=-1;
-    if (calibDB->GetCalibRaw())
-      nFailL1=calibDB->GetCalibRaw()->GetNFailL1Phase();
-    //production information
-    Int_t nalien=0,nRawAlien=0,nlocal=0,nRawLocal=0;
-//     GetProductionInfo(irun, nalien, nRawAlien, nlocal,nRawLocal);
-    //run type
-    TObjString runType(AliTPCcalibDB::GetRunType(irun).Data());
     
-    for (Int_t itime=startTime; itime<endTime; itime+=dtime){
+    //tempMap->GetLinearFitter(0,0,itime);
+    (*pcstream)<<"dcs"<<
+      "run="<<irun<<
+      "time="<<itime<<
+      "startTimeGRP="<<startTimeGRP<<
+      "stopTimeGRP="<<stopTimeGRP<<
+      //run type
+      "runType.="<<&runType<<
+      // voltage setting
+      "VIROC.="<<&voltagesIROC<<
+      "VOROC.="<<&voltagesOROC<<
+      "medianVIROC="<<voltIROC<<
+      "medianVOROC="<<voltOROC<<
+      "coverIA=" << coverIA <<
+      "coverIC=" << coverIC <<
+      "coverOA=" << coverOA <<
+      "coverOC=" << coverOC <<
+      "skirtA=" << skirtA <<
+      "skirtC=" << skirtC <<
+      "ggOffA=" << ggOffA <<
+      "ggOffC=" << ggOffC <<
       //
-      TTimeStamp tstamp(itime);
-      Float_t valuePressure  = calibDB->GetPressure(tstamp,irun,0);
-      Float_t valuePressure2 = calibDB->GetPressure(tstamp,irun,1);
-      //temperature fits
-      TLinearFitter * fitter = 0;
-      TVectorD vecTemp[10];
-      for (Int_t itype=0; itype<5; itype++)
-        for (Int_t iside=0; iside<2; iside++){
-          fitter= tempMap->GetLinearFitter(itype,iside,tstamp);          
-          if (!fitter) continue;
-          fitter->Eval();
-          fitter->GetParameters(vecTemp[itype+iside*5]);
-          delete fitter;
-        }
-      //measured skirt temperatures
-      TVectorD vecSkirtTempA(18);
-      TVectorD vecSkirtTempC(18);
-      Int_t nsenTemp=tempArray->NumSensors();
-      for (Int_t isenTemp=0;isenTemp<nsenTemp;++isenTemp){
-        AliTPCSensorTemp *senTemp=(AliTPCSensorTemp*)tempArray->GetSensorNum(isenTemp);
-        if (senTemp->GetType()!=3) continue;
-        if (TMath::Sqrt(senTemp->GetX()*senTemp->GetX()+senTemp->GetY()*senTemp->GetY())<100) continue; //only skirt, outer FC vessel
-        Double_t val=senTemp->GetValue(tstamp);
-        if (senTemp->GetSide()==0)
-          vecSkirtTempA[senTemp->GetSector()]=val;
-        else
-          vecSkirtTempC[senTemp->GetSector()]=val;
-      }
-      //goofie data
-      TVectorD vecGoofie;
-      if (goofieArray){
-        vecGoofie.ResizeTo(goofieArray->NumSensors());
-        for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
-          AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
-          if (gsensor){
-            vecGoofie[isensor] = gsensor->GetValue(tstamp);
-          }
-        }
-      } else {
-        vecGoofie.ResizeTo(19);
-      }
-      Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative((UInt_t)itime,irun,0);
-      Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative((UInt_t)itime,irun,1);
+      "ptrel0="<<ptrelative0<<  // deltaTP/TP  - A side
+      "ptrel1="<<ptrelative1<<  // deltaTP/TPC - C side
+      "goofie.="<<&vecGoofie<<
+      "goofieE.="<<&vecEntries<<
+      "goofieMean.="<<&vecMean<<
+      "goofieMedian.="<<&vecMedian<<
+      "goofieRMS.="<<&vecRMS<<
       //
-      TVectorD voltagesIROC(36);
-      TVectorD voltagesOROC(36);
-      for(Int_t j=1; j<36; j++) voltagesIROC[j-1] = AliTPCcalibDB::Instance()->GetChamberHighVoltage(irun, j,itime);
-      for(Int_t j=36; j<72; j++) voltagesOROC[j-36] = AliTPCcalibDB::Instance()->GetChamberHighVoltage(irun, j,itime);
-      Double_t voltIROC = TMath::Median(36, voltagesIROC.GetMatrixArray());
-      Double_t voltOROC = TMath::Median(36, voltagesOROC.GetMatrixArray());
+      "press="<<valuePressure<<
+      "press2="<<valuePressure2<<
+      "temp00.="<<&vecTemp[0]<<
+      "temp10.="<<&vecTemp[1]<<
+      "temp20.="<<&vecTemp[2]<<
+      "temp30.="<<&vecTemp[3]<<
+      "temp40.="<<&vecTemp[4]<<
+      "temp01.="<<&vecTemp[5]<<
+      "temp11.="<<&vecTemp[6]<<
+      "temp21.="<<&vecTemp[7]<<
+      "temp31.="<<&vecTemp[8]<<
+      "temp41.="<<&vecTemp[9]<<
+      "tempSkirtA.="<<&vecSkirtTempA<<
+      "tempSkirtC.="<<&vecSkirtTempC;
+    
+    (*pcstream)<<"dcs"<<       
+      //noise data
+      "meanNoise.="<<&vNoiseMean<<
+      "meanNoiseSen.="<<&vNoiseMeanSenRegions<<
+      "rmsNoise.="<<&vNoiseRMS<<
+      "rmsNoiseSen.="<<&vNoiseRMSSenRegions<<
+      "zeroNoise="<<nonMaskedZero<<
+      //pulser data
+      "timePulser.=" << &vTimePulser <<
+      "nOffPulser="<<nOffChannels<<
+      //altro data
+      "nMasked="<<nMasked<<
+      //ce data -Jens version
+      "CEfitA.="<<&fitResultsA<<
+      "CEfitC.="<<&fitResultsC<<
+      "nmaskedCE="<<nmaskedCE<<
+      "chi2ACE="<<chi2ACE<<
+      "chi2CCE="<<chi2CCE<<
+      //ce data new - MI version
+      "CEfitAMI.="<<&fitCEResultsA<<
+      "CEfitCMI.="<<&fitCEResultsC<<
+      "chi2CEA="<<chi2CEA<<
+      "chi2CEC="<<chi2CEC<<
       //
-      Float_t  coverIA=AliTPCcalibDB::GetCoverVoltage(irun,0,itime);
-      Float_t  coverIC=AliTPCcalibDB::GetCoverVoltage(irun,18,itime);
-      Float_t  coverOA=AliTPCcalibDB::GetCoverVoltage(irun,36,itime);
-      Float_t  coverOC=AliTPCcalibDB::GetCoverVoltage(irun,54,itime);
-      Float_t  skirtA=AliTPCcalibDB::GetSkirtVoltage(irun,0,itime);
-      Float_t  skirtC=AliTPCcalibDB::GetSkirtVoltage(irun,18,itime);
-      Float_t  ggOffA=AliTPCcalibDB::GetGGoffsetVoltage(irun,0,itime);
-      Float_t  ggOffC=AliTPCcalibDB::GetGGoffsetVoltage(irun,18,itime);
-      //drift velocity
-      Float_t dvCorr=-5;
-      if (fitVdrift) dvCorr=fitVdrift->Eval(itime);
-      
-      
-      //tempMap->GetLinearFitter(0,0,itime);
-      (*pcstream)<<"dcs"<<
-        "run="<<irun<<
-        "time="<<itime<<
-        "startTimeGRP="<<startTimeGRP<<
-        "stopTimeGRP="<<stopTimeGRP<<
-        //run type
-        "runType.="<<&runType<<
-        // voltage setting
-        "VIROC.="<<&voltagesIROC<<
-        "VOROC.="<<&voltagesOROC<<
-        "medianVIROC="<<voltIROC<<
-        "medianVOROC="<<voltOROC<<
-        "coverIA=" << coverIA <<
-        "coverIC=" << coverIC <<
-        "coverOA=" << coverOA <<
-        "coverOC=" << coverOC <<
-        "skirtA=" << skirtA <<
-        "skirtC=" << skirtC <<
-        "ggOffA=" << ggOffA <<
-        "ggOffC=" << ggOffC <<
-        //
-        "ptrel0="<<ptrelative0<<  // deltaTP/TP  - A side
-        "ptrel1="<<ptrelative1<<  // deltaTP/TPC - C side
-        "goofie.="<<&vecGoofie<<
-        "goofieE.="<<&vecEntries<<
-        "goofieMean.="<<&vecMean<<
-        "goofieMedian.="<<&vecMedian<<
-        "goofieRMS.="<<&vecRMS<<
-        //
-        "press="<<valuePressure<<
-        "press2="<<valuePressure2<<
-        "temp00.="<<&vecTemp[0]<<
-        "temp10.="<<&vecTemp[1]<<
-        "temp20.="<<&vecTemp[2]<<
-        "temp30.="<<&vecTemp[3]<<
-        "temp40.="<<&vecTemp[4]<<
-        "temp01.="<<&vecTemp[5]<<
-        "temp11.="<<&vecTemp[6]<<
-        "temp21.="<<&vecTemp[7]<<
-        "temp31.="<<&vecTemp[8]<<
-        "temp41.="<<&vecTemp[9]<<
-        "tempSkirtA.="<<&vecSkirtTempA<<
-        "tempSkirtC.="<<&vecSkirtTempC;
-      
-      (*pcstream)<<"dcs"<<     
-        //noise data
-        "meanNoise.="<<&vNoiseMean<<
-        "meanNoiseSen.="<<&vNoiseMeanSenRegions<<
-        "rmsNoise.="<<&vNoiseRMS<<
-        "rmsNoiseSen.="<<&vNoiseRMSSenRegions<<
-        "zeroNoise="<<nonMaskedZero<<
-        //pulser data
-        "timePulser.=" << &vTimePulser <<
-        "nOffPulser="<<nOffChannels<<
-        //altro data
-        "nMasked="<<nMasked<<
-        //ce data
-        "CEfitA.="<<&fitResultsA<<
-        "CEfitC.="<<&fitResultsC<<
-        "nmaskedCE="<<nmaskedCE<<
-       "chi2ACE="<<chi2ACE<<
-       "chi2CCE="<<chi2CCE<<
-        //ce graph data
-        "CEgrTEntries.="<<&vecTEntries<<
-        "CEgrTMean.="<<&vecTMean<<
-        "CEgrTRMS.="<<&vecTRMS<<
-        "CEgrTMedian.="<<&vecTMedian<<
-        "CEgrQEntries.="<<&vecQEntries<<
-        "CEgrQMean.="<<&vecQMean<<
-        "CEgrQRMS.="<<&vecQRMS<<
-        "CEgrQMedian.="<<&vecQMedian<<
-        "CEgrDriftA="<<driftTimeA<<
-        "CEgrDriftC="<<driftTimeC<<
-        //calib raw data
-        "nFailL1="<<nFailL1<<
-        // b field
-        "Bz="<< bz <<
-        "L3polarity="<<l3pol<<
-        // production information
-        "nalien="<<nalien<<
-        "nRawAlien="<<nRawAlien<<
-        "nlocal="<<nlocal<<
-        "nRawLocal="<<nRawLocal<<
-        //comparisons with ref data
-        "pedestalDeviations.="<<&pedestalDeviations<<
-        "noiseDeviations.="<<&noiseDeviations<<
-        "pulserQdeviations.="<<&pulserQdeviations<<
-//         "pulserVarQMean="<<varQMean<<
-        "pulserNpadsOutOneTB="<<npadsOutOneTB<<
-        "pulserNpadsOffAdd="<<npadsOffAdd<<
-        "driftCorrCosmAll="<<dvCorr<<
-      
-        "\n";
-    }//end run loop
-//     delete fitVdrift;
-//     fitVdrift=0;
-  }
-  delete pcstream;
+      //ce graph data
+      "CEgrTEntries.="<<&vecTEntries<<
+      "CEgrTMean.="<<&vecTMean<<
+      "CEgrTRMS.="<<&vecTRMS<<
+      "CEgrTMedian.="<<&vecTMedian<<
+      "CEgrQEntries.="<<&vecQEntries<<
+      "CEgrQMean.="<<&vecQMean<<
+      "CEgrQRMS.="<<&vecQRMS<<
+      "CEgrQMedian.="<<&vecQMedian<<
+      "CEgrDriftA="<<driftTimeA<<
+      "CEgrDriftC="<<driftTimeC<<
+      //calib raw data
+      "nFailL1="<<nFailL1<<
+      // b field
+      "Bz="<< bz <<
+      "L3polarity="<<l3pol<<
+      // production information
+      "nalien="<<nalien<<
+      "nRawAlien="<<nRawAlien<<
+      "nlocal="<<nlocal<<
+      "nRawLocal="<<nRawLocal<<
+      //comparisons with ref data
+      "pedestalDeviations.="<<&pedestalDeviations<<
+      "noiseDeviations.="<<&noiseDeviations<<
+      "pulserQdeviations.="<<&pulserQdeviations<<
+      //         "pulserVarQMean="<<varQMean<<
+      "pulserNpadsOutOneTB="<<npadsOutOneTB<<
+      "pulserNpadsOffAdd="<<npadsOffAdd<<
+      "driftCorrCosmAll="<<dvCorr<<
+      // time dependence of gain
+      "gainMIP="<<gainMIP<<
+      "gainCosmic="<<gainCosmic<<      
+      "\n";
+  }//end run loop
+  //     delete fitVdrift;
+  //     fitVdrift=0; 
 }
 
 
@@ -377,30 +437,3 @@ void GetProductionInfo(Int_t run, Int_t &nalien, Int_t &nRawAlien, Int_t &nlocal
   nRawLocal=sNlines.Atoi();
 }
 
-void FilterMag(const char * runList){
-  //
-  //
-  //
-  //  AliTPCcalibDB * calibDB = AliTPCcalibDB::Instance();
-  ifstream in;
-  in.open(runList);
-  Int_t irun=0;
-  while(in.good()) {
-    in >> irun;
-    if (irun==0) continue;
-    AliGRPObject *grp = AliTPCcalibDB::GetGRP(irun);
-    Float_t current = -1;
-    Float_t bz      = -1;
-//     Float_t press   =  0;
-    if (grp){
-      current = grp->GetL3Current((AliGRPObject::Stats)0);
-      bz = 5*current/30000.;
-      printf("Run%d\tL3 current%f\tBz\t%f\n",irun,current,bz);
-    }
-    else{
-      printf("Run%d\tL3 current%f\tBz\t%f\n",irun,current,bz);
-    }
-  }
-  
-}
-