MinuitMigrad.cxx

Go to the documentation of this file.

#include "MinuitMigrad.h"

#include "StatedFCN.h"

#include "FitterException.h"

#include <Minuit/FunctionMinimum.h>
#include <Minuit/MinuitParameter.h>
#include <Minuit/MnMigrad.h>
#include <Minuit/MnUserParameters.h>

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

MinuitMigrad::
MinuitMigrad ( const char * name )
  : Fitter ( name ),
    m_minimizer ( 0 ) // null pointer
{
}

MinuitMigrad::
MinuitMigrad ( const MinuitMigrad & mm )
  : Fitter ( mm ),
    m_minimizer ( 0 )
{
  if ( mm.m_minimizer != 0 ) {
    m_minimizer = new MnMigrad ( * mm.m_minimizer );
  }
}

Fitter *
MinuitMigrad::
clone ( ) const
{
  return new MinuitMigrad ( *this );
}

bool
MinuitMigrad::
needsDerivatives () const
{
  return false;
}

void
MinuitMigrad::
initialize ()
{
  m_minimizer = 0;
  bool yes = m_fcn -> hasFunction ();
  if ( yes ) {
    const vector < string > & names = m_fcn -> getParmNames ();
    unsigned int size = names.size ();
    const vector < double > & parms = m_fcn -> getParameters ();
    MnUserParameters mn_parms;

    for ( unsigned int i = 0; i < size; i++ ) {
      const string name = names [ i ];
      double value = parms [ i ];
      mn_parms.add ( name.c_str(), value, 0.1 );
    }
    m_minimizer = new MnMigrad ( *m_fcn, mn_parms );
  }
}

void
MinuitMigrad::
checkIndex ( unsigned int i )
{
  if ( m_minimizer == 0 ) initialize ();

  if ( m_minimizer == 0 ) {
    string what ( m_name );
    what += ": model function no yet set";
    throw FitterException ( what );
  }
  const vector < double > & parms = m_fcn -> getParameters ();
  unsigned int size = parms.size();
  if ( i < size == false ) {
    string what ( m_name );
    what += ": index to parameter out of range";
    throw FitterException ( what );
  }
}

void
MinuitMigrad::
setLimits ( unsigned int i, double lower, double upper )
{
  checkIndex ( i );
  m_minimizer -> setLimits ( i, lower, upper );
}

void
MinuitMigrad::
setStepSize ( unsigned int i, double size )
{
  checkIndex ( i );
  m_minimizer -> setError ( i, size );
}

bool
MinuitMigrad::
calcBestFit ()
{
  if ( m_minimizer == 0 ) initialize ();

  const vector < double > & init_parms = m_fcn -> getParameters ();
  const vector < int > & fixes = m_fcn -> getFixedFlags ();
  const vector < MinuitParameter> & mparms 
    = m_minimizer -> minuitParameters ();
  unsigned int size = fixes.size();

  for ( unsigned int i = 0; i < size; i++ ) {
    m_minimizer -> setValue ( i, init_parms[i] );
    bool yes = fixes [i] != 0;
    const MinuitParameter & mparm = mparms [i];
    if ( yes && mparm.isFixed() == false ) {
      m_minimizer -> fix ( i );
    }
    else if ( mparm.isFixed() == true ) {
      m_minimizer -> release ( i );
    }
  }
  
  FunctionMinimum fun_min = m_minimizer -> operator () ();

  // in the above call, hidden parameters are maxfcn=0, edm=1.e-5

  std::vector < double > cur_parms = m_minimizer -> params();
  m_fcn -> setParameters ( cur_parms );

  return fun_min.isValid ();
}

int 
MinuitMigrad::
calcCovariance ( std::vector < std::vector < double > >& covar )
{ 
//   assert ( m_minimizer != 0 );
  if ( m_minimizer == 0 ) initialize ();

  const MnUserCovariance & covar_m = m_minimizer ->covariance ();

  unsigned int size = covar_m.nrow ();
  covar.resize ( size );
  for ( unsigned int i = 0; i < size; i++ ) {
    covar[i].resize ( size, 0.0 );
  }

  for ( unsigned int i = 0; i < size; i++ ) {
    for ( unsigned int j = 0; j < size; j++ ) {
      covar[i][j] = covar_m ( i, j );
    }
  }

  return EXIT_SUCCESS;
}

---