#include <nagdmc.h>

/*
  handle_info() prints to screen information based on the value of the info
  parameter.
*/
int
handle_info(const char func[], int info);

int
main(void) {
    const char    file[] = {"iris.dat"};
    long          rec1 = 0;
    long          nvar = 5;
    long          nrec = 150;
    long          dblk = 150;
    double       *data = 0;
    long          nxvar = 0;
    long         *xvar = 0;
    long          yvar = 4;
    long          nrbf = 3;
    double        cen[] = {59.36,27.7,42.6,13.26,
                           50.06,34.28,14.62,2.46,
                           65.88,29.74,55.52,20.26};
    int           dtype = 0;
    int           mtype = 1;
    double       *m = 0;
    int           rbftype = 4;
    double        alpha = 1.0;
    double        tau = 1.0e-6;
    double        lambda = 0.0;
    int           info = 0;
    double        fv;
    double       *w = 0;
    double       *yhat = 0;
    double       *model = 0;
    double       *rwork = 0;
    FILE         *fp = 0;
    long          i, j, k;

    /*
      Memory allocation.
    */
    if (mtype == 1)
        k = nvar;
    else if (mtype == 3)
        k = nvar * (1 + nvar) / 2;
    else
        k = 0;

    k += 6 + (nvar+1)*(nrbf+1);

    if (!(rwork = (double *)malloc((1+nrbf+nrec+k) * sizeof(double))) ||
        !(data = (double *)malloc(dblk*nvar * sizeof(double)))) {
        printf(" Memory allocation failure.\n\n");
        return 2;
    }
    w = rwork;
    yhat = w + nrbf + 1;
    model = yhat + nrec;

    /*
      Read data values.
    */
    if ((fp = fopen(file,"r")) == 0) {
        printf(" Data file named %s was not found.\n\n",file);
        free(data);
        free(rwork);
        return 2;
    }
    
    for (i=0; i<dblk; ++i) {
        for (j=0; j<nvar; ++j) 
            fscanf(fp,"%lf ",&data[i*nvar+j]);
    }

    fclose(fp);

    /*
      Summary information.
    */

    /*
      Fit RBF model.
    */
    nagdmc_rbf(rec1,nvar,nrec,dblk,data,nxvar,xvar,yvar,nrbf,cen,dtype,
               mtype,m,rbftype,alpha,tau,lambda,w,yhat,model,&info);

    if (handle_info("nagdmc_rbf",info)) {
        free(rwork);
        return 2;
    }

    /*
      Use fitted RBF model to predict data.
    */
    printf("\n Observed\tPrediction\tResidual \n\n");
    for (i=0; i<nrec; ++i) {
        nagdmc_predict_rbf(&data[i*nvar],model,&fv,&info);

        if (handle_info("nagdmc_predict_rbf",info)) {
            free(rwork);
            return 2;
        }
        
        printf(" %-8li\t%-8.4f\t%-8.4f\n",(long)data[i*nvar+yvar],fv,
               data[i*nvar+yvar]-fv);
    }

    /*
      Return allocated memory to the operating system.
    */
    free(rwork);

    return 0;
}

int
handle_info(const char func[], int info) {
    if (info == -999) {
        printf("\n Invalid licence, please contact NAG.\n\n");
        return 2;
    }
    else if (info > 0) {
        printf("\n Error code %i from %s\n\n",info,func);
        return 1;
    }
    else if (info < 0)
        printf("\n Information code %i from %s\n\n",info,func);

    return 0;
}