// @(#)root/mathcore:$Id$
// Author: L. Moneta Thu Sep 21 16:21:29 2006

/**********************************************************************
 *                                                                    *
 * Copyright (c) 2006  LCG ROOT Math Team, CERN/PH-SFT                *
 *                                                                    *
 *                                                                    *
 **********************************************************************/

// Implementation file for class FitConfig

#include "Fit/FitConfig.h"

#include "Fit/FitResult.h"

#include "Math/IParamFunction.h"
#include "Math/Util.h"

#include "Math/Minimizer.h"
#include "Math/Factory.h"

#include <cmath> 

#include <string> 
#include <sstream> 

#include "Math/Error.h"

//#define DEBUG
#ifdef DEBUG
#endif
#include <iostream>

namespace ROOT { 

namespace Fit { 



FitConfig::FitConfig(unsigned int npar) : 
   fNormErrors(false),
   fParabErrors(false), // ensure that in any case correct parabolic errors are estimated
   fMinosErrors(false),    // do full Minos error analysis for all parameters
   fUpdateAfterFit(true),    // update after fit
   fWeightCorr(false),
   fSettings(std::vector<ParameterSettings>(npar) )  
{
   // constructor implementation
}


FitConfig::~FitConfig() 
{
   // destructor implementation. No Operations
}

FitConfig::FitConfig(const FitConfig &rhs) { 
   // Implementation of copy constructor
   (*this) = rhs; 
}

FitConfig & FitConfig::operator = (const FitConfig &rhs) { 
   // Implementation of assignment operator.
   if (this == &rhs) return *this;  // time saving self-test

   fNormErrors = rhs.fNormErrors; 
   fParabErrors = rhs.fParabErrors; 
   fMinosErrors = rhs.fMinosErrors; 
   fUpdateAfterFit = rhs.fUpdateAfterFit;
   fWeightCorr     = rhs.fWeightCorr;

   fSettings = rhs.fSettings; 
   fMinosParams = rhs.fMinosParams; 

   fMinimizerOpts = rhs.fMinimizerOpts;

   return *this;
}

void FitConfig::SetFromFitResult(const FitResult &result) { 
   // Implementation of setting of parameters from the result of the fit
   // all the other options will stay the same. 
   // If the size of parameters do not match they will be re-created
   // but in that case the bound on the parameter will be lost

   unsigned int npar = result.NPar();
   if (fSettings.size() !=  npar) {
      fSettings.clear();
      fSettings.resize(npar);
   }
   // fill the parameter settings 
   for (unsigned int i = 0; i < npar; ++i) {
      if (result.IsParameterFixed(i) )
         fSettings[i].Set(result.ParName(i), result.Value(i) ); 
      else { 
         if (result.IsParameterBound(i) && !fSettings[i].IsBound() ) {
            // bound on parameters will be lost- must be done by hand by user 
            std::string msg = "Bound on parameter " + result.ParName(i) + " is lost; it must be set again by the user";
            MATH_WARN_MSG("FitConfig::SetFromResult",msg.c_str() );
         }
         fSettings[i].Set( result.ParName(i), result.Value(i), result.Error(i) ); 

         // query if parameter needs to run Minos
         if (result.HasMinosError(i) ) {
            if (fMinosParams.size() == 0) { 
               fMinosErrors = true; 
               fMinosParams.reserve(npar-i);
            }
            fMinosParams.push_back(i);
         }
      }
   }

   // set information about errors 
   SetNormErrors( result.NormalizedErrors() );

   // set also minimizer type 
   // algorithm is after " / "
   const std::string & minname = result.MinimizerType();
   size_t pos = minname.find(" / ");
   if (pos != std::string::npos) { 
      std::string minimType = minname.substr(0,pos);
      std::string algoType = minname.substr(pos+3,minname.length() );
      SetMinimizer(minimType.c_str(), algoType.c_str() );
   }
   else { 
      SetMinimizer(minname.c_str());
   }      
}


void FitConfig::SetParamsSettings(unsigned int npar, const double *params, const double * vstep ) { 
   // initialize FitConfig from given parameter values and step sizes
   // if npar different than existing one - clear old one and create new ones
   if (params == 0) { 
      fSettings =  std::vector<ParameterSettings>(npar); 
      return; 
   }
   // if a vector of parameters is given and parameters are not existing or are of different size
   bool createNew = false; 
   if (npar != fSettings.size() ) { 
      fSettings.clear(); 
      fSettings.reserve(npar); 
      createNew = true; 
   }
   unsigned int i = 0; 
   const double * end = params+npar;
   for (const double * ipar = params; ipar !=  end; ++ipar) {  
      double val = *ipar;       
      double step = 0; 
      if (vstep == 0) {  
         step = 0.3*std::fabs(val);   // step size is 30% of par value
         //double step = 2.0*std::fabs(val);   // step size is 30% of par value
         if (val ==  0) step  =  0.3; 
      }
      else 
         step = vstep[i]; 

      if (createNew) 
         fSettings.push_back( ParameterSettings("Par_" + ROOT::Math::Util::ToString(i), val, step ) ); 
      else {
         fSettings[i].SetValue(val); 
         fSettings[i].SetStepSize(step); 
      }

      i++;
   }
}

void FitConfig::CreateParamsSettings(const ROOT::Math::IParamMultiFunction & func) { 
   // initialize from model function
   // set the parameters values from the function
   unsigned int npar = func.NPar(); 
   const double * begin = func.Parameters(); 
   if (begin == 0) { 
      fSettings =  std::vector<ParameterSettings>(npar); 
      return; 
   }

   fSettings.clear(); 
   fSettings.reserve(npar); 
   const double * end =  begin+npar; 
   unsigned int i = 0; 
   for (const double * ipar = begin; ipar !=  end; ++ipar) {  
      double val = *ipar; 
      double step = 0.3*std::fabs(val);   // step size is 30% of par value
      //double step = 2.0*std::fabs(val);   // step size is 30% of par value
      if (val ==  0) step  =  0.3; 
      
      fSettings.push_back( ParameterSettings(func.ParameterName(i), val, step ) ); 
#ifdef DEBUG
      std::cout << "FitConfig: add parameter " <<  func.ParameterName(i) << " val = " << val << std::endl;
#endif
      i++;
   } 

}

ROOT::Math::Minimizer * FitConfig::CreateMinimizer() { 
   // create minimizer according to the chosen configuration using the 
   // plug-in manager

   const std::string & minimType = fMinimizerOpts.MinimizerType(); 
   const std::string & algoType  = fMinimizerOpts.MinimizerAlgorithm(); 

   std::string  defaultMinim = ROOT::Math::MinimizerOptions::DefaultMinimizerType(); 
   ROOT::Math::Minimizer * min = ROOT::Math::Factory::CreateMinimizer(minimType, algoType); 
   // check if a different minimizer is used (in case a default value is passed, then set correctly in FitConfig)
   const std::string & minim_newDefault = ROOT::Math::MinimizerOptions::DefaultMinimizerType();
   if (defaultMinim != minim_newDefault )  fMinimizerOpts.SetMinimizerType(minim_newDefault.c_str());
      
   if (min == 0) { 
      // if creation of minimizer failed force the use by default of Minuit
      std::string minim2 = "Minuit"; 
      if (minimType == "Minuit") minim2 = "Minuit2";
      if (minimType != minim2 ) {
         std::string msg = "Could not create the " + minimType + " minimizer. Try using the minimizer " + minim2; 
         MATH_WARN_MSG("FitConfig::CreateMinimizer",msg.c_str());
         min = ROOT::Math::Factory::CreateMinimizer(minim2,"Migrad"); 
         if (min == 0) { 
            MATH_ERROR_MSG("FitConfig::CreateMinimizer","Could not create the Minuit2 minimizer");
            return 0; 
         }
         SetMinimizer( minim2.c_str(),"Migrad"); 
      }
      else {
         std::string msg = "Could not create the Minimizer " + minimType; 
         MATH_ERROR_MSG("FitConfig::CreateMinimizer",msg.c_str());
         return 0;
      }
   } 

   // set default max of function calls according to the number of parameters
   // formula from Minuit2 (adapted)
   if (fMinimizerOpts.MaxFunctionCalls() == 0) {  
      unsigned int npar =  fSettings.size();      
      int maxfcn = 1000 + 100*npar + 5*npar*npar;
      fMinimizerOpts.SetMaxFunctionCalls(maxfcn); 
   }


   // set default minimizer control parameters 
   min->SetPrintLevel( fMinimizerOpts.PrintLevel() ); 
   min->SetMaxFunctionCalls( fMinimizerOpts.MaxFunctionCalls() ); 
   min->SetMaxIterations( fMinimizerOpts.MaxIterations() ); 
   min->SetTolerance( fMinimizerOpts.Tolerance() ); 
   min->SetPrecision( fMinimizerOpts.Precision() ); 
   min->SetValidError( fParabErrors );
   min->SetStrategy( fMinimizerOpts.Strategy() );
   min->SetErrorDef( fMinimizerOpts.ErrorDef() );


   return min; 
} 

void FitConfig::SetDefaultMinimizer(const char * type, const char *algo ) { 
   // set the default minimizer type and algorithms
   ROOT::Math::MinimizerOptions::SetDefaultMinimizer(type, algo); 
} 

void FitConfig::SetMinimizerOptions(const ROOT::Math::MinimizerOptions & minopt) {  
   // set all the minimizer options
   fMinimizerOpts = minopt; 
}


   } // end namespace Fit

} // end namespace ROOT