/**
   \file


   \author Maximo Ave
   \author Javier Gonzalez
   \version $Id$
   \date 09 Mar 2014
*/

#ifndef __DetectorResponse_h_
#define __DetectorResponse_h_

static const char CVSId__DetectorResponse[] =
  "$Id$";

#include <vector>
#include <map>
#include <cmath>
#include <string>

#define NTYPES 6
#define NCDFMAX 300
#ifndef UNDEF
#define UNDEF -100
#endif
#define DEBUG_RESP 0

// NOTE : no distinction between muon decay in flight and at rest.


namespace TabulatedTankSimulatorNS {

  //
  const double UnitPerPE_wcd = 1. / 244.3;
  const double UnitPerPE_scin = 1. / 15.65;
  const double UnitPerPE_md = 1.;

  // Enlarge Area factor  ( WCD -> radius*f ,  ASCII -> lengh*f, widht*f )
  const double fEnlargeArea = 2.;

  // Auger Water Cherenkov
  const double tankh = 120.;
  const double tankr = 180.;

  // Scintillator 2 m^2
  const double scinlength = 180.;
  const double scinwidth = 108.;
  const double scinheight = 1.;

  // Muon detector:  64 strips in 2 sides with 5.12 m^2 each
  const double nMD = 2.;
  const double MDlength = 400.;
  const double MDwidth = 128.;
  const double MDheight = 1.;
  const double SatLevelMD = 20.; // In a 25 [ns] window for a 64 strips detector

  const double MDBoundary = 50.; // In G4 particles are injected in a larger area to account for MS
  const double MDAreaRatio = (MDwidth + 2 * MDBoundary)* (MDlength + 2 * MDBoundary) / (MDlength* MDwidth);

  /*CDF for one Specie and (to,ke)*/
  class PeCDF {
  public:
    int itype;  //  [0] Muon [1] Electron [2] Gamma [3] MuonLow [4] MichelElectron [5] MichelElectronLow
    double lke, theta, theta_deg; //--KE in GeV, theta in rad, theta in deg
    double mean, rms, median;
    int nentries, nentriesCut;
    int PeCut;
    double probZero, probCut, probDecay, probHit;
    std::vector <std::pair<int, double> > prob;         //--  Differential Probabilities (removed after CDFs are created)
    std::vector <std::pair<int, double> > cdf;          //-- Cummulative Probability Distribution
  };


  /*Tank Response Class*/

  class DetectorResponse {

  private:
    // WCD [VEM] / ASCII [VEM 1 GeV] / AMIGA [muon reaching bottom part of the detector]
    double UnitPerPE;
    int NPECUT;

    std::vector<PeCDF> fMuonCDFs, fElectronCDFs, fGammaCDFs, fMuonLowCDFs , fMichelElectronCDFs,  fMichelElectronLowCDFs ;

    std::string fDataDir;
    int fDetectorType;  // [0] Auger Water Cherenkov [1]  ASCII [2] AMIGA
    int iGrid;          // [0] Auger Water Cherenkov [1]  ASCII [2] AMIGA
    bool fUseEnlargeArea;

  public:
    DetectorResponse(bool flag = false, int DetectorType = 0 , bool UseEnlargeArea = false);
    ~DetectorResponse() {}

    double GetUnitPerPE()
    {
      return UnitPerPE;
    }
    int GetiGrid()
    {
      return iGrid;
    }
    double Interpolate(double* x, double* y, double f);
    double GetSignalBin(double f, PeCDF& fPeCDF);

    bool IsHit(double ke, double theta, int itype , double f);
    double GetPe(double ke, double theta, int itype , double f); // f= -100 Mean Interpolation | f=(0,1) CDF Interpolation

    int GetBin(double lkeo, double to, int itype, int* ibins, bool& flag);
    PeCDF GetCDF(int itype, int ibin);
    const char* GetParticleName(int itype);

    bool ReadG4File(PeCDF& pCDF, PeCDF& pCDF_dec);
    void SetUpCDFs(PeCDF* pCDF);

    bool ReadCDFFile(std::vector<PeCDF>* fCDFs, int size);
    bool WriteCDFFile(std::vector<PeCDF>* fCDFs);

    double GetAreaProj(int itype, double theta, double ke);      // Projected area on ground
    double GetAreaPerp(int itype, double theta, double ke);      // Perpendicular Area
    double GetRelativeTrack(double theta);            // Relative track of a throughgoing particle

  };


} //


#endif // __DetectorResponse_h_

// Configure (x)emacs for this file ...
// Local Variables:
// mode:c++
// compile-command: "make -C .. -k"
// End: