train.c

SOM-sd 是现在国外非常流行的一个神经网络的结构自组织特征映射算法

/*
  Contents: Routines used to train the som-sd.

  Author: Markus Hagenbuchner

  Comments and questions concerning this program package may be sent
  to '[e-mail]'


  ChangeLog:
    03/10/2006:
    - BugFix: File name for interrupted training runs was not unique if several
              processes terminated within one second.
    30/03/2006:
    - Bug removed: ComputeHexaDistance(.) now computes the correct distance.
    - Speed improvement: The Distance(.) functions now return a squared value
      which results in a speed improvement of about 4-8% depending on which
      Distance(.) function is used.
    - Speed improvement: Forced inlining of function AdaptVector(.) resulted in
      a speed gain of approx. 8%.
    20/05/2005:
    - Bug removed: Confusion of topology and neighborhood in TrainMap(.)

 */


/************/
/* Includes */
/************/
#include <ctype.h>
#include <float.h>
#include <limits.h>
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <signal.h>
#include <sys/types.h>
#include <time.h>
#include <unistd.h>
#include "common.h"
#include "data.h"
#include "fileio.h"
#include "system.h"
#include "train.h"
#include "utils.h"


/********************/
/* Global variables */
/********************/
int _save_then_exit_ = 0; /* Indicate if an interrupt was caught */



/* Begin functions... */

/******************************************************************************
Description: Function is used to compute the current alpha value (the
             learning rate) for the current training iteration, given the
             maximum number of training iterations, and given initial value
             of alpha at iter=0, and assuming a linearly decrease of alpha.

Return value: The alpha value for the given itaration assuming a linearly
              decreasing learning rate.
******************************************************************************/
FLOAT LinearDecrease(UNSIGNED iter, UNSIGNED length, FLOAT alphastart)
{
  return (alphastart * (FLOAT) (length - iter) / (FLOAT) length);
}

/******************************************************************************
Description: Function the current alpha value. Other parameters are ignored.
             It is a dummy function which is required to maintain compatibility
             with other alpha functions.

Return value: The alpha value for the given itaration which is assumed to
              remain unchanged.
******************************************************************************/
FLOAT ConstantAlpha(UNSIGNED iter, UNSIGNED length, FLOAT alphastart)
{
  return alphastart;
}

/******************************************************************************
Description: Function is used to compute the current alpha value (the
             learning rate) for the current training iteration, given the
             maximum number of training iterations, and given initial value
             of alpha at iter=0, and assuming an exponential decrease of alpha.
             The Contant INV_ALPHA_CONSTANT specifies the "steepness" of the
             decrease, and the lower limit of alpha.

Return value: The alpha value for the given itaration assuming a exponential
              decreasing learning rate.
******************************************************************************/
#define INV_ALPHA_CONSTANT 100.0
FLOAT ExponentialAlpha(UNSIGNED iter, UNSIGNED length, FLOAT alphastart)
{
  FLOAT c;
  
  c = length / INV_ALPHA_CONSTANT;
  
  return (alphastart * c / (c + iter));
}

/******************************************************************************
Description: Function is used to compute the current alpha value (the
             learning rate) for the current training iteration, given the
             maximum number of training iterations, and given initial value
             of alpha at iter=0, and assuming an sigmoidal decrease of alpha.
             The constant GAUSS_ALPHA_CONSTANT specifies the "steepness" of the
             decrease, and influences the lower limit of alpha.

Return value: The alpha value for the given itaration assuming a sigmoidal
              decreasing learning rate.
******************************************************************************/
#define GAUSS_ALPHA_CONSTANT 4.0
float SigmoidalAlpha(UNSIGNED iter, UNSIGNED length, FLOAT alphastart)
{
  return alphastart*exp(-((FLOAT)iter*iter*GAUSS_ALPHA_CONSTANT)/((FLOAT)length*length));
}

/******************************************************************************
Description: Compute the geographical distance between two codebook coordinates
             given a rectagonal neighborhood relationship.

Return value: The squared!! Eucledian distance between two points.
******************************************************************************/
FLOAT ComputeRectDistance(int bx, int by, int tx, int ty)
{
  FLOAT ret, diff;

  diff = bx - tx;
  ret = diff * diff;
  diff = by - ty;
  ret += diff * diff;

  return ret;
  //  return (FLOAT)sqrtf(ret);
}

/******************************************************************************
Description: Compute the geographical distance between two codebook coordinates
             given a hexagonal neighborhood relationship.

Return value: The squared!!  eucledian distance between two points.
******************************************************************************/
FLOAT ComputeHexaDistance(int bx, int by, int tx, int ty)
{
  FLOAT ret, dy;
  int dx;

  dx = bx - tx;
  dy = by - ty;

  if (dx & 1){
    if (tx & 1)
      dy += 0.5;
    else
      dy -= 0.5;
  }

  ret = dy * dy;
  ret += 0.75 * dx * dx;

  return ret;
  //  return (FLOAT)sqrtf(ret);
}

/******************************************************************************
Description: Compute the geographical distance between two codebook coordinates
             given a octagonal neighborhood relationship.

Return value: The squared distance between two points.
******************************************************************************/
FLOAT ComputeOctDistance(int bx, int by, int tx, int ty)
{
  FLOAT ret;

  ret = (FLOAT)max(abs(bx-tx), abs(by-ty));

  return ret*ret;
}

/******************************************************************************
Description: move a codebook vector towards another vector

Return value: 
******************************************************************************/
inline void AdaptVector(FLOAT *codebook, FLOAT *sample, UNSIGNED dim, FLOAT alpha)
{
  UNSIGNED i;

  for (i = 0; i < dim; i++)
    codebook[i] += alpha * (sample[i] - codebook[i]);
}

/******************************************************************************
Description: Find best matching codebook using the Eucledian distance meassure.

Return value: The best matching codebook is returned to parameter "winner".
******************************************************************************/
void FindWinnerEucledian(struct Map *map, struct Node *node, struct Graph *gptr, struct Winner *winner)
{
  FLOAT *mu;
  UNSIGNED vdim;
  UNSIGNED noc;  /* Number of codebooks in the map */
  FLOAT *codebook, *sample;
  UNSIGNED n, i;
  FLOAT diffsf, diff, difference;

  vdim = gptr->dimension;
  mu = node->mu;
  noc = map->xdim * map->ydim;
  diffsf = FLT_MAX;
  sample = node->points;
  for (n = 0; n < noc; n++){  /* For every codebook of the map */
    codebook = map->codes[n].points;
    difference = 0.0;

    /* Compute the difference between codebook and input entry */
    for (i = 0; i < vdim; i++){
      diff = codebook[i] - sample[i];
      difference += diff * diff * mu[i];
      if (difference > diffsf)
	break;
    }
    /* If distance is smaller than previous distances */
    if (difference < diffsf){
      winner->codeno = n;
      diffsf         = difference;
    }
  }
  winner->diff   = diffsf;

  return;
}

/******************************************************************************
Description: 

Return value: 
******************************************************************************/
void VQFindWinnerEucledian(struct Map *map, struct Node *node, struct Graph *gptr, struct Winner *winner)
{
  FLOAT *mu;
  UNSIGNED ldim, fanout, fanin, tend;
  UNSIGNED noc;  /* Number of codebooks in the map */
  FLOAT *codebook, *sample;
  UNSIGNED n, i;
  int id;

  FLOAT diffsf, diff, difference;

  ldim = gptr->ldim;            /* Offset for label component         */
  fanout = gptr->FanOut;        /* Offset for child state component   */
  fanin = gptr->FanIn;          /* Offset for parent state component  */
  tend = ldim+2*(fanin+fanout)+gptr->tdim; /* End of target vector component */

  mu = node->mu;

  noc = map->xdim * map->ydim;
  diffsf = FLT_MAX;
  sample = node->points;
  for (n = 0; n < noc; n++){  /* For every codebook of the map */
    codebook = map->codes[n].points;
    difference = 0.0;

    /* Compute the difference between codebook and input entry label */
    for (i = 0; i < ldim; i++){
      diff = codebook[i] - sample[i];
      difference += diff * diff * mu[i];
      if (difference >= diffsf)
      	goto big_difference;
    }

    /* Consider children coordinate vector */
    diff = map->codes[n].a;
    for (i = 0; i < fanout; i++){
      id = (int)sample[ldim + i*2];
      if (id >= 0)
	diff += (1.0 - 2 * codebook[ldim+noc*i+id]);
    }
    difference += diff * mu[i-1];
    if (difference >= diffsf)
      goto big_difference;

    /* Difference to parent coordinate vector */
    diff = map->codes[n].b;
    for (i = 0; i < fanin; i++){  
      id = (int)sample[ldim + 2+ fanout + i*2];
      if (id >= 0)
	diff += (1.0 - 2 * codebook[ldim+noc*fanout+noc*i+id]);
    }
    difference += diff * mu[i-1];
    if (difference >= diffsf)
      goto big_difference;

    /* Difference to target vector component */
    for (i = ldim + 2*fanin + 2*fanout; i < tend; i++){
      diff = codebook[i] - sample[i];
      difference += diff * diff * mu[i];
      if (difference >= diffsf)
      	goto big_difference;
    }

    /* Distance is smaller than previous distances */
    winner->codeno = n;
    diffsf         = difference;
  big_difference:
    continue;
  }
  winner->diff   = diffsf;

  return;
}

/******************************************************************************
Description: Adapt all codebook vectors which are located within a fixed
             radius around the winning codebook.

Return value: none
******************************************************************************/
void BubbleAdapt(struct Graph *gptr,struct Map *map, struct Node *node, struct Winner *winner, FLOAT radius, FLOAT alpha)
{
  UNSIGNED n, noc;
  FLOAT dist;
  FLOAT (*ComputeDistance)(int bx, int by, int tx, int ty);