FitObject.cpp

#include "FitObject.h"
#include "Sct/OutOfRangeError.h"

#include <TF1.h>
#include <iomanip>
#include <cfloat>
#include <algorithm>

using namespace std;

using namespace Sct;

namespace SctData {

FitObject::FitObject(string formula, unsigned int nPar, const vector<string>& parNames, unsigned int nDim) : 
    m_parameter(nPar, 0), m_parError(nPar, 0), m_parName(parNames), m_parFixed(nPar, false), m_varMax(nDim, DBL_MAX), 
    m_varMin(nDim, -DBL_MAX), m_chiSquared(0), m_nDF(0), m_formula(formula) {
    
}
    
// Create a FitObject from a root TF1 object
// uses FitObject::operator= from below
FitObject::FitObject(string formula, unsigned int nPar, const vector<string>& parNames, const TF1& rootfunc) : 
    m_parameter(nPar, 0), m_parError(nPar, 0), m_parName(parNames), m_parFixed(nPar, false), m_varMax(1, DBL_MAX), 
    m_varMin(1, -DBL_MAX), m_chiSquared(0), m_nDF(0), m_formula(formula) {
    (*this) = rootfunc;
}

void FitObject::reset() throw() {
    //not much to be done here
    m_nDF = 0;
    m_chiSquared = 0;
    fill(m_parameter.begin(), m_parameter.end(), 0);
    fill(m_parError.begin(), m_parError.end(), 0);
    fill(m_parFixed.begin(), m_parFixed.end(), false);
    fill(m_varMax.begin(), m_varMax.end(), DBL_MAX);
    fill(m_varMin.begin(), m_varMin.end(), -DBL_MAX);
}

char* FitObject::getParName(const int ipar) const throw(LogicError) {
#ifndef NDEBUG
    if ( ipar<0 || ipar>=this->getNPar() ) {
        throw OutOfRangeError<int>("FitObject::getParName ", __FILE__, __LINE__, ipar,0,this->getNPar()-1);
    }
#endif
    return const_cast <char*> ( (m_parName[ipar]).c_str() ) ;
}

int FitObject::getParIndex(const string& name) const throw(LogicError) {
    for (unsigned int i=0; i<m_parName.size(); ++i) {
        if (m_parName[i] == name)
            return i;
    }
    throw InvalidArgumentError("FitObject::getParIndex Didn't find " + name, __FILE__, __LINE__);
}


double FitObject::getParameter(const int ipar) const throw(LogicError) {
#ifndef NDEBUG
    if ( ipar<0 || ipar>=this->getNPar() ) {
        throw OutOfRangeError<int>("FitObject::getParError ",__FILE__, __LINE__, ipar,0,this->getNPar()-1);
    }
#endif
    return m_parameter[ipar];
}

double FitObject::getParError(const int ipar) const throw(LogicError) {
#ifndef NDEBUG
    if ( ipar<0 || ipar>=this->getNPar() ) {
        throw OutOfRangeError<int>("FitObject::getParError ", __FILE__, __LINE__, ipar,0,this->getNPar()-1);
    }
#endif
    return m_parError[ipar];
}

void FitObject::setParameter(const int ipar, const double value)
throw(LogicError) {
#ifndef NDEBUG
    if ( ipar<0 || ipar>=this->getNPar() ) {
        throw OutOfRangeError<int>("FitObject::setParameter", __FILE__, __LINE__, ipar,0,this->getNPar()-1);
    }
#endif
    m_parameter[ipar]=value;
}

void FitObject::setParError(const int ipar, const double value)
throw(LogicError) {
#ifndef NDEBUG
    if ( ipar<0 || ipar>=this->getNPar() ) {
        throw OutOfRangeError<int>("FitObject::setParError", __FILE__, __LINE__, ipar,0,this->getNPar()-1);
    }
#endif
    m_parError[ipar]=value;
}

void FitObject::setParName(int ipar, const char* name) throw(LogicError) {
#ifndef NDEBUG
    if ( ipar<0 || ipar>=this->getNPar() ) {
        throw OutOfRangeError<int>("FitObject::setParName", __FILE__, __LINE__, ipar,0,this->getNPar()-1);
    }
#endif
    m_parName[ipar]=name;
}

double FitObject::getVarMax(const int ivar) const throw (LogicError) {
#ifndef NDEBUG
    if ( ivar<0||ivar>=this->getNDim() ) {
        throw OutOfRangeError<int>("FitObject::getVarMax", __FILE__, __LINE__, ivar,0,this->getNDim()-1);
    }
#endif
    return m_varMax[ivar];
}

double FitObject::getVarMin(const int ivar) const throw (LogicError) {
#ifndef NDEBUG
    if ( ivar<0||ivar>=this->getNDim() ) {
        throw OutOfRangeError<int>("FitObject::getVarMin", __FILE__, __LINE__, ivar,0,this->getNDim()-1);
    }
#endif
    return m_varMin[ivar];
}

void FitObject::setVarMax(const int ivar, const double value) throw(LogicError) {
#ifndef NDEBUG
    if ( ivar<0||ivar>=this->getNDim() ) {
        throw OutOfRangeError<int>("FitObject::setVarMax", __FILE__, __LINE__, ivar,0,this->getNDim()-1);
    }
#endif
    m_varMax[ivar]=value;
}

void FitObject::setVarMin(const int ivar, const double value) throw(LogicError) {
#ifndef NDEBUG
    if ( ivar<0||ivar>=this->getNDim() ) {
        throw OutOfRangeError<int>("FitObject::setVarMin", __FILE__, __LINE__, ivar,0,this->getNDim()-1);
    }
#endif
    m_varMin[ivar]=value;
}

void FitObject::fixParameter(int ipar, bool fix) {
#ifndef NDEBUG
    if ( ipar<0||ipar>=getNPar() ) {
        throw OutOfRangeError<int>("FitObject::isFixed", __FILE__, __LINE__, ipar,0,getNPar()-1);
    }
#endif    
    m_parFixed[ipar] = fix;
}

bool FitObject::isFixed(int ipar) const {
#ifndef NDEBUG
    if ( ipar<0||ipar>=getNPar() ) {
        throw OutOfRangeError<int>("FitObject::isFixed", __FILE__, __LINE__, ipar,0,getNPar()-1);
    }
#endif    
    return m_parFixed[ipar];
}


auto_ptr<TF1> FitObject::makeRootTF1() const throw (LogicError) {

    // Check this is a one-dimensional function, otherwise
    // it dosen't make sense
#ifndef NDEBUG
    if ( this->getNDim()!=1 ) {
        throw OutOfRangeError<int>("FitObject::makeRootTF1 nDim", __FILE__, __LINE__, this->getNDim(),1,1);
    }
#endif

    auto_ptr<TF1> rootfunc  = this->makeBasicRootTF1();

    for ( int ipar=0 ; ipar<this->getNPar() ; ipar++ ) {
        // set the parameters
        rootfunc->SetParameter(ipar, this->getParameter(ipar) ) ;
        // set the errors
        rootfunc->SetParError(ipar, this->getParError(ipar) );
        // set the names
        rootfunc->SetParName(ipar, this->getParName(ipar) );
    if (isFixed(ipar)) rootfunc->SetParLimits(ipar, this->getParameter(ipar), this->getParameter(ipar));
    }

    // set the number of degrees of freedom and the chisquared
    rootfunc->SetNDF( this->getNDF() );
    rootfunc->SetChisquare( this->getChiSquared() );

    // set the minimum and maximum
    if (this->getNDim()==1) {
        rootfunc->SetRange ( this->getVarMin(0) , this->getVarMax(0) );
    } else if (this->getNDim()==2) {
        rootfunc->SetRange ( this->getVarMin(0) ,this->getVarMin(1) , this->getVarMax(0) , this->getVarMax(1) );
    } else if (this->getNDim()==3) {
        rootfunc->SetRange ( this->getVarMin(0) , this->getVarMin(1) ,this->getVarMin(2) ,
                             this->getVarMax(0) , this->getVarMax(1) ,this->getVarMax(3) );
    } else {
        throw OutOfRangeError<int>("FitObject::makeRootTF1 ", __FILE__, __LINE__, this->getNDim(),1,3);
    }


    //Set title
    rootfunc->SetTitle(this->getFormula());
    rootfunc->SetName(this->getFormula());
    return rootfunc;
}

// make this FitObject have the same parameter as a ROOT TF1 function
FitObject& FitObject::operator=(const TF1& rootfunc) throw(LogicError) {
    // get the number of parameters of the root function
    int npar=rootfunc.GetNpar();
#ifndef NDEBUG
    if ( this->getNPar()!=npar ) {
        throw OutOfRangeError<int>("FitObject::operator=(const TF1& rootfunc) - different numbers of parameters", __FILE__, __LINE__, npar, getNPar(),getNPar());
    }
#endif

    for (int ipar=0 ; ipar<getNPar() ; ipar++ ) {
        this->setParameter(ipar, rootfunc.GetParameter(ipar) );
        this->setParError(ipar, rootfunc.GetParError(ipar) );
        // leave parameter names as in constructor.
    }

    // Set the chisquared and the number of degrees of freedom:
    this->setChiSquared( rootfunc.GetChisquare() );
    this->setNDF( rootfunc.GetNDF() );

    // Set the X max and X min
    this->setVarMax(0, rootfunc.GetXmax() );
    this->setVarMin(0, rootfunc.GetXmin() );

    // Don't set the formula: leave the one created by the constructor!
    // this->setFormula(rootfunc.GetTitle() );
    return (*this);
}
void FitObject::print() const throw(LogicError) {
    cout <<*this;
}

ostream& operator << (ostream & os, const SctData::FitObject& f) throw (LogicError) {
    os <<"========================================"<<endl;
    os << "FitObject, " << f.getFormula() << ", (" << f.getNDim() << " dimensional)" << endl;
    os << "----variables----------range-----------"<<endl;
    for (int idim=0; idim<f.getNDim(); idim++ ) {
        os << " var:"<< idim << " : " << setw(10) << f.getVarMin(idim)
        << " -> " << setw(10) << f.getVarMax(idim)<<endl;
    }
    os << "----parameters--------values----errors--"<<endl;
    for (int ipar=0; ipar<f.getNPar(); ipar++) {
        os << " par:"<<ipar<<" "<< setw(10) << f.getParName(ipar) << " = " << setw(10) << f.getParameter(ipar)
        << " +- " << f.getParError(ipar)<<endl;
    }
    os << " Chisq / NDF    =    " << f.getChiSquared() << " / " << f.getNDF()<<endl;

    return os;
}

} // end of namespace SctData

---