#include <RAT/Log.hh>
#include <RAT/Minuit.hh>
#include <RAT/OptimisedComponent.hh>
#include <RAT/DS/FitResult.hh>
#include <RAT/Processor.hh>
using namespace RAT::Optimisers;

#include <Minuit2/MnMigrad.h>
#include <Minuit2/MnMinimize.h>
#include <Minuit2/MnSimplex.h>
#include <Minuit2/MnUserParameters.h>
#include <Minuit2/FunctionMinimum.h>
using namespace ROOT::Minuit2;

#include <TError.h>

#include <vector>
using namespace std;

void
Minuit::Initialise( const std::string& param )
{
  // Default values (can be changed by macro file inputs):
  fTol = 0.1;
  fMaxCalls = 0;
  fMethod = "Migrad";

  // Limit ROOT Minuit2 verbosity (global variable)
  gErrorIgnoreLevel = 1001;
}

void
Minuit::BeginOfRun( DS::Run& )
{
}

void
Minuit::SetD( const std::string& param,
              double value )
{
  if( param == string( "fitTolerance" ) )
    fTol = value;
  else
    throw Processor::ParamUnknown( param );
}

void
Minuit::SetI( const std::string& param,
              int value )
{
  if( param == string( "maxCalls" ) )
    fMaxCalls = value;
  else
    throw Processor::ParamUnknown( param );
}

void
Minuit::SetS( const string& param,
              const string& value )
{
  if( param == string( "method" ) )
    fMethod = value;
  else
    throw Processor::ParamUnknown( param );

  if(fMethod!="Migrad" && fMethod!="Minimize" && fMethod!="Simplex") {
    throw Processor::ParamUnknown( param+" not recognized: "+fMethod );
    fMethod = "Migrad";
  }
}

double
Minuit::Minimise()
{
  fMinFactor = 1.0;
  fUp = 0.5;
  return fMinFactor * MinuitOptimise(); // Ensure Minimum is returned
}

double
Minuit::Maximise()
{
  fMinFactor = -1.0;
  fUp = 0.5;
  return fMinFactor * MinuitOptimise(); // Ensure Maximum is returned
}

double
Minuit::MinuitOptimise()
{
  // First create the params
  vector<double> params = fComponent->GetParams();
  vector<double> posErrors = fComponent->GetPositiveErrors();
  vector<double> negErrors = fComponent->GetNegativeErrors();

  static bool warned = false;
  if( !warned && !equal( posErrors.begin(), posErrors.end(), negErrors.begin() ) )
    {
      warned = true;
      warn << "Minuit::MinuitOptimise: Negative errors differ, but negative errors are ignored." << newline;
    }
  MnUserParameters mnParams = MnUserParameters( params, posErrors );

  // Following needed because of scope problems.  Have to create theMin before
  // the if{...} logic; otherwise theMin goes out of scope when leaving the
  // if{...} logic.
  MnUserTransformation trafo;
  MinimumState minState(params.size());
  MinimumSeed minSeed(minState, trafo);
  FunctionMinimum theMin(minSeed, 1.);

  if(fMethod=="Migrad") {
      MnMigrad migrad( *this, mnParams );
      theMin = migrad(fMaxCalls,fTol);
  }else if(fMethod=="Minimize") {
      MnMinimize minimize( *this, mnParams );
      theMin = minimize(fMaxCalls,fTol);
  }else{
      MnSimplex simplex ( *this, mnParams );
      theMin = simplex(fMaxCalls,fTol);
  }

  debug << "Minuit::MinuitOptimise: -- IsValid, Fval, Edm, NFcn = "<< theMin.IsValid() <<" "<< theMin.Fval()<<" "<< theMin.Edm()<<" "<< theMin.NFcn() << endl;

  // Now update the FitResult with answer
  MnUserParameters result = theMin.UserParameters();
  fParams = result.Params();
  fNegativeErrors = fPositiveErrors = result.Errors();
  fValid = theMin.IsValid();

  fComponent->SetFOM( "EDM", theMin.Edm() );
  fComponent->SetFOM( "NFnCalls", (double)theMin.NFcn() );

  return theMin.Fval();
}

double
Minuit::operator()( const std::vector<double>& lParams ) const
{
  double returnVal = fMinFactor * (*fComponent)( lParams ); // Call the component operator()
  return returnVal;
}