/******************************************************************************

                      ===================
        Network:      Self-Organizing Map
                      ===================

        Application:  Control
                      Pole Balancing Problem

        Author:       Karsten Kutza
        Date:         6.6.96

        Reference:    T. Kohonen	
                      Self-Organized Formation
                      of Topologically Correct Feature Maps
                      Biological Cybernetics, 43, pp. 59-69, 1982

 ******************************************************************************/




/******************************************************************************
                            D E C L A R A T I O N S
 ******************************************************************************/


#include <stdlib.h>
#include <stdio.h>
#include <math.h>


typedef int           BOOL;
typedef int           INT;
typedef double        REAL;

#define FALSE         0
#define TRUE          1
#define NOT           !
#define AND           &&
#define OR            ||

#define MIN_REAL      -HUGE_VAL
#define MAX_REAL      +HUGE_VAL
#define MIN(x,y)      ((x)<(y) ? (x) : (y))
#define MAX(x,y)      ((x)>(y) ? (x) : (y))

#define PI            (2*asin(1))
#define sqr(x)        ((x)*(x))


typedef struct {                     /* A LAYER OF A NET:                     */
        INT           Units;         /* - number of units in this layer       */
        REAL*         Output;        /* - output of ith unit                  */
        REAL**        Weight;        /* - connection weights to ith unit      */
        REAL*         StepSize;      /* - size of search steps of ith unit    */
        REAL*         dScoreMean;    /* - mean score delta of ith unit        */
} LAYER;

typedef struct {                     /* A NET:                                */
        LAYER*        InputLayer;    /* - input layer                         */
        LAYER*        KohonenLayer;  /* - Kohonen layer                       */
        LAYER*        OutputLayer;   /* - output layer                        */
        INT           Winner;        /* - last winner in Kohonen layer        */
        REAL          Alpha;         /* - learning rate for Kohonen layer     */
        REAL          Alpha_;        /* - learning rate for output layer      */
        REAL          Alpha__;       /* - learning rate for step sizes        */
        REAL          Gamma;         /* - smoothing factor for score deltas   */
        REAL          Sigma;         /* - parameter for width of neighborhood */
} NET;


/******************************************************************************
        R A N D O M S   D R A W N   F R O M   D I S T R I B U T I O N S
 ******************************************************************************/


void InitializeRandoms()
{
  srand(4715);
}      


REAL RandomEqualREAL(REAL Low, REAL High)
{
  return ((REAL) rand() / RAND_MAX) * (High-Low) + Low;
}      


REAL RandomNormalREAL(REAL Mu, REAL Sigma)
{
  REAL x,fx;

  do {
    x = RandomEqualREAL(Mu-3*Sigma, Mu+3*Sigma);
    fx = (1 / (sqrt(2*PI)*Sigma)) * exp(-sqr(x-Mu) / (2*sqr(Sigma)));
  } while (fx < RandomEqualREAL(0, 1));
  return x;
}      


/******************************************************************************
               A P P L I C A T I O N - S P E C I F I C   C O D E
 ******************************************************************************/


#define ROWS          25
#define COLS          25

#define N             2
#define C             (ROWS * COLS)
#define M             1

#define TRAIN_STEPS   10000
#define BALANCED      100

FILE*                 f;


void InitializeApplication(NET* Net)
{
  INT i;
   
  for (i=0; i<Net->KohonenLayer->Units; i++) {
    Net->KohonenLayer->StepSize[i] = 1;
    Net->KohonenLayer->dScoreMean[i] = 0;
  }
  f = fopen("SOM.txt", "w");
}


void WriteNet(NET* Net)
{
  INT  r,c;
  REAL x,y,z;

  fprintf(f, "\n\n\n");
  for (r=0; r<ROWS; r++) {
    for (c=0; c<COLS; c++) {
      x = Net->KohonenLayer->Weight[r*ROWS+c][0];
      y = Net->KohonenLayer->Weight[r*ROWS+c][1];
      z = Net->OutputLayer->Weight[0][r*ROWS+c];
      fprintf(f, "([%5.1f