#include <ctype.h>
#include <string.h>
#include <stdio.h>
#ifndef NOMALLOC_H
#include <malloc.h>
#endif
#ifndef NOUNISTD_H
#include <unistd.h>
#endif

#include "fudgit.h"
#include "symbol.h"
#include "macro.h"
#include "code.h"
#include "math.tab.h"
#include "functions.h"
#include "head.h"

/* All the global internal variables  */
int Ft_Iter = 0;
int Ft_Mlist = 0;
int Ft_Mode = 0;
double *Ft_A = NULL;
double *Ft_DA = NULL;
double *Ft_Mfparx1 = NULL;
int *Ft_Miparx1 = NULL;
double **Ft_M1parxpar = NULL;
double **Ft_M2parxpar = NULL;
double **Ft_Mparxsamp = NULL;
double Ft_Q = 0.0;
double Ft_Cortest[3];
int Ft_Samples = 0;
int Ft_Debug = 0;
int Ft_Check = INF_CHK | NAN_CHK | EDOM_CHK | ERANGE_CHK;
int Ft_Expandhist = 1;
int Ft_Dolevel = 0;
double *Ft_X2;
double *Ft_Data;
double *Ft_If_value;
double *Ft_Param;

FILE *Ft_Inread;
FILE *Ft_Outprint;
char Ft_Outname[TOKENSIZE+8];
char Ft_Inname[TOKENSIZE+8];
char Ft_Home[PATH_MAXIM+4];
char Ft_Shell[TOKENSIZE+8];
char Ft_Pager[TOKENSIZE+8];
char Ft_Prompt_cm[MAXPROMPT+4];
char Ft_Prompt_fm[MAXPROMPT+4];
char Ft_Prompt_pm[MAXPROMPT+4];
char Ft_Format[TOKENSIZE+8];
char Ft_TFormat[TOKENSIZE+8];
char Ft_Vformat[TOKENSIZE+8];
char Ft_UFunction[MAXMACRO+8] = { '\0' };
char Ft_Pname[TOKENSIZE+8] = { '\0' };
char Ft_ReadFile[TOKENSIZE+8];
char Ft_Cwd[PATH_MAXIM];
char Ft_Tmp[TOKENSIZE+8];
char *Ft_Plotting[MAXPARG];
char Ft_Comchar;
int Ft_Methi;
int Ft_Funci;

Meth Ft_Method[METHNUM] = {
	{"none", "none"},
	{"ls_r!eg", "least square linear regression"},
	{"lad!_reg", "least absolute deviation linear regression"},
	{"ls_f!it", "least square linear fit"},
	{"sv!d_fit", "singular value decomposition linear fit"},
	{"ml!_fit", "Marquardt-Levenberg non-linear fit"}
};

Func Ft_Function[FUNCNUM] = {
	{"none", "none"},
	{"str!aight", "straight line"},
	{"po!lynomial", "polynomial"},
	{"leg!endre", "Legendre polynomial"},
	{"si!ne", "sine series"},
	{"cos!ine", "cosine series"},
	{"ex!ponential", "exponential series"},
	{"gau!ssian", "gaussian series"},
	{"us!er", "user defined function"}
};

extern double *Ft_dvector(int nl, int nh);
extern double **Ft_dmatrix(int nrl, int nrh, int ncl, int nch);
extern void Ft_free_dvector(double *v, int nl, int nh);
extern void Ft_free_dmatrix(double **m, int nrl, int nrh, int ncl, int nch);
extern void Ft_free_ivector(int *v, int nl, int nh);
extern int *Ft_ivector(int nl, int nh);
extern int Ft_exit(int);
extern int Ft_symremove (char *name, int verb);

int Ft_initsetup(void)
{
    Symbol *sym, *Ft_lookup(char *);
    char *cp, *getenv(const char *);
    int i;

    sym = Ft_lookup("Cwd");
	sym->u.str = Ft_Cwd;
    sym = Ft_lookup("ReadFile");
	sym->u.str = Ft_ReadFile;
    sprintf(Ft_ReadFile, "none");
    sym = Ft_lookup("Tmp");
#ifdef AMIGA
    sprintf(Ft_Tmp, "t:fudgit%d", getpid());
#else
    sprintf(Ft_Tmp, "/tmp/fudgit%d", getpid());
#endif
	sym->u.str = Ft_Tmp;
    if ((cp = getenv("PAGER")))
        sprintf(Ft_Pager, "%s", cp);
    else
        sprintf(Ft_Pager, "%s", DEFPAGER);
    if ((cp = getenv("SHELL")))
        sprintf(Ft_Shell, "%s", cp);
    else
        sprintf(Ft_Shell, "%s", DEFSHELL);
    if ((cp = getenv("HOME")))
        sprintf(Ft_Home, "%s", cp);
    else {
#ifdef AMIGA
        sprintf(Ft_Home, "./");
#else
        fputs("Fatal: Could not find home directory!\n", stderr);
        Ft_exit(1);
#endif
    }
    sprintf(Ft_Format, "%s", FORMAT);
    sprintf(Ft_TFormat, "\t%s", FORMAT);
    sprintf(Ft_Vformat, "%s", VFORMAT);
    sprintf(Ft_Prompt_cm, "%s", PROMPT_CM);
    sprintf(Ft_Prompt_fm, "%s", PROMPT_FM);
    sprintf(Ft_Prompt_pm, "%s", PROMPT_PM);
    for (i=0;i<MAXPARG-1;i++) {
        if ((Ft_Plotting[i] = (char *)calloc(TOKENSIZE+1, 1)) == NULL) {
            fputs("Fatal: Allocation error.\n", stderr);
            Ft_exit(1);
        }
    }
    sprintf(Ft_Plotting[0], "%s", PLOTTING);
    Ft_Plotting[1][0] = '\0';
    Ft_Comchar = '#';
    Ft_Samples = MAXPTS;
    sym = Ft_lookup("data");
	Ft_Data = &(sym->u.val);
    sym = Ft_lookup("chi2");
	Ft_X2 = &(sym->u.val);
    sym = Ft_lookup("param");
    Ft_Param = &(sym->u.val);
    sym = Ft_lookup("if_value");
    Ft_If_value = &(sym->u.val);
    Ft_Methi = 0;
    Ft_Funci = 0;
    Ft_Iter = ITER;
	Ft_Outprint = stdout;
	Ft_Inread = stdin;
	strcpy(Ft_Outname, "stdout");
	strcpy(Ft_Inname, "stdin");

    return(0);
}

/* defines the name and number of parameters  */
int Ft_setparam(char *name, int n)
{
    int i;
    char dname[TOKENSIZE+6];
    Symbol *sym;
    extern Symbol *Ft_install(char *, int, int);

    if (!isupper((int)*name)) {
        fprintf(stderr, "%s: Illegal vector name.\n", name);
        return(ERRR);
    }
    for (i=1;i<strlen(name);i++) {
        if (!isupper((int)name[i]) || !isdigit((int)name[i]))  {
            fprintf(stderr, "%s: Illegal vector name.\n", name);
            return(ERRR);
        }
    }
    if (strlen(Ft_Pname)) {
        sprintf(dname, "D%s", Ft_Pname);
        Ft_symremove(Ft_Pname, 0);
        Ft_symremove(dname, 0);
    }
    sprintf(Ft_Pname, "%s", name);
    sym = Ft_install(Ft_Pname, PARAM, n);
    Ft_A = sym->u.vec;
    sprintf(dname, "D%s", Ft_Pname);
    sym = Ft_install(dname, PARAM, n);
    Ft_DA = sym->u.vec;
    /* Allocate internal matrices */
    if (Ft_Mfparx1 != (double *)NULL) {
        Ft_free_dvector(Ft_Mfparx1, 1, (int) *Ft_Param);
    }
    if ((Ft_Mfparx1 = Ft_dvector(1, n)) == (double *)NULL) {
        return(ERRR);
    }
    if (Ft_Miparx1 != (int *)NULL) {
        Ft_free_ivector(Ft_Miparx1, 1, (int) *Ft_Param);
    }
    if ((Ft_Miparx1 = Ft_ivector(1, n)) == (int *)NULL) {
        return(ERRR);
    }
    if (Ft_M1parxpar != (double **)NULL) {
        Ft_free_dmatrix(Ft_M1parxpar, 1, (int) *Ft_Param, 1, (int) *Ft_Param);
    }
    if ((Ft_M1parxpar = Ft_dmatrix(1, n, 1, n)) == (double**)NULL) {
        return(ERRR);
    }
    if (Ft_M2parxpar != (double **)NULL) {
        Ft_free_dmatrix(Ft_M2parxpar, 1, (int) *Ft_Param, 1, (int) *Ft_Param);
    }
    if ((Ft_M2parxpar = Ft_dmatrix(1, n, 1, n)) == (double**)NULL) {
        return(ERRR);
    }
    if (Ft_Mparxsamp != (double **)NULL) {
        free(Ft_Mparxsamp+1);
    }
    /* Make my own matrix skeleton */
    Ft_Mparxsamp = (double **)malloc((unsigned)n*sizeof(double *));
    if (Ft_Mparxsamp == (double **)NULL) {
        fputs("set parameters: Allocation error.\n", stderr);
        return(ERRR);
    }
    Ft_Mparxsamp--;
       
    *Ft_Param = n;
    Ft_Mlist = 0;
    return(0);
}

int Ft_showsetup(void)
{
    int i = 0;

    fprintf(stdout, "%28s: \"%s\"\n", "ReadFile", Ft_ReadFile);
    fprintf(stdout, "%28s: %s\n",
    "Fitting method", Ft_Method[Ft_Methi].name);
    fprintf(stdout, "%28s: \"%d\"\n", "Iteration number", Ft_Iter);
    fprintf(stdout, "%28s: %s\n",
    "Function to fit", Ft_Function[Ft_Funci].name);
    fprintf(stdout, "%28s: %d\n", "Number of parameters", (int) *Ft_Param);
    fprintf(stdout, "%28s: %d points\n", "Current capacity", Ft_Samples);
    fprintf(stdout, "%28s: %d\n", "Number of data points", (int) *Ft_Data);
    fprintf(stdout, "%28s: ", "Plotting program");
    while (Ft_Plotting[i][0]) {
        fprintf(stdout, "%s ", Ft_Plotting[i]);
        i++;
    }
    fputc('\n', stdout);
    fprintf(stdout, "%28s: %s\n", "Pager program", Ft_Pager);
    fprintf(stdout, "%28s: \"%s\"\n", "Output format", Ft_Format);
    fprintf(stdout, "%28s: '%c'\n", "Comment character", Ft_Comchar);
    fprintf(stdout, "%28s: \"%s\"\n", "Temporary file", Ft_Tmp);
    return(0);
}

int Ft_showfit(void)
{
    int i, j;

    if ((int) *Ft_Param == 0) {
        fprintf(stderr, "No parameter!\n");
        return(ERRR);
    }
    for (i=1;i <= (int) *Ft_Param;i++) {
        fprintf(stdout, "\t%s[%d]: ", Ft_Pname, i);
        fprintf(stdout, Ft_Format, Ft_A[i]);
        fputs("\t +/- ", stdout);
        fprintf(stdout, Ft_Format, Ft_DA[i]);
        fputc('\n', stdout);
    }
    if (Ft_Methi == LA_REG) {
        fputs("Mean absolute deviation: ", stdout);
        fprintf(stdout, Ft_Format, *Ft_X2);
        fputc('\n', stdout);
    }
    else {
        fputs("Chi 2: ", stdout);
        fprintf(stdout, Ft_Format, *Ft_X2);
        fputc('\n', stdout);
    }
    if (Ft_Methi == LS_REG) {
        fputs("Goodness-of-fit probability: ", stdout);
        fprintf(stdout, Ft_Format, Ft_Q);
        fputc('\n', stdout);
    }
    if (Ft_Mlist) {
        fputs("Adjusting:", stdout);
        for (i=1; i<= Ft_Mlist; i++) {
            fprintf(stdout, " %d", Ft_Miparx1[i]);
        }
        fputc('\n', stdout);
    }
    if (Ft_Methi == ML_FIT || Ft_Methi == SVD_FIT || Ft_Methi == LS_FIT) {
        fputs("Covariance matrix:\n", stdout);
        for (i=1;i<= (int) *Ft_Param;i++) {
            fputs(" |", stdout);
            for (j=1;j<= (int) *Ft_Param;j++) {
                fprintf(stdout, "\t% 10.8e", Ft_M1parxpar[i][j]);
            }
            fputs("\t |\n", stdout);
        }
    }
    if (Ft_Methi == ML_FIT) {
        fputs("Curvature matrix:\n", stdout);
        for (i=1;i<= (int) *Ft_Param;i++) {
            fputs(" |", stdout);
            for (j=1;j<= (int) *Ft_Param;j++) {
                fprintf(stdout, "\t% 10.8e", Ft_M2parxpar[i][j]);
            }
            fputs("\t |\n", stdout);
        }
    }
    if (Ft_Methi == LS_REG || Ft_Methi == LA_REG) {
        fputs("Linear correlation tests\n", stdout);
        fprintf(stdout, "Correlation coefficient: %g\n", Ft_Cortest[0]);
        fprintf(stdout, "Fisher's `z' coefficient: %g\n", Ft_Cortest[2]);
        fprintf(stdout, "Significance: %g\n", Ft_Cortest[1]);
    }

    fputc('\n', stdout);
    return(0);
}