#include "OfflineMon/TPhysMon.hh"

ClassImp(TPhysMon)


//_____________________________________________________________________________
double f_ep(double* X, double* Par) {
  // straw-man parametrization of E/P distribution - function taken from
  // the real histogram

  double x[] = {
    0.025,0.075, 	0.125, 	0.175, 	0.225, 0.275, 	0.325, 	0.375, 	0.425, 	0.475,
    0.525,0.575, 	0.625, 	0.675, 	0.725, 0.775,	0.825, 	0.875, 	0.925, 	0.975,
    1.025,1.075, 	1.125, 	1.175, 	1.225, 1.275,	1.325, 	1.375, 	1.425, 	1.475,
    1.525,1.575, 	1.625, 	1.675, 	1.725, 1.775,	1.825, 	1.875, 	1.925, 	1.975,
    2.025,2.075, 	2.125, 	2.175
  };

  double f[] = {  
    0.000,   0.000,    0.000,    0.000,    0.000,
    0.000,   0.000,    0.000,    0.000,    0.000,
    3.000,   4.000,    6.000,    9.000,   15.000,
    40.000, 102.000,  305.000,  862.000, 1391.000,
    1212.000, 805.000,  494.000,  364.000,  283.000,
    211.000, 188.000,  170.000,  136.000,  134.000,
    114.000,  85.000,   82.000,   78.000,   73.000,
    63.000,  62.000,   62.000,   48.000,   44.000,
    10.000,   0.000,    0.000,    0.000
  };

  // 2 parameters: Par[0] : overall scale, Par[1]: offset in X
  // Assumes the same binning!
 
  double ftouse = 0;
  double bin    = 0.05;
  double xx     = X[0]+Par[1];
  
  for (int i=0; i<43; i++) {
    if ( (x[i] <= xx) && (x[i+1] >= xx)) {
      ftouse = f[i]+(xx-x[i])/bin*(f[i+1]-f[i]);
      break;
    }
  }

  double dxx = 0;
  //  if (xx > 1) dxx = 5*(xx-1);
  if (xx > Par[3]) dxx = 7*(xx-Par[3]);

  return Par[0]*(ftouse + Par[2]*(1-TMath::Exp(-dxx)));

}

//_____________________________________________________________________________
TPhysMon::TPhysMon() {
  fFuncEp = new TF1("func_ep",f_ep,0.85,1.9,4);
  fFuncEp->SetParNames("scale","offset","constant","thr");
  Clear();
}

//_____________________________________________________________________________
TPhysMon::~TPhysMon() {
  delete fFuncEp;
}

//_____________________________________________________________________________
void TPhysMon::Clear(const char* Option) {


  fEOverP        = -1.;
  fSigEOverP     = -1;

  fEoverPincMean = -1.;
  fEoverPincN    = -1.;
  fEoverPincErr  = -1.;

  fEoverPNEMean = -1.;
  fEoverPNEN    = -1.;
  fEoverPNEErr  = -1.;

  fEoverPNWMean = -1.;
  fEoverPNWN    = -1.;
  fEoverPNWErr  = -1.;
}