/* GAUSS FIT
   uses: nonlfit.c
   need: mrqmin.c, gaussj.c covsrt.c
*/
#include <math.h>

void fgauss_c(x,a,y,dyda)
float x,a[5],*y,dyda[5];
{
  double fac,ex,arg;

  *y=0.0;
  arg=(x-a[2])/a[3];
  ex=exp(-arg*arg);
  fac=a[1]*ex*2.0*arg;
  *y = a[1]*ex + a[4];
  dyda[1]=ex;
  dyda[2]=fac/a[3];
  dyda[3]=fac*arg/a[3];
  dyda[4]=1.0;
}

/* necessary for FORTRAN call in MIDAS */
/* CALL GAUSS_FIT(...) */
/* esoext.exe will interpret this call as GAUSST !! */
int gausst_(y,sig,ndata,param,npar,istat)
float *y, *sig, *param;
int   *ndata, *npar, *istat;
{
  gauss_fit(y-1,sig-1,*ndata,param-1,*npar,istat);
}

int gausst(y,sig,ndata,param,npar,istat)
float *y, *sig, *param;
int   *ndata, *npar, *istat;
{
  gauss_fit(y-1,sig-1,*ndata,param-1,*npar,istat);
}


int gauss_fit(y,sig,ndata,param,npar,istat)
float *y, *sig, *param;
int   ndata, npar, *istat;
{
/*  void  fgauss_c();		because of VMS... */
  int   listp[4], i;

  for (i=0;i<4;i++) listp[i] = i+1;
  *istat = non_lfit(y,sig,ndata,param,npar,listp-1,4,fgauss_c);
}