NonLinearLSQ.h

#ifndef NonLinearLSQ_h
#define NonLinearLSQ_h
 
/* SPDX-License-Identifier: CC0-1.0 */                                                                    
 
#include <iostream>
#include <sstream>
#include <string>
#include <vector>
#include <algorithm>
 
#include "tnt_array1d.h"
#include "tnt_array1d_utils.h"
#include "tnt_array2d.h"
#include "tnt_array2d_utils.h"
 
#include <gsl/gsl_vector.h>
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_blas.h>
#include <gsl/gsl_multifit_nlin.h>
 
#include "IException.h"
 
namespace Isis {
 
  class NonLinearLSQ {
    public:
      typedef TNT::Array1D<double> NLVector;
      typedef TNT::Array2D<double> NLMatrix;
 
      //  Constructors and Destructor
      NonLinearLSQ() : _fitParms(), _uncert(), _nIters(0), _maxIters(50),
                       _status(0), _userTerminated(false), _userMessage() {}
 
      virtual ~NonLinearLSQ() { }
 
      virtual int nSize() const = 0;
      virtual int nParms() const = 0;
 
      void setMaxIters(int m) { _maxIters = m; }
 
      int maxIters() const { return (_maxIters); }
 
      virtual NLVector guess()  = 0;
      virtual NLVector f_x(const NLVector &x) = 0;
      virtual NLMatrix df_x(const NLVector &x) = 0;
 
      virtual double absErr() const { return (1.0E-4); }
      virtual double relErr() const { return (1.0E-4); }
 
      int curvefit();
      inline int status() const { return (_status); }
      inline bool success() const { return (_status == GSL_SUCCESS); }
      inline bool success(int status) const { return (status == GSL_SUCCESS); }
      inline std::string statusstr() const {
        return (std::string(gsl_strerror(_status)));
      }
      inline std::string statusstr(int status) const {
        return (std::string(gsl_strerror(status)));
      }
 
      virtual int checkIteration(const int Iter, const NLVector &fitcoefs,
                                 const NLVector &uncerts, double cplxconj,
                                 int Istatus) {
        return (Istatus);
      }
 
      inline NLVector coefs() const { return (_fitParms); }
      inline NLVector uncert() const { return (_uncert); }
      inline int nIterations() const { return (_nIters); }
 
    protected:
      void Terminate(const std::string &message = "");
      void Abort(const std::string &reason = "");
 
      bool doContinue() const { return (!_userTerminated); }
 
 
    private:
      NLVector    _fitParms;
      NLVector    _uncert;
      int         _nIters;
      int         _maxIters;
      int         _status;
      bool        _userTerminated;
      std::string _userMessage;
 
      struct _nlsqPointer {
          NonLinearLSQ *nlsq;
      };
      static int f(const gsl_vector *x, void *params, gsl_vector *fx);
      static int df(const gsl_vector *x, void *params, gsl_matrix *J);
      static int fdf(const gsl_vector *x, void *params, gsl_vector *fx,
                     gsl_matrix *J);
 
      NLVector gslToNlsq(const gsl_vector *v) const;
      NLMatrix gslToNlsq(const gsl_matrix *m) const;
      gsl_vector *NlsqTogsl(const NLVector &v, gsl_vector *gv = 0) const;
      gsl_matrix *NlsqTogsl(const NLMatrix &m, gsl_matrix *gm = 0) const;
      NLVector getUncertainty(const gsl_matrix *m) const;
  };
 
}     // namespace Isis
#endif