/*                  Zadanie z predmetu: NEURONOVE SIETE
                    Projekt c. 1: Back-Propagation momentum
                    Riesitel: Tomas Smajda, 4. roc., 19/11/2001
                                                                             */
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define	N1	     2	   // pocet neuronov vo vstupnej vrstve
#define	N2     11	   // pocet neuronov v 1. skrytej vrstve
#define	N3	     7	   // pocet neuronov v 2. skrytej vrstve
#define	N4	     1	   // pocet neuronov vo vystupnej vrstve
#define PRAH     1      // hodnota prahu

#define aktiv_fcia(x)       ( 1/(1+exp(-1*x)) )
#define deriv_aktiv_fcia(x) ( exp(-1*x)/((1+exp(-1*x))*(1+exp(-1*x))))

double in2[N2], in3[N3], in4[N4];                   // hodnoty neuronov pred zaktivovanim
double ou1[N1], ou2[N2], ou3[N3], ou4[N4];          // hodnoty neuronov po zaktivovani
double sv1[N1][N2], sv2[N2][N3], sv3[N3][N4];       // jednotlive vahy (k)
double d2[N2], d3[N3], d4[N4];                      // jednotlive delty
double svv1[N1][N2], svv2[N2][N3], svv3[N3][N4];    // jednotlive vahy (k-1)
double svvv1[N1][N2], svvv2[N2][N3], svvv3[N3][N4]; // jednotlive vahy (k-2)

double ev, Gama, Alpha;
FILE *train_file, *test_file, *error_file, *vahy;
double x[1251], y[1251];                            // polia na ulozenie vstupov
int  tr_value1[1251], tr_value2, te_value1, te_value2;
int stvorec_zle=0, spirala_zle=0, stvorec_dobre=0, spirala_dobre=0;

void zivot(void)
{
    int i,j;

    for(i=0;i<N2;i++){
       in2[i]=0;
       for(j=0;j<N1;j++){
	       in2[i]=in2[i]+(ou1[j]*sv1[j][i]);
       }
       in2[i]=in2[i]-PRAH;
       ou2[i]=aktiv_fcia(in2[i]);
    }
    for(i=0;i<N3;i++){
       in3[i]=0;
       for(j=0;j<N2;j++){
	       in3[i]=in3[i]+(ou2[j]*sv2[j][i]);
       }
       in3[i]=in3[i]-PRAH;
       ou3[i]=aktiv_fcia(in3[i]);
    }
    for(i=0;i<N4;i++){
       in4[i]=0;
       for(j=0;j<N3;j++){
	       in4[i]=in4[i]+(ou3[j]*sv3[j][i]);
       }
       in4[i]=in4[i]-PRAH;
       ou4[i]=aktiv_fcia(in4[i]);
    }
}

void Back_Propagation(void)
{
    int i,j;

    // vypocet delty
    for(i=0;i<N4;i++){
       d4[i]=(ev-ou4[i])*deriv_aktiv_fcia(in4[i]);
    }
    for(i=0;i<N3;i++){
       d3[i]=0;
       for(j=0;j<N4;j++){
	       d3[i]=d3[i]+(d4[j]*sv3[i][j]);
       }
       d3[i]=d3[i]*deriv_aktiv_fcia(in3[i]);
    }
    for(i=0;i<N2;i++){
       d2[i]=0;
       for(j=0;j<N3;j++){
	       d2[i]=d2[i]+(d3[j]*sv2[i][j]);
       }
	    d2[i]=d2[i]*deriv_aktiv_fcia(in2[i]);
    }

    // vypocet novych vah
    for(i=0;i<N3;i++){
       for(j=0;j<N4;j++){
          svvv3[i][j]=svv3[i][j];
          svv3[i][j]=sv3[i][j];
          sv3[i][j]+=(Gama*d4[j]*ou3[i])+(Alpha*(fabs(svvv3[i][j]-svv3[i][j])));
       }
    }
    for(i=0;i<N2;i++){
       for(j=0;j<N3;j++){
          svvv2[i][j]=svv2[i][j];
          svv2[i][j]=sv2[i][j];
	       sv2[i][j]+=(Gama*d3[j]*ou2[i])+(Alpha*(fabs(svvv2[i][j]-svv2[i][j])));
       }
    }
    for(i=0;i<N1;i++){
       for(j=0;j<N2;j++){
          svvv1[i][j]=svv1[i][j];
          svv1[i][j]=sv1[i][j];
          sv1[i][j]+=(Gama*d2[j]*ou1[i])+(Alpha*(fabs(svvv1[i][j]-svv1[i][j])));
       }
    }
}

int main(void)
{
   int i, j, epocha=0;
   double J, uspesnost;
   Gama=0.1;
   Alpha=0;

   if ((train_file = fopen("spiral_train.pat", "rt")) == NULL){
       fprintf(stderr, "Trenovaci subor sa neda otvorit.\n");
       printf("Stlacte ENTER a aplikacia sa skonci...");
       getchar();
       exit(1);
   }
   if ((error_file = fopen("error_file.txt", "a")) == NULL){
       fprintf(stderr, "Vystupny - chybovy subor sa neda otvorit.\n");
       printf("Stlacte ENTER a aplikacia sa skonci...");
       getchar();
       exit(1);
   }
   if ((test_file = fopen("spiral_test.pat", "r")) == NULL){
       fprintf(stderr, "Testovaci subor sa neda otvorit.\n");
       printf("Stlacte ENTER a aplikacia sa skonci...");
       getchar();
       exit(1);
   }
   if ((vahy = fopen("vahy.txt", "a")) == NULL){
       fprintf(stderr, "Nedaju sa zapisat vahy do vystupneho suboru.\n");
       printf("Stlacte ENTER a aplikacia sa skonci...");
       getchar();
       exit(1);
   }

   // prvotna inicializacia vah
   printf("Boli inicializovane nasledovne vahy...\n");
   for(i=0;i<N1;i++){
       for(j=0;j<N2;j++){
          sv1[i][j]=pow(-1,j)*rand()/RAND_MAX;
          printf("sv1[%d][%d]=%lf\n", i+1 , j+1, sv1[i][j]);
       }
   }
   getchar();
   for(i=0;i<N2;i++){
       for(j=0;j<N3;j++){
          sv2[i][j]=pow(-1,j)*rand()/RAND_MAX;
          printf("sv2[%d][%d]=%lf\n", i+1 , j+1, sv2[i][j]);
       }
   }
   getchar();
   for(i=0;i<N3;i++){
       for(j=0;j<N4;j++){
          sv3[i][j]=pow(-1,j)*rand()/RAND_MAX;
          printf("sv3[%d][%d]=%lf\n", i+1 , j+1, sv3[i][j]);
       }
   }

   // naplni sa uciace pole
   for(i=0;i<1251;i++){
      fscanf(train_file,"%lf %lf %d %d\n", &x[i], &y[i], &tr_value1[i], &tr_value2);
   }

   // faza ucenia
   printf("Stlacte ENTER a prebehne ucenie Neuronovej siete.");
   getchar();
   do{
      J=0;      for(i=0;i<1251;i++){
         ou1[0]=x[i];
         ou1[1]=y[i];
         ev=(double)tr_value1[i];
         zivot();                           // jeden cyklus ucenia
         Back_Propagation();
         J = J + ((ev-ou4[0]) * (ev-ou4[0]));
      }
      J = 0.5 * (J/1251.0);
      epocha++;
      if(J<0.10)  Gama=0.09;
      if(J<0.097) Gama=0.05;
      if(J<0.09)  Gama=0.01;
      if(J<0.082) Gama=0.005;
      //if (!(epocha%100)){
      //   printf(".");
      //}
      printf("J(t)=%lf\t Epocha=%d\n", J, epocha);
      fprintf(error_file,"J(t)=%lf\t Epocha=%d\n", J, epocha); // zapis chyby pri uceni do suboru
   }while(epocha<20000);

   // ulozenie konecnych vah do suboru "vahy.txt"
   fprintf(vahy,"A takto vyzeraju jednotlive vahy...\n");
   for(i=0;i<N1;i++){
       for(j=0;j<N2;j++){
          fprintf(vahy,"sv1[%d][%d]=%lf\n", j+1, i+1, sv1[i][j]);
       }
   }
   fprintf(vahy,"\n");
   for(i=0;i<N2;i++){
       for(j=0;j<N3;j++){
          fprintf(vahy,"sv2[%d][%d]=%lf\n", j+1, i+1, sv2[i][j]);
       }
   }
   fprintf(vahy,"\n");
   for(i=0;i<N3;i++){
       for(j=0;j<N4;j++){
          fprintf(vahy,"sv3[%d][%d]=%lf\n", j+1, i+1, sv3[i][j]);
       }
   }

   fclose(error_file);
   fclose(train_file);
   fclose(vahy);

   // nasleduje faza zivota
   printf("\nPrebieha testovanie Neuronovej siete...\n");
   for(i=0;i<10002;i++){
      fscanf(test_file,"%lf %lf %d %d\n", &ou1[0], &ou1[1], &te_value1, &te_value2);
      zivot();
      if((ou4[0]<=0.5) && (te_value1==1)) stvorec_zle++;
      if((ou4[0]>=0.5) && (te_value1==0)) spirala_zle++;
      if((ou4[0]<0.5) && (te_value1==0))  stvorec_dobre++;
      if((ou4[0]>0.5) && (te_value1==1))  spirala_dobre++;
   }
   uspesnost=(stvorec_dobre+spirala_dobre)*100.00/10002.00;

   printf("\n%d krat bol zle identifikovany stvorec.\n", stvorec_zle);
   printf("%d krat bola zle identifikovana spirala.\n", spirala_zle);
   printf("%d krat bol dobre identifikovany stvorec.\n", stvorec_dobre);
   printf("%d krat bola dobre identifikovana spirala.\n", spirala_dobre);
   printf("Celkova uspesnost je %.2lf%.\n\n", uspesnost);

   fclose(test_file);
   printf("Stlacte ENTER a aplikacia sa skonci...");
   getchar();
}

