#ifndef EC_H
#define EC_H

#include <unistd.h>
#include <vector>
#include <iomanip>
#include <fstream>
#include <sstream>
#include <algorithm>
#include <math.h>
#include <string>

#include "TMinuit.h"
#include "TF1.h"
#include "IoSd.h"
#include "IoAuger.h"
#include "LatLong-UTMConversion.h"
#include "EcConstants.h"
#include <map>
#include <iostream>

using namespace std;


namespace kEc {
  void Debug(int);
  void StrictCalib(int);
  
  extern int gDebug;
}

// to switch to old calib
int Ec_UseOldCalib(int i=-1);

//Fit function for parabolic decay time
int fit_exp_parab(double *,int ,double * ,double *);


/// Describes a bit of FADC trace with signal
class TSignal {
public:
  int fGain;            ///< If was from low or high gain
  int fStart;           ///< In time bin
  int fEnd;             ///< In time bin
  double StartTime();   ///< StartTime with respect of EOT
  double fSignal;
  double fPeak;
  double fNbBin;
  int fT10;
  //double fNoise;

  //bool operator <(const TSignal& s) const;
};


/// Basic property of a trace, for PMTs or for Station
class TTraceData {
public:
  float fSigInVEM;            ///< Integrated signal in VEM
  float fPeakInVEM;           ///< Peak signal in VEM
  mutable float fT90;
  mutable float fT70;
  mutable float fT50;
  mutable float fT10;
  int fStartIntegration;
  int fEndIntegration;
  TTraceData() : fSigInVEM(-1), fPeakInVEM(-1), fT90(-1), fT70(-1), fT50(-1),fT10(-1.),fStartIntegration(0),fEndIntegration(0){}
  float fTrace[kIoSd::MAXSAMPLE]; /// fadc trace of each PMT or sum if taken from TCalibStation
  float Tx(float) const;
  float T90() const;
  float T70() const;
  float T50() const;
  float T10() const;
  float SigInVEM();
  float PeakInVEM();
  
  int StartIntegration() {return fStartIntegration;}
  int EndIntegration() {return fEndIntegration;}

  // new calibration
  vector<TSignal> fRawSignals;
  vector<TSignal> fSignals;
};

/// Calibrated PMT trace
class TCalibSdPmt:public TTraceData {
public:

  TCalibSdPmt();
  int fHighGainSat;           ///< Number of channels of High gain saturated (0 = no saturation)
  int fLowGainSat;            ///< Number of channels of Low gain saturated

  int fCalibratedState;       ///< How well is the pmt calibrated: 0=no calib, 1=rough, 2=online, 3=histo, 4=extra water info
  
  float fVemPeak;             ///< Peak VEM values, from online calibration
  float fVemCharge;           ///< Charge VEM values, from online calibration or histograms (best available)
  float fVemChargePB;           ///< Charge VEM values, as fitted by Pierre Billoir's code
  float fBaseline[2];         ///< Baselines, from offline calibration
  float fPmtBaselineSigma;    ///< Baselines sigma, from offline calibration, for the PMT trace
  float fDA;                  ///< Dynode/Anode ratio, offline one
  float fPulseSlope;
  float fHistSlope;
  float fChi;

  // New calibration
  float fADCBaseline[2][kIoSd::MAXSAMPLE];  ///< Baselines, from offline calibration, with undershoot correction
  float fBaselineSigma[2];    ///< Baselines sigma, from offline calibration, for the fADC

  float VemPeak() {
    if (fVemPeak < 0 && kEc::gDebug)cerr <<"VemPeak is not calculated, please call TCalibStation::Signal()."<< endl;
    return fVemPeak;
  }
  float VemCharge() {
    if (fVemCharge < 0 && kEc::gDebug)cerr <<"VemCharge is not calculated, please call TCalibStation::Signal()."<< endl;
    return fVemCharge;
  }
  float PulseSlope() {
    //    if (!fPulseSlopeComputed) cerr<<"PulseSlope is not calculated"<<endl;
    return fPulseSlope;
  }
  float HistSlope() {
    //    if (!fHistSlopeComputed) cerr<<"HistSlope is not calculated"<<endl;
    return fHistSlope;
  }
};

/// Calibrated station (has x,y,z, and Pmt traces in VEM)
class TCalibStation:public IoSdStation, public TTraceData {
public:

  TCalibStation();
  TCalibStation(const IoSdStation &);

  TCalibSdPmt fPmt[kIoSd::NPMT];

  unsigned int fRefSecond;    ///< reference second for StartTime() like functions
  double fX;                  ///< x cartesian coordinate centered on last Northing/Easting/Altitude reference
  double fY;                  ///< y cartesian coordinate
  double fZ;                  ///< z cartesian coordinate
  double fLatitude;           ///< latitude of the Tank (using WGS-84)
  double fLongitude;          ///< longitude of the Tank (using WGS-84)

  int fPmTStatus;             ///< is equal to the number of PMTs working well in the tank
  int fSat;                   ///< 1 or more if station saturated (with low gain)
  int fPLDOffset;
  int fSeveralSignals;
  int fIsInSeed;
  unsigned int fStartBin;
   bool IsLowGainSaturated()  const {
    return (fPmt[0].fLowGainSat  || fPmt[1].fLowGainSat  || fPmt[2].fLowGainSat);
  }
  bool IsHighGainSaturated() const {
    return (fPmt[0].fHighGainSat || fPmt[1].fHighGainSat || fPmt[2].fHighGainSat);
  }

  int Anomaly();
  bool IsLightningEvent();    //should be in Es later
  bool OldIsLightningEvent();    //should be in Es later
  double ShapeParameter(double dist);
  void OldChargeHistogramFit(int,int);
  void ChargeHistogramFit(int,int);
  double DoTheDip(int);//NEW VEM calibration , GAP 2015-050
  int fit_parab(double *, int , int , double& , double& , double& , double& );
  int compute_dip_hump(int, int, double *, double , double &, double &, double &, double& , double &, double &, double& , double&);

  void WaterEstimatesFit(int);
  void OldCalibrate(int, int);
  void CorrectPosition();
  bool CheckPLDOffset();
  void PLDOffset();
  void Calibrate();
  void EstimateCalibration(int);
  void PmtSignal(int);
  void OldSignal(int loadCalib = 1);
  void Signal();
  int CleanFilter(int &, int &, double*);
  void CorrectTraceTime();
  double fDist;               ///< distance to the axis if defined
  double fDcore;              ///< distance to the core if defined
  double fCosZeta;            ///< angle to shower axis
  double fDt;                 ///< dt station to shower front if defined
  double fdt;                 ///< time interval between trigger times
  double fWgt;                ///< Weight of station in time fit. Normally 1.0
  double DecayTime(int); /// returns fitted decay time (by Pierre Billoir)
  double DecayTime(int pmt,int which); /// returns fitted decay time (by Pierre Billoir) which=0 for linear fit 1 for parabolic DT, 2 for parabolic curvature
  bool LoadCalib();
  void CalibrateASCII();       ///< ASCII scintillator calibration (in PMT 1/2)

  float GetBaseline(int, int, int, int);
  float SigOverPeak() {
    return SigInVEM() / PeakInVEM();
  }
  void FindOtherSignal(int, int, float *);
  double StartTime() const {
  
    return 1e-9 * (Gps->NanoSecond -(int)(kEc::FADCDepth - fStartBin) * kEc::TimeSlot + fPLDOffset) + (int)(Gps->Second - fRefSecond);
  }
      
  TCalibSdPmt & operator[](int i) {return fPmt[i];}
  void Print(ostream&);

  // new offline calibration, in EcOff.cc
  int ComputeBaselines(int pmt,int gain);
  
  int ComputeBaselines();
  int BuildSignals(int pmt);
  int BuildSignals();
  int MergeSignals();
  int SelectSignal();
  int TreatSaturated();
};

inline ostream& operator<<(ostream& s, TCalibStation& stat) {
  stat.Print(s);
  return s;
}

/// Calibrated event. Has a list of calibrated stations and a reference Northing, Easting for coordinates transformation
class TEcEvent:public IoSdEvent {
public:
  TEcEvent();
  TEcEvent(EventPos);
  TEcEvent(EventId);
  TEcEvent(const IoSdEvent & sdevent);
  ~TEcEvent();

  vector < TCalibStation > fCalibStations;    ///< stations
  unsigned int fNCalibStations;
  double fRefNorthing;
  double fRefEasting;
  unsigned int RefSecond() {
    return Trigger.Second;
  }
  virtual bool EnoughTanks();
  void Signals();
  void Init();
  void Print();
  void SaveCalib();
  void SetReference(double northing, double easting);
  void ReadPldVersions();
};

#endif