#ifndef TPCCLUSTERERRORSCALER_H
#define TPCCLUSTERERRORSCALER_H

#include "TString.h"
#include "TVector3.h"
#include "TpcDigiPar.h"

#include "TpcGas.h"
#include "TMatrixDSym.h"
#include "TClonesArray.h"
#include "TMatrixD.h"

class TpcClusterErrorScaler
{
 private:
  static TpcClusterErrorScaler *finstance;
  TpcDigiPar * fdigiPar;
  TString fparFile;
  Double_t fpar[3][8];
  Bool_t fInit;
  Bool_t fNoFairRun;
  Bool_t fManualDiffusion;
  Int_t ffunc;
  Double_t fdx;
  Double_t fdy;
  Double_t fdz;
  Double_t fwx;
  Double_t fwy;
  Double_t fdecayX;
  Double_t fdecayY;
  Double_t fdecayZ;
  Double_t fDl;
  Double_t fDt;
  std::vector< std::vector< std::vector<double> > > flookupTable;
  int nZ,dZ;
  int nTheta,dTheta;

  TClonesArray* fClusterArray;
  TString fClusterBranchName;

  virtual ~TpcClusterErrorScaler();  
  TpcClusterErrorScaler();
  void readParFile();
  void readLookupFile();
  
  Double_t scaleFunc_phys   (Double_t driftl, Int_t i);
  Double_t scaleFunc_pol    (Double_t driftl, Int_t i);
  Double_t scaleFunc_diff   (Double_t driftl, Int_t i, TVector3 ax);
  Double_t scaleFunc_mean   (Double_t driftl, Int_t opt);
  Double_t scaleFunc_lookup (Double_t driftl, Double_t theta, Int_t i);
  Double_t scaleFunc_lookup2(Double_t driftl, Double_t theta, Int_t i);


  Double_t linInterpol(double v1,double v2,double x);
  Double_t biLinInterpol(double v1, double v2, double v3, double v4, double x, double y);



 public:
  void init(TString filename,Int_t func);
  void postinit();
  void SetParFile(TString f){fparFile=f;}
  void SetNoFairRun(Bool_t noFair=kTRUE){fNoFairRun=noFair;}
  void SetManualDiffusion(double dl, double dt){fDl=dl;fDt=dt;fManualDiffusion=kTRUE;}

  bool IsInit() {return fInit;}
  Int_t getScaleMode(){return ffunc;}
  void getPRFWeight(double& wx, double& wy){postinit();wx=fwx;wy=fwy;}
  double getZJitter(){postinit();return fdz;}

  Double_t scaleFunc(Double_t driftl, Int_t i, TVector3 ax, TVector3 diff);
  Double_t scaleFunc(Double_t driftl, Double_t theta, Int_t i, bool mirr=false);


  void paramVariance(double& varX, double& varY, double& varZ, Double_t driftl, TVector3 dir);
  double paramAmp();
  void setParameter(int parset, int parnum, double par);
  double getParameter(int parset, int parnum){return fpar[parset][parnum];}
  double getFitFractionPar(){return fpar[2][7];}
  TMatrixDSym getErrorCov(double driftl, TVector3 dir);
  TMatrixDSym getErrorCovOnPlane(int clusterID, TMatrixD plane, TMatrixDSym _shapeCovPlane=TMatrixDSym(2), double _zpos=-1);

  inline static TpcClusterErrorScaler* getInstance() {
    if(finstance == NULL) {
      finstance = new TpcClusterErrorScaler();
    }
    return finstance;
  }
  
};
#endif