#include <RAT/ModeCut.hh>
#include <RAT/FitterPMT.hh>
#include <RAT/Processor.hh>
#include <RAT/Log.hh>
#include <RAT/PhysicsUtil.hh>

using namespace RAT;
using namespace RAT::PMTSelectors;

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

using std::vector;
using std::string;

void ModeCut::Initialise( const std::string& ) {
  // Default values
  fLowCut  = -50.0;
  fHighCut =  50.0;
  // -- This time domain for the in-time window was discussed in the reconstruction call
  // of 17/05/2017 and agreed in these values.
  // The values can however be changed through macro
  fLowDomainLimit = -50;
  fHighDomainLimit = 600;
  detail << "ModeCut::Initialise: Initializing mode cut. Default domain [" << fLowDomainLimit << "," << fHighDomainLimit << "]" << newline;
}

void ModeCut::SetD( const std::string& param, double value ) {
  if( param == string( "lowLimit" ) ) {
    fLowCut = value;
  }else if( param == string( "highLimit" ) ) {
    fHighCut = value;
  }else{
    throw Processor::ParamUnknown( param );
  }
}
void ModeCut::SetI( const std::string& param, int value) {
  if( param == string( "lowLimit" ) ) {
    fLowCut =  static_cast<double>(value);
  }else if( param == string( "highLimit" ) ) {
    fHighCut =  static_cast<double>(value);
  } else if (param == string("low_domain_limit")) {
    fLowDomainLimit = value;
  } else if (param == string("high_domain_limit")) {
    fHighDomainLimit = value;
  } else {
    throw Processor::ParamUnknown( param );
  }
}


vector<FitterPMT> ModeCut::GetSelectedPMTs( const std::vector<FitterPMT>& data,
    const DS::FitVertex& ) {

  // First find the mode
  double modeTime = GetModeTime( data );

  // Now select the PMTs
  vector<FitterPMT> selectedPMTs;
  for( vector<FitterPMT>::const_iterator iPMT = data.begin(); iPMT != data.end(); ++iPMT ){
    const double hitTime = iPMT->GetTime();
    if( hitTime > (modeTime+fLowCut) && hitTime < (modeTime+fHighCut) ){
      selectedPMTs.push_back( *iPMT );
    }
  }
  return selectedPMTs;
}

double ModeCut::GetModeTime( const std::vector<FitterPMT>& data )  {
  std::vector<unsigned int> timeHist((fHighDomainLimit-fLowDomainLimit),0);
  for( vector<FitterPMT>::const_iterator iPMT = data.begin(); iPMT != data.end(); ++iPMT ) {
    // find the bin to increment
    // -- If the time is outside the window, throw an exception? Discard it?
    double pmt_time = iPMT->GetTime();
    if ((util::round_2_int(pmt_time) < fLowDomainLimit) || (util::round_2_int(pmt_time) >= fHighDomainLimit)) {
      static bool first = true;
      if (first) {
        warn << BRED << "ModeCut::GetModeTime : WARNING::Found a PMT time outside the expected window. This should be investigated." << CLR << newline;
        first = false;
      }
      continue; // bypass. Don't count to calculate the mode
    }
    // Apply an offset the size of the lower domain limit
    size_t bin = util::round_2_int(pmt_time) - fLowDomainLimit;
    timeHist.at(bin)++;
  }
  // Now find the mode
  return std::distance(timeHist.begin(),std::max_element(timeHist.begin(), timeHist.end())) + static_cast<double>(fLowDomainLimit);

}