/******
* The program is used for training a neural network according to the art algorithm
* usage:
* a) art
*		or
* b) art -l input_file data_file weight_file
*		or
* c) art -o weight_file data_file output_file
*
* a) -- interactive user querying
* b) -- learning. input_file contains data on structure and training parameters.
*	data_file is training data. weight_file is output weights.
* c) -- output generation. weight_file is generated by learning. data_file is
*	test data. output_file will contain test data with associated class.
******/
#include <stdio.h>
#include <sys/time.h>

#define READFILE 0
#define WRITEFILE 1
#define TASK_NONE 0
#define TASK_LEARNING 1
#define TASK_OUTPUT 2
#define MAXRAND (1<<16)

struct node
{	int *t;			/* Top-down weights of the node */
	double *b;		/* Bottom-up weights of the node */
	int active;		/* Is the node active or not */
	struct node *next;	/* Link to next node */
};

struct input_s
{	int *input;		/* input vector */
	double *output;		/* output bottom-up weights */
};

struct network
{	/* The structure */
	int indim; 		/* #input dimensions */
	int num_nodes;		/* Superfluous. Can derive from list */
	struct node *list;	/* List of nodes in the network */

	/* The training parameters */
	double vigilance;	/* Vigilance parameter */
	int learning;		/* Is the network learning? */

	/* The inputs */
	int num_inputs;		/* #inputs */
	struct input_s *inputs;	/* Storage for input vectors */
};


struct network *readinput(), *readfromfile(), *readfromkeyboard();
struct network *getnetwork();
void openfile(), initrandom(), deallocate_network();
struct node *best_node(), *get_node(), *largest_output();
int update_best(), satisfactory();
void new_node(), activate_all();
void art(), trainnetwork(), get_outputs();
void writeoutputweights(), writeoutputdata();

/* For learning:
	argv[1] - "-l"
	argv[2] - input structure file
	argv[3] - training data file 
	argv[4] - output weights file
   For output generation:
	argv[1] - "-o"
	argv[2] - weights file
	argv[3] - test data file
	argv[4] - output data file   
*/
void main(argc, argv)
int argc;
char *argv[];
{
	if (argc == 1)
		art(TASK_NONE, NULL, NULL, NULL);
	else if( argv[1][0] == '-' && argv[1][1] == 'l')
		art(TASK_LEARNING,  argv[2], argv[3], argv[4]);
	else if( argv[1][0] == '-' && argv[1][1] == 'o')
		art(TASK_OUTPUT, argv[2], argv[3], argv[4]);
	else
	{	printf("Usage\n");
		printf("For interactive mode: \n");
		printf("\tart\n");
		printf("For learning\n");
		printf("\tart -l input_file data_file weights_file\n");
		printf("For output generation\n");
		printf("\tart -o weights_file data_file output_file\n");
	}
}

/* Initialize random generator with number of microseconds since
   January 1, 1970, which is almost a random value in itself.
*/
void initrandom()
{	struct timeval tp;
	struct timezone tzp;
	if( gettimeofday(&tp,&tzp) == -1 )
	{	printf("Could not initialize random number generator\n");
		_exit(0);
	}
	srandom((int)tp.tv_usec);
}


/* Sets all nodes in the list of nodes to active */
void activate_all(np)
struct node *np;
{
	while(np)
	{   np->active = 1;
	    np = np->next;
	}
}


/* Calculates a scalar product of two vectors */
double scalar_product(x, y, n)
double *x;
int *y;
int n;
{	int i;
	double sum = 0;

	for(i=0; i<n; i++)
	    sum += x[i] * y[i];

	return(sum);
}


/* Returns the node with the largest output. Only 'active' nodes are considered. 
*/
struct node *largest_output(np, input, indim)
struct node *np;
int *input;
int indim;
{	int first;
	double sum, maxsum;
	struct node *best = NULL;

	first = 1;
	while(np)
	{   if(np->active)
	    {	sum = scalar_product(np->b, input, indim);
		if( first || sum > maxsum)
		{   maxsum = sum;
		    best = np;
		    first = 0;
		}
	    }
	    np = np->next;
	}

	return(best);
}


/* Tests whether the chosen node satisfies the vigilance criterion. Returns 1 if it 
   does, and 0 if it does not. 
*/
int satisfactory(np, input, indim, vigilance)
struct node *np;
int *input;
int indim;
double vigilance;
{	int i;
	int count, sum;

	sum = count = 0;
	for(i=0; i<indim; i++)
	{   sum += input[i] * np->t[i];
	    count += input[i];
	}

	return(sum/(double)count > vigilance);
}


/* Returns a pointer to the best matching node, subject to the vigilance criterium. 
   If a matched node does not conform with the vigilance standard, it is 
   deactivated, and the next largest is found. 
*/
struct node *best_node(p, input)
struct network *p;
int *input;
{	int found = 0;
	struct node *best;

	/* All are active at the beginning */
	activate_all(p->list);
	while(!found)
	{   best = largest_output(p->list, input, p->indim);
	    if(best == NULL)
		found = 1;
	    else if(satisfactory(best, input, p->indim, p->vigilance))
		found = 1;
	    else
		best->active = 0;
	}

	return(best);
}


/* Update the node 'np', as to match the input sample 'input'. */
int update_best(np, input, indim)
struct node *np;
int *input;
int indim;
{	int i;
	int change = 0;
	int count = 0;

	for(i=0; i<indim; i++)
	{   if(np->t[i] != np->t[i] * input[i])
	    {	np->t[i] = np->t[i] * input[i];
		change = 1;
	    }
	    count += np->t[i];
	}
	for(i=0; i<indim; i++)
	    np->b[i] = np->t[i]/(0.5 + count);

	return(change);
}


/* Allocates a node. Allocates internal arrays in the node. Returns pointer to the 
   node. 
*/
struct node *get_node(indim)
int indim;
{	struct node *np;

	np = (struct node *)malloc(sizeof(struct node));
	if(np == NULL)
	    return(NULL);
	np->t = NULL;
	np->b = NULL;
	np->next = NULL;

	np->t = (int *)malloc(indim*sizeof(int));
	if(np->t == NULL)
	{   free(np);
	    return(NULL);
	}
	np->b = (double *)malloc(indim*sizeof(double));
	if(np->b == NULL)
	{   free(np->t);
	    free(np);
	    return(NULL);
	}
	np->active = 0;

	return(np);
}


/* Allocates new node. If not possible, terminates program after deallocating 
   machine. After allocating node, initializes node to the input sample 'input'. 
*/
void new_node(p, input)
struct network *p;
int *input;
{	struct node *np;
	int i, count;

	np = get_node(p->indim);
	if(np == NULL)
	{   printf("Cannot allocate new node. Out of memory\n");
	    deallocate_network(p);
	    _exit(0);
	}
	count = 0;
	for(i=0; i<p->indim; i++)
	{   np->t[i] = input[i];
	    count += np->t[i];
	}
	for(i=0; i<p->indim; i++)
	    np->b[i] = np->t[i] / (0.5 + count);
	np->next = p->list;
	p->list = np;
	p->num_nodes++;
}


/* The ART training algorithm. While changes are made to the nodes, for each input, 
   find the best matching neuron, under the vigilance constraint. If such a node 
   exists, its position is updated. If such a node does not exist, a new node is 
   created to conform with the input. 
*/
void trainnetwork(p)
struct network *p;
{	int change, i;
	struct node *best;

	change = 1;
	while (change)
	{   change = 0;
	    for(i=0; i<p->num_inputs; i++)
	    {	best = best_node(p, p->inputs[i].input);
		if(best)
		{   if(update_best(best, p->inputs[i].input, p->indim))
			change = 1;
		}
		else
		{   new_node(p, p->inputs[i].input);
		    change = 1;
		}
	    }
	}
}


/* For each output, find the winner node, and record it's bottom-up weights in the 
   structure input_s of the current input 
*/
void get_outputs(p)
struct network *p;
{	int i, j;
	struct node *best;

	activate_all(p->list);
	for(i=0; i<p->num_inputs; i++)
	{   best = largest_output(p->list, p->inputs[i].input, p->indim);
	    for(j=0; j<p->indim; j++)
		p->inputs[i].output[j] = best->b[j];
	}
}


/* Allocates the structure network, initializes all pointers to NULL, and all 
   numbers to 0
*/
struct network *getnetwork()
{	struct network *p;

	p = (struct network *)malloc(sizeof(struct network));
	if (p == NULL)
		return(NULL);

	p->indim = p->num_nodes = 0;
	p->vigilance = 0.0;
	p->list = NULL;
	p->inputs = NULL;
}


/* Deallocates all already allocated arrays and structures in the network
*/
void deallocate_network(p)
struct network *p;
{	int i;
	struct node *n;

	if(p->list)
	{   n = p->list;
	    while(n)
	    {	p->list = n->next;
		if(n->t)
		    free(n->t);
		if(n->b)
		    free(n->b);
		free(n);
		n = p->list;
	    }
	}

	if(p->inputs)
	{   for(i=0; i<p->num_inputs; i++)   
		if (p->inputs[i].input)
		    free(p->inputs[i].input);
	    free(p->inputs);
	}

	if (p)
	    free(p);
}


/* Allocate one node, read data for one node, place it in p. Return 0 if out of 
   memory, 1 if successfull
*/ 
int read_one_node(in, p)
FILE *in;
struct network *p;
{	struct node *np;
	int i;

	np = (struct node *)malloc(sizeof(struct node));
	if(!np)	
	    return(0);
	np->t = (int *)malloc(p->indim * sizeof(int));
	if(!np->t)
	{   free(np);
	    return(0);
	}
	np->b = (double *)malloc(p->indim * sizeof(double));
	if(!np->b)
	{   free(np->b);
	    free(np);
	    return(0);
	}

	for(i=0; i<p->indim; i++)
	    fscanf(in, "%d ", np->t + i);
	fscanf(in, "\n");
	for(i=0; i<p->indim; i++)
	    fscanf(in, "%lf ", np->b + i);
	fscanf(in, "\n");

	np->next = p->list;
	p->list = np;
	return(1);
}


/* If in == NULL, the user will be interactively queried for all the values. One of
   the values: the weights file, will give pout its value.
*/
struct network *readinput(task, in, data, pout)
int task;
FILE *in, *data;
FILE **pout;
{	struct network *p;

	if (task == TASK_NONE)
		p = readfromkeyboard(pout);
	else	p = readfromfile(task, in, data);

	return(p);
}


/* Prompts user for parameters of the machine, as well as names of data and output
   files. Allocates the structure network and places all values of interest in it.
   Returns the pointer towards the structure, and opens the weights file for output
   and places it's handler in pout.
*/      
struct network *readfromkeyboard(pout)
FILE **pout;
{	struct network *p;
	char buf[80];
	char ch;
	int ok, i, j;
	FILE *data, *in;

	/* Allocate the structure network */
	if( (p = getnetwork()) == NULL)
		return( NULL );

	ok = 0;
	while (!ok)
	{	printf("** Select L(earning) or O(utput generation) **\n");
		scanf("%s", buf);
		ch = buf[0];
		if( ch == 'l' || ch == 'L' )
		{	p->learning = 1;
			ok = 1;
		}
		else if( ch == 'o' || ch == 'O' )
		{	p->learning = 0;
			ok = 1;
		}
		else printf("Incorrect option entered. Please re-enter\n");
	}
	if( p->learning )
	{   printf("How many dimensions does the input pattern have?: ");
	    scanf("%d", &p->indim);
	    printf("Enter the name of the training data file: ");
	    scanf("%s", buf);  
	    openfile(buf, &data, READFILE);
	    printf("Enter the name of the output weight file: ");
	    scanf("%s", buf);
	    openfile(buf, pout, WRITEFILE);

	    printf("Vigilance parameter (number between 0 and 1): ");
	    scanf("%lf", &p->vigilance);
	    printf("Total number of input samples?: ");
	    scanf("%d", &p->num_inputs);
	}
	else
	{
	    printf("Enter the name of the weights file: ");
	    scanf("%s", buf);
	    openfile(buf, &in, READFILE);
	    printf("Enter the name of the test data file: ");
	    scanf("%s", buf);  
	    openfile(buf, &data, READFILE);
	    printf("Enter the name of the output file: ");
	    scanf("%s", buf);
	    openfile(buf, pout, WRITEFILE);

	    /* Get the number of input dimensions */	
	    fscanf(in, "%d\n", &p->indim);

	    /* Get the number of nodes */
	    fscanf(in, "%d\n", &p->num_nodes);

	    /* For each node, allocate the space, and read in the weights */
	    for(i=0; i<p->num_nodes; i++)
	    	if( !read_one_node(in, p))
		{   deallocate_network(p);
		    return(NULL);
		}

	    /* Get number of inputs */
	    fscanf(data, "%d\n", &p->num_inputs);
	}
	/* Allocate the input storage */
	p->inputs = (struct input_s *)malloc(p->num_inputs*sizeof(struct input_s));
	if( p->inputs == NULL )
	{	deallocate_network(p);
		return(NULL);
	}
	for( i=0; i<p->num_inputs; i++)
	{   p->inputs[i].input = NULL;
	    p->inputs[i].output = NULL;
	}
	for( i=0; i<p->num_inputs; i++ )
	{   p->inputs[i].input = (int *)malloc(p->indim * sizeof(int));
	    if( p->inputs[i].input == NULL )
	    {	deallocate_network(p);
		return(NULL);
	    }
	    p->inputs[i].output = (double *)malloc(p->indim * sizeof(double));
	    if( p->inputs[i].output == NULL )
	    {	deallocate_network(p);
		return(NULL);
	    }
	}
	/* Read the inputs from the data file. */
	for( i=0; i<p->num_inputs; i++)
	{   for( j=0; j<p->indim; j++)
		fscanf(data, "%d ", p->inputs[i].input + j);
	    fscanf(data, "\n");
	}

	return(p);
}


/* Reads the input files and initializes the structure networs. Places all values
   of interest in the structure. If one of the arrays in the structure cannot be
   allocated, it deallocates the whole structure and returns NULL
*/
struct network *readfromfile(task, in, data)
int task;
FILE *in, *data;
{	struct network *p = NULL;
	int i, j;

	/* Allocate the structure network */
	if( (p = getnetwork()) == NULL)
		return(NULL);

	if( task == TASK_LEARNING )
	    p->learning = 1;
	else
	    p->learning = 0;

	/* Get the # of input dimensions */
	fscanf(in, "%d\n", &p->indim);

	if( p->learning )
	{
	    /* Get the vigilance parameter */
	    fscanf(in, "%lf\n", &p->vigilance);

	    /* Get number of inputs */
	    fscanf(in, "%d\n", &p->num_inputs);
	}
	else
	{
	    /* Get the number of nodes */
	    fscanf(in, "%d\n", &p->num_nodes);

	    /* For each node, allocate the space, and read in the weights */
	    for(i=0; i<p->num_nodes; i++)
	    	if( !read_one_node(in, p))
		{   deallocate_network(p);
		    return(NULL);
		}

	    /* Get number of inputs */
	    fscanf(data, "%d\n", &p->num_inputs);
	}

	/* Allocate the input storage */
	p->inputs = (struct input_s *)malloc(p->num_inputs*sizeof(struct input_s));
	if( p->inputs == NULL )
	{	deallocate_network(p);
		return(NULL);
	}
	for( i=0; i<p->num_inputs; i++)
	{   p->inputs[i].input = NULL;
	    p->inputs[i].output = NULL;
	}
	for( i=0; i<p->num_inputs; i++ )
	{   p->inputs[i].input = (int *)malloc(p->indim * sizeof(int));
	    if( p->inputs[i].input == NULL )
	    {	deallocate_network(p);
		return(NULL);
	    }
	    p->inputs[i].output = (double *)malloc(p->indim * sizeof(double));
	    if( p->inputs[i].output == NULL )
	    {	deallocate_network(p);
		return(NULL);
	    }
	}
	/* Read the inputs from the data file. */
	for( i=0; i<p->num_inputs; i++)
	{   for( j=0; j<p->indim; j++)
		fscanf(data, "%d ", p->inputs[i].input + j);
	    fscanf(data, "\n");
	}

	return(p);
}

/* Write the weights, with the structural parameters of the LVQ network, to the file
   out
*/
void writeoutputweights(p, out)
struct network *p;
FILE *out;
{	struct node *np;
	int i;

	/* Write the number of dimensions and the number of nodes */
	fprintf(out, "%d\n", p->indim);
	fprintf(out, "%d\n", p->num_nodes);

	/* For each node, write the t-weights, then the b-weights */
	np = p->list;
	while(np)
	{   for(i=0; i<p->indim; i++)
		fprintf(out, "%d ", np->t[i]);
	    fprintf(out, "\n");

	    for(i=0; i<p->indim; i++)
		fprintf(out, "%lf ", np->b[i]);
	    fprintf(out, "\n");

	    np = np->next;
	}
}

/* For every input sample, write the input sample, followed by the output generated 
   for that sample. Output will contain both top-down and bottom-up weights for the 
   winner node 
*/
void writeoutputdata(p, out)
struct network *p;
FILE *out;
{	int i, j;

	for(i=0; i<p->num_inputs; i++)
	{   /* The input */
	    for(j=0; j<p->indim; j++)
		fprintf(out, "%d ", p->inputs[i].input[j]);
	    fprintf(out, "\n");
	    /* The top-down weights */
	    for(j=0; j<p->indim; j++)
		fprintf(out, "%d ", (p->inputs[i].output[j] != 0) );
	    fprintf(out, "\n");
	    /* The bottom-up weights */
	    for(j=0; j<p->indim; j++)
		fprintf(out, "%lf ", p->inputs[i].output[j]);
	    fprintf(out, "\n");
	}
}


/* Opens the file filename for reading or writing (depending on type).
   Terminates the program on error.
*/
void openfile( filename, fp, type)
char *filename;
FILE **fp;
int type;
{
	if( filename )
	{   if( type == READFILE )
		*fp = fopen(filename, "rt");
	    else
		*fp = fopen(filename, "wt");

	    if( *fp == NULL )
	    {	printf("Unable to open file %s for %s\n",filename,
				(type == READFILE) ? "reading" : "writing");
		_exit(0);
	    }
	}
}

/* The ART program. Opens all necessary files, reads inputs, trains the network,
   (or tests it) and writes the derived weights (or outputs).
*/
void art(task, inputfile, datafile, outfile)
int task;
char *inputfile, *datafile, *outfile;
{	FILE *in, *data, *out; 
	struct network *p;
	long t1, t2;
        
	t1 = clock();
	in = out = data = NULL;
	openfile( inputfile, &in, READFILE);
	openfile( datafile, &data, READFILE);
	openfile( outfile, &out, WRITEFILE);

	if( (p = readinput(task, in, data, &out)) == NULL)
	{	printf("Not enough memory\n");
		_exit(0);
	}
	initrandom();
	if(p->learning)
	{   trainnetwork(p);
	    writeoutputweights(p, out);
	}
	else
	{
	    get_outputs(p);
	    writeoutputdata(p, out);
	}
	deallocate_network(p);

	t2 = clock();
	printf("Elapsed time %ld milliseconds\n",(t2-t1)/1000);
}