#include "objfitter.h"

TMatrix<r_4>
ObjectFitter::FitResidus(TMatrix<r_4> const & mtx, GeneralFit& gfit,
			 double xorg, double yorg, double dx, double dy)
{
  return( ArrayFitter<r_4>::FitResidus(mtx, gfit, xorg, yorg, dx, dy) );
}

TMatrix<r_8>
ObjectFitter::FitResidus(TMatrix<r_8> const & mtx, GeneralFit& gfit,
			 double xorg, double yorg, double dx, double dy)
{
  return( ArrayFitter<r_8>::FitResidus(mtx, gfit, xorg, yorg, dx, dy) );
}

TMatrix<r_4>
ObjectFitter::FitFunction(TMatrix<r_4> const & mtx, GeneralFit& gfit,
			 double xorg, double yorg, double dx, double dy)
{
  return( ArrayFitter<r_4>::FitFunction(mtx, gfit, xorg, yorg, dx, dy) );
}

TMatrix<r_8>
ObjectFitter::FitFunction(TMatrix<r_8> const & mtx, GeneralFit& gfit,
			 double xorg, double yorg, double dx, double dy)
{
  return( ArrayFitter<r_8>::FitFunction(mtx, gfit, xorg, yorg, dx, dy) );
}

TVector<r_4>
ObjectFitter::FitResidus(TVector<r_4> const & vec, GeneralFit& gfit,
			       double xorg, double dx)
{  
  return( ArrayFitter<r_4>::FitResidus(vec, gfit, xorg, dx) );
}

TVector<r_8>
ObjectFitter::FitResidus(TVector<r_8> const & vec, GeneralFit& gfit,
			       double xorg, double dx)
{  
  return( ArrayFitter<r_8>::FitResidus(vec, gfit, xorg, dx) );
}

TVector<r_4>
ObjectFitter::FitFunction(TVector<r_4> const & vec, GeneralFit& gfit,
			       double xorg, double dx)
{  
  return( ArrayFitter<r_4>::FitFunction(vec, gfit, xorg, dx) );
}

TVector<r_8>
ObjectFitter::FitFunction(TVector<r_8> const & vec, GeneralFit& gfit,
			       double xorg, double dx)
{  
  return( ArrayFitter<r_8>::FitFunction(vec, gfit, xorg, dx) );
}


// ---------------------------------------------------------------------------

template <class T>
TMatrix<T> 
ArrayFitter<T>::FitResidus(TMatrix<T> const & mtx, GeneralFit& gfit,
			 double xorg=0.,double yorg=0.,double dx=1.,double dy=1.)
// Retourne une classe contenant les residus du fit ``gfit''.
// On suppose que x=j (colonnes) et y=i (lignes) pour m(i,j).
// Les coordonnees de l'element (i,j) sont :
//    (i,j) -> x = xorg+j*dx , y = yorg+i*dy
{
if(mtx.NCols()<=0||mtx.NRows()<=0)
  throw(SzMismatchError("ArrayFitter::FitResidus(TMatrix<T>...) size mismatch\n"));
GeneralFunction* f = gfit.GetFunction();
if(f==NULL)
  throw(NullPtrError("ArrayFitter::FitResidus(TMatrix<T>...) GeneraFit==NULL\n"));
int npar =  gfit.GetNPar();
if(npar==0)
  throw(SzMismatchError("ArrayFitter::FitResidus(TMatrix<T>...) GeneraFit 0 parametre\n"));
double* par = new double[npar];
{for(int i=0;i<npar;i++) par[i] = gfit.GetParm(i);}
TMatrix<T> m(mtx);
m.SetTemp(true);
for(uint_4 i=0;i<mtx.NRows();i++) for(uint_4 j=0;j<mtx.NCols();j++) {
  double x[2] = {xorg+j*dx,yorg+i*dy};
  m(i,j) -= f->Value(x,par);
}
delete [] par;
return m;
}


template <class T>
TMatrix<T> 
ArrayFitter<T>::FitFunction(TMatrix<T> const & mtx, GeneralFit& gfit,
			  double xorg=0.,double yorg=0.,double dx=1.,double dy=1.)

// Retourne une classe contenant la fonction du fit ``gfit''.
// On suppose que x=j (colonnes) et y=i (lignes) pour m(i,j).
// Les coordonnees de l'element (i,j) sont :
//    (i,j) -> x = xorg + j*dx , y = yorg + i*dy

{
if(mtx.NCols()<=0||mtx.NRows()<=0)
  throw(SzMismatchError("ArrayFitter::FitFunction(TMatrix<T>...) size mismatch\n"));
GeneralFunction* f = gfit.GetFunction();
if(f==NULL)
  throw(NullPtrError("ArrayFitter::FitFunction(TMatrix<T>...) GeneraFit==NULL\n"));
int npar =  gfit.GetNPar();
if(npar==0)
  throw(SzMismatchError("ArrayFitter::FitFunction(TMatrix<T>...) GeneraFit 0 parametre\n"));
double* par = new double[npar];
{for(int i=0;i<npar;i++) par[i] = gfit.GetParm(i);}
TMatrix<T> m(mtx.NRows(), mtx.NCols());
m.SetTemp(true);
for(uint_4 i=0;i<mtx.NRows();i++) for(uint_4 j=0;j<mtx.NCols();j++) {
  double x[2] = {xorg+j*dx,yorg+i*dy};
  m(i,j) = f->Value(x,par);
}
delete [] par;
return m;

}

template <class T>
TVector<T> 
ArrayFitter<T>::FitResidus(TVector<T> const & vec, GeneralFit& gfit,
			 double xorg=0.,double dx=1.)
// Retourne une classe contenant les residus du fit ``gfit''.
// La coordonnee de l'element (i) est -> x = xorg + i*dx
{
if(vec.NElts()<=0)
  throw(SzMismatchError("ArrayFitter::FitResidus(TVector<T>...)  size mismatch\n"));
GeneralFunction* f = gfit.GetFunction();
if(f==NULL)
  throw(NullPtrError("ArrayFitter::FitResidus(TVector<T>...) GeneraFit==NULL\n"));
int npar =  gfit.GetNPar();
if(npar==0)
  throw(SzMismatchError("ArrayFitter::FitResidus(TVector<T>...) GeneraFit 0 parametre\n"));
double* par = new double[npar];
{for(int i=0;i<npar;i++) par[i] = gfit.GetParm(i);}
TVector<T> v(vec);
v.SetTemp(true);
for(uint_4 i=0;i<vec.NElts();i++) {
  double x = xorg+i*dx;
  v(i) -= f->Value(&x,par);
}
delete [] par;
return v;
}

template <class T>
TVector<T> 
ArrayFitter<T>::FitFunction(TVector<T> const & vec, GeneralFit& gfit,
			  double xorg=0.,double dx=1.)
// Retourne une classe contenant les residus du fit ``gfit''.
// La coordonnee de l'element (i) est -> x = xorg + i*dx
{
if(vec.NElts()<=0)
  throw(SzMismatchError("ArrayFitter::FitResidus(TVector<T>...)  size mismatch\n"));
GeneralFunction* f = gfit.GetFunction();
if(f==NULL)
  throw(NullPtrError("ArrayFitter::FitResidus(TVector<T>...) GeneraFit==NULL\n"));
int npar =  gfit.GetNPar();
if(npar==0)
  throw(SzMismatchError("ArrayFitter::FitResidus(TVector<T>...) GeneraFit 0 parametre\n"));
double* par = new double[npar];
{for(int i=0;i<npar;i++) par[i] = gfit.GetParm(i);}
TVector<T> v(vec.NElts());
v.SetTemp(true);
for(uint_4 i=0;i<vec.NElts();i++) {
  double x = xorg+i*dx;
  v(i) = f->Value(&x,par);
}
delete [] par;
return v;
}


///////////////////////////////////////////////////////////////
#ifdef __CXX_PRAGMA_TEMPLATES__
#pragma define_template ArrayFitter<r_4>
#pragma define_template ArrayFitter<r_8>
#endif

#if defined(ANSI_TEMPLATES) || defined(GNU_TEMPLATES)
template class ArrayFitter<r_4>;
template class ArrayFitter<r_8>;
#endif