#include <stdio.h>
#include <math.h>
#include <errno.h>
#include <stdlib.h>
#include <string.h>
#include <memory.h>
#include "neuron_switch.h"

/* -------------------------------------------------------------------------------
 Lecture des differrentes Look Up Tables
 -------------------------------------------------------------------------------  */
void neuron_lect_LUTswitch()
{
  FILE *fic;
  int i,j,poub;
  char *ligne=malloc(sizeof(char)*150);
  float fpoub;

  /* ----- LUTs for Rrs490/Rrs555 > 1.5 ------ */
  if( (fic=fopen(rsup_LUT_POIDS,"r")) == NULL) {perror(rsup_LUT_POIDS); exit(-1);}
  fgets(ligne,150,fic);
  for(i=0; i<rsup_NC1; i++)
    fscanf(fic,"%d %d %f",&poub,&poub,&rsup_b1[i]);
  for(i=0; i<rsup_NC2; i++)
    fscanf(fic,"%d %d %f",&poub,&poub,&rsup_b2[i]);
  fscanf(fic,"%d %d %f",&poub,&poub,&rsup_b3);
  for(j=0; j<rsup_NE; j++){
   for(i=0; i<rsup_NC1; i++)
     fscanf(fic,"%d %d %f",&poub,&poub,&rsup_w1[j][i]);
  }
  for(j=0; j<rsup_NC1; j++){
    for(i=0; i<rsup_NC2; i++)
      fscanf(fic,"%d %d %f",&poub,&poub,&rsup_w2[j][i]);
  }
  for(i=0; i<rsup_NC2; i++)
    fscanf(fic,"%d %d %f",&poub,&poub,&rsup_w3[i]);
  fclose(fic);

  if( (fic=fopen(rsup_LUT_MOY,"r")) == NULL) {perror(rsup_LUT_MOY); exit(-1);}
  for(i=0; i<rsup_NE-1; i++)
    fscanf(fic,"%f",&rsup_moy[i]);
  fscanf(fic,"%f",&fpoub);
  fscanf(fic,"%f",&rsup_moy[rsup_NES-2]);
  fscanf(fic,"%f",&fpoub);
  fscanf(fic,"%f",&rsup_moy[rsup_NES-1]);
  fclose(fic);
  if( (fic=fopen(rsup_LUT_ECART,"r")) == NULL) {perror(rsup_LUT_ECART); exit(-1);}
  for(i=0; i<rsup_NE-1; i++)
    fscanf(fic,"%f",&rsup_ecart[i]);
  fscanf(fic,"%f",&fpoub);
  fscanf(fic,"%f",&rsup_ecart[rsup_NES-2]);
  fscanf(fic,"%f",&fpoub);
  fscanf(fic,"%f",&rsup_ecart[rsup_NES-1]);
  fclose(fic);


/* ----- LUTs for Rrs490/Rrs555 <= 1.5 ------ */
  if( (fic=fopen(rinf_LUT_POIDS,"r")) == NULL) {perror(rinf_LUT_POIDS); exit(-1);}
  fgets(ligne,150,fic);
  for(i=0; i<rinf_NC1; i++)
    fscanf(fic,"%d %d %f",&poub,&poub,&rinf_b1[i]);
  for(i=0; i<rinf_NC2; i++)
    fscanf(fic,"%d %d %f",&poub,&poub,&rinf_b2[i]);
  fscanf(fic,"%d %d %f",&poub,&poub,&rinf_b3);
  for(j=0; j<rinf_NE; j++){
   for(i=0; i<rinf_NC1; i++)
     fscanf(fic,"%d %d %f",&poub,&poub,&rinf_w1[j][i]);
  }
  for(j=0; j<rinf_NC1; j++){
    for(i=0; i<rinf_NC2; i++)
      fscanf(fic,"%d %d %f",&poub,&poub,&rinf_w2[j][i]);
  }
  for(i=0; i<rinf_NC2; i++)
    fscanf(fic,"%d %d %f",&poub,&poub,&rinf_w3[i]);
  fclose(fic);

  if( (fic=fopen(rinf_LUT_MOY,"r")) == NULL) {perror(rinf_LUT_MOY); exit(-1);}
  for(i=0; i<rinf_NES-1; i++)
    fscanf(fic,"%f",&rinf_moy[i]);
  fscanf(fic,"%f",&fpoub);
  fscanf(fic,"%f",&rinf_moy[rinf_NES-1]);
  fclose(fic);
  if( (fic=fopen(rinf_LUT_ECART,"r")) == NULL) {perror(rinf_LUT_ECART); exit(-1);}
  for(i=0; i<rinf_NES-1; i++)
    fscanf(fic,"%f",&rinf_ecart[i]);
  fscanf(fic,"%f",&fpoub);
  fscanf(fic,"%f",&rinf_ecart[rinf_NES-1]);
  fclose(fic);

}

/* -------------------------------------------------------------------------------
 Calcul du Kd a partir des poids
 - Input:
  input[NE] = Rrs 412 443 490 510 555 Lambda
 ------------------------------------------------------------------------------- */
float rsup_neuron_passe_avant(float input[rsup_NE+1])
{
  float a[rsup_NC1], b[rsup_NC2], y=0.0, x[rsup_NE];
  int i,j;

  /* Normalisation */
  for(i=0; i<rsup_NE; i++){
    x[i] = ((2./3.)*(input[i]-rsup_moy[i]))/rsup_ecart[i];
  }

  for(i=0;i<rsup_NC1;i++){
    a[i] = 0.0;
    for(j=0;j<rsup_NE;j++){
      a[i] += (x[j]*rsup_w1[j][i]);
    }
    a[i] = 1.715905*(float)tanh((2./3.)*(double)(a[i] + rsup_b1[i]));
  }
  for(i=0;i<rsup_NC2;i++){
    b[i] = 0.0;
    for(j=0;j<rsup_NC1;j++){
      b[i] += (a[j]*rsup_w2[j][i]);
    }
    b[i] = 1.715905*(float)tanh((2./3.)*(double)(b[i] + rsup_b2[i]));
  }
  for(j=0;j<rsup_NC2;j++){
    y += (b[j]*rsup_w3[j]);
  }

  /* Denormalisation */
  y = 1.5*(y + rsup_b3)*rsup_ecart[rsup_NES-1] + rsup_moy[rsup_NES-1];
  y = (float)pow(10.,y);
  return(y);
}

/* -------------------------------------------------------------------------------
 Calcul du Kd a partir des poids
 - Input:
  input[NE] = Rrs 412 443 490 510 555 670 Lambda
 ------------------------------------------------------------------------------- */
float rinf_neuron_passe_avant(float input[rinf_NE])
{
  float a[rinf_NC1], b[rinf_NC2], y=0.0, x[rinf_NE];
  int i,j;

  /* Normalisation */
  for(i=0; i<rinf_NE; i++){
    x[i] = ((2./3.)*(input[i]-rinf_moy[i]))/rinf_ecart[i];
  }

  for(i=0;i<rinf_NC1;i++){
    a[i] = 0.0;
    for(j=0;j<rinf_NE;j++){
      a[i] += (x[j]*rinf_w1[j][i]);
    }
    a[i] = 1.715905*(float)tanh((2./3.)*(double)(a[i] + rinf_b1[i]));
  }
  for(i=0;i<rinf_NC2;i++){
    b[i] = 0.0;
    for(j=0;j<rinf_NC1;j++){
      b[i] += (a[j]*rinf_w2[j][i]);
    }
    b[i] = 1.715905*(float)tanh((2./3.)*(double)(b[i] + rinf_b2[i]));
  }
  for(j=0;j<rinf_NC2;j++){
    y += (b[j]*rinf_w3[j]);
  }

  /* Denormalisation */
  y = 1.5*(y + rinf_b3)*rinf_ecart[rinf_NES-1] + rinf_moy[rinf_NES-1];
  y = (float)pow(10.,y);
  return(y);
}


/*int main (int argc, char *argv[])
{
  FILE *fic;
  float data_in[rinf_NE],result_in,result_out,res_norm, sumInf_rms=0., sumInf_rel=0., sumSup_rms=0., sumSup_rel=0.;
  int i,nb,lu;

  if( (fic=fopen("../LUTS/base_KdSeaWiFS_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_670_log_Kd_lambda_merge_seuil_15_publi_091213_test.dat","r")) == NULL) {perror("input"); exit(-1);}
  neuron_lect_LUTs();

  nb = 0;
  while((lu=fscanf(fic,"%f",&data_in[0])) == 1){
    for(i=1; i<rinf_NE; i++)
      fscanf(fic,"%f",&data_in[i]);
    fscanf(fic,"%f",&result_in);
    result_out = rinf_neuron_passe_avant(data_in);
    result_in = (float)pow(10.,result_in);
    sumInf_rms += (float)pow((double)(result_in-result_out),2.);
    sumInf_rel +=  (float)sqrt((double)(((result_in-result_out)/result_in) * ((result_in-result_out)/result_in)));
    printf("%f %f",result_in,result_out);

    data_in[5] = data_in[6];
    result_out = rsup_neuron_passe_avant(data_in);
    sumSup_rms += (float)pow((double)(result_in-result_out),2.);
    sumSup_rel +=  (float)sqrt((double)(((result_in-result_out)/result_in) * ((result_in-result_out)/result_in)));
    printf(" %f\n",result_out);
    nb++;
  }
  printf("rmseInf = %f, rel errInf= %f,  rmseSup = %f, rel errSup= %f  nb = %d\n",(float)sqrt((double)sumInf_rms/nb),sumInf_rel/nb, (float)sqrt((double)sumSup_rms/nb),sumSup_rel/nb, nb);



  fclose(fic);
  exit(1);
}*/