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

	IMAGE PROCESSING AND NEURAL NETWORKS LABORATORY
		The University of Texas at Arlington


PROGRAM		:	MLPTEST.C
AUTHOR		:	Hema Chandrasekaran

SYNOPSIS 	:	This program reads weights from user chosen file and applies the
				weights to a trained network.

INPUTS		:	data file name and the weight filename. The weight file is in 
				text form and begins with the topology of the MLP.

OUTPUTS		:	Outputs saved to a file, if user chooses that option.


*******************************************************************************/
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <math.h>


#define	 MAXFEATURES		100
#define  MAXHIDDENLAYERS	2
#define	 SIGDIGITS			9
#define	 FILENAMELEN		50
#define	 TRUE				1
#define	 FALSE				0



struct  UserEntry
{
	int 		Features;
	int 		DesiredOutputs;
	int 		MaxPatterns;
	char		SaveOutput;
	FILE		*fpwt, *fpfeat, *fplog, *fplist, *fpout;
	char		*WtFile, *ListFile, *FeatFile, *LogFile, *OutFile;
};

struct TopologyInfo
{
	FILE		*fpTopology;
	char		*TopologyFile;
	int			oldTopology;
	int			HiddenLayers;
	int			HiddenLayerUnits[MAXHIDDENLAYERS];
	int 		LayerUnits[MAXHIDDENLAYERS+2];
};

struct MiscArrays
{
	double 		*Inputs, *Outputs;
	double  	****Wts;
	double  	**Thresholds;
	double  	**Nets;
	double  	**LayerOutputs;
	double		MSError[MAXFEATURES];
	double		TotalError;
};

struct Origin
{
	struct UserEntry 	*userptr;
	struct TopologyInfo  *tptr;
	struct MiscArrays    *arrayptr;
};



FILE *OpenFile(char *, char *);
double  *AllocateMemory(int );
double  **AllocateMatrixMemory(int, int);
char *get_string(char *);

void main(void)
{
	struct  	UserEntry	 *userptr = NULL;
	struct  	TopologyInfo *tptr = NULL;
	struct  	MiscArrays	 *arrayptr = NULL;
	struct		Origin	*sosptr;
	double		TotalError = 0.0;
	int			Nv = 0, IsNotEof = 1;

	struct Origin *InitializeStructures(struct UserEntry *,struct TopologyInfo *,
													struct MiscArrays *);

	void 	 WriteLogandGetUserInput(struct UserEntry *);

	int	 ReadOnePattern(struct UserEntry *,	 struct MiscArrays *);

	void	 CalculateNet(struct TopologyInfo *, struct MiscArrays *);

	void	 ReadWeights(struct UserEntry *,struct TopologyInfo *,
				struct MiscArrays * );

	double	 CalculateError(struct UserEntry *,struct TopologyInfo *,
				struct MiscArrays *);

	void	 WriteErrorLog(struct UserEntry *,struct TopologyInfo *,
				struct MiscArrays *);

	sosptr = InitializeStructures(userptr,tptr, arrayptr);

	userptr = sosptr->userptr;
	tptr = sosptr->tptr ;
	arrayptr = sosptr->arrayptr;


	WriteLogandGetUserInput(userptr);
	ReadWeights(userptr, tptr, arrayptr);

	fseek(userptr->fpfeat, (long)0, 0);
	while((IsNotEof = ReadOnePattern(userptr, arrayptr)) != EOF)
	{
		CalculateNet(tptr, arrayptr);
		TotalError += CalculateError(userptr, tptr, arrayptr);
		if (IsNotEof != EOF)
		{
			Nv++;
		}
	}

	userptr->MaxPatterns = Nv;
	arrayptr->TotalError = TotalError;
	WriteErrorLog(userptr, tptr, arrayptr);
	fcloseall();
}

/***************************************************************************/
struct Origin *InitializeStructures(struct UserEntry *userptr,
								struct TopologyInfo *tptr,struct MiscArrays *arrayptr)
{
	/* declare a structure of all structures	*/
	struct Origin  *sosptr;

	sosptr = (struct  Origin*) malloc(sizeof(struct Origin));

	userptr = (struct UserEntry *) malloc(sizeof(struct UserEntry));

	tptr = (struct TopologyInfo *) malloc(sizeof(struct TopologyInfo));

	arrayptr = (struct MiscArrays *) malloc(sizeof(struct MiscArrays));


	if ((userptr == NULL)||(tptr == NULL)||(sosptr == NULL)||(arrayptr == NULL))
	{
		printf("unable to allocate memory ...\n");
		exit(1);
	}

	sosptr->userptr =	userptr ;
	sosptr->tptr = 	tptr;
	sosptr->arrayptr = arrayptr;

	userptr->ListFile = "mlp.lst";

	userptr->OutFile = (char *)calloc(FILENAMELEN,sizeof(char));
	userptr->LogFile = (char *)calloc(FILENAMELEN,sizeof(char));

	userptr->fplist = NULL;
	userptr->fpfeat = NULL;
	userptr->fplog = NULL;
	userptr->fpout = NULL;

	return(sosptr);
}
/*****************************************************************************/
void ReadWeights(struct UserEntry *userptr,
							struct TopologyInfo *tptr,struct MiscArrays *arrayptr)
{
	int i,j,k,l, TotalLayers, *LUnits, Lhid, Nh, len;
	double 	**Thresholds;
	char    *str, *str2, *str3;


	void InitializeArrays(struct UserEntry *,struct TopologyInfo *,struct MiscArrays *);
	str = (char *) calloc(FILENAMELEN, sizeof(char));
	str2 = (char *) calloc(FILENAMELEN,  sizeof(char));
	str3 = (char *) calloc(FILENAMELEN,  sizeof(char));


	rewind(userptr->fpwt);

	printf("\nReading weights from %s \n", userptr->WtFile);

	fscanf(userptr->fpwt,"%d  ", &Lhid); /*number of hidden layers first */
	tptr->HiddenLayers = Lhid;
	TotalLayers = Lhid+2;


	/* initialize	*/
	for ( i= 0; i < TotalLayers; i++)
	{
		tptr->LayerUnits[i] = 0;
	}

	for ( i= 0; i < TotalLayers-2; i++)
	{
		tptr->HiddenLayerUnits[i] = 0;
	}


	for(i =0; i < MAXFEATURES; i++)
	{
		arrayptr->MSError[i] = 0.0;
	}


	/* inpus, hidden units, outputs  */

	for ( i= 0; i < TotalLayers; i++)
	{
		fscanf(userptr->fpwt,"%d ",&tptr->LayerUnits[i]);
	}
	for ( i= 1; i < TotalLayers-1; i++)
	{
		tptr->HiddenLayerUnits[i-1] =  tptr->LayerUnits[i];

	}

	LUnits = tptr->LayerUnits;
	userptr->Features = tptr->LayerUnits[0];
	userptr->DesiredOutputs = tptr->LayerUnits[TotalLayers-1];

	InitializeArrays(userptr,tptr,arrayptr);

	Thresholds = 	arrayptr->Thresholds;

	/* restore weights	*/
	for ( i= 1; i < TotalLayers; i++)
	{
		for ( j= 0; j < i; j++)
		{
			for ( k= 0; k < LUnits[i]; k++)
			{
				for ( l = 0; l < LUnits[j]; l++)
				{
					fscanf(userptr->fpwt,"%s",str3);
					arrayptr->Wts[i][j][k][l] = atof(str3);
				}
			}
		}
	}
	/* restore thresholds too ..	*/
	for ( i= 1; i < TotalLayers; i++)
	{
		for ( j = 0; j < LUnits[i]; j++)
		{
			fscanf(userptr->fpwt,"%s",str3);
			Thresholds[i][j] = atof(str3);
		}
	}


	Nh = tptr->HiddenLayerUnits[0];
	itoa(Nh,str2,10);
	len = strlen(userptr->FeatFile);	

	/* identify data file */
	strncpy(str,userptr->FeatFile,len-4);	
	strcat(str,"_");	
	strcat(str,str2);	
	
	strcat(userptr->OutFile,str);
	strcat(userptr->OutFile,".out");
	
	strcat(userptr->LogFile,str);
	strcat(userptr->LogFile,".log");
	
	if (userptr->SaveOutput)
	{
		do
		{
			 userptr->fpout = OpenFile(userptr->OutFile,"wt");
		}	
		while (userptr->fpout == NULL);
	}
	do
	{
		 userptr->fplog = OpenFile(userptr->LogFile,"wt");
	}
	while (userptr->fplog == NULL);
	fprintf(userptr->fplog," Image Processing and Neural Networks laboratory\n");
	fprintf(userptr->fplog,"          University of Texas at Arlington \n");
	fprintf(userptr->fplog,"\nMLP Testing Program\n");

	fprintf(userptr->fplog,"\nReading weights from %s \n", userptr->WtFile);

}
/***************************************************************************
Function ReadOnePattern
***************************************************************************/
int	ReadOnePattern(struct UserEntry *userptr,struct MiscArrays *arrayptr)
{
	char str[FILENAMELEN];
	int i,j, Nf, Nmax, Nout, retval;
	double  *Inputs, *Outputs,  **Hout;
	FILE *fp;


	Nf = userptr->Features;
	Nout = userptr->DesiredOutputs;
	Nmax = Nf+Nout;
	fp = userptr->fpfeat;
	Inputs = arrayptr->Inputs;
	Outputs = arrayptr->Outputs;
	Hout = arrayptr->LayerOutputs;

	if (!feof(fp))
	{
		retval = 1;
		for ( i = 0; i < Nmax; i++)
		{
			if (fscanf(fp, "%s",str ))
			{
				if (feof(fp))
				{
					retval = EOF;
					break;
				}
				if ( i < Nf)
				{
					Inputs[i] = atof(str);
					Hout[0][i] = Inputs[i];
				}
				else
				{
					j = i-Nf;
					Outputs[j] = atof(str);
				}
			}
			else
			{
				break;
			}
		}
	}
	else
	{
		retval = EOF;
	}
	return(retval);
}
/***************************************************************************/
void	 CalculateNet(struct TopologyInfo *tptr,struct MiscArrays *arrayptr)
{
	int i,j,k,layernum,  TotalLayers,  Lhid, *LUnits, OutLayer;
	double net,   **Nets, **Thresholds, **Hout;

	LUnits = tptr->LayerUnits;
	Lhid = tptr->HiddenLayers;
	OutLayer = Lhid+1;
	TotalLayers = Lhid+2;

	Nets = arrayptr->Nets;
	Hout = arrayptr->LayerOutputs;
	Thresholds = arrayptr->Thresholds;



	/* need to reinitialize Nets[i][j] = 0.0 for each pattern	*/

	for ( i= 1; i < TotalLayers; i++)
	{
		for ( j= 0; j < LUnits[i]; j++)
		{
			Nets[i][j] = 0.0;
		}
	}

	for ( i= 1; i < TotalLayers; i++)
	{
		for ( j= 0; j < LUnits[i]; j++)
		{
			Nets[i][j] += Thresholds[i][j];

			for ( layernum = 0; layernum < i; layernum++)
			{
				for ( k = 0; k < LUnits[layernum]; k++)
				{
					Nets[i][j] += arrayptr->Wts[i][layernum][j][k]*Hout[layernum][k];
				}
			}
			/* once net function for each unit has been calculated, f(net) = output
			for that unit should be calculated and written back to the array
			to be ready for the next iteration*/
			net =  Nets[i][j] ;
			/* output layer has  linear activation	*/
			if ( i == OutLayer)
			{
				Hout[i][j] = net;
			}
			else
			{
				/* to take care of underflow, overflow problems	*/
				if (net >= 25.0)
				{
					Hout[i][j] = 1.0;
				}
				else if (net <= -25.0)
				{
					Hout[i][j] = 0.0;
				}
				else
				{
					Hout[i][j] = (1.0 /(1.0 + exp(-net)));
				}
			}
		}
	}
}
/*****************************************************************************
Function CalculateError calculates error in the k th output
***************************************************************************/
double CalculateError(struct UserEntry *userptr,struct TopologyInfo *tptr,
							struct MiscArrays *arrayptr)
{
	int 		Nout,i,Lhid,	OutLayer;
	double   *Outputs,  **Hout, error, TotalMse = 0.0;



	Nout = userptr->DesiredOutputs;
	Lhid = tptr->HiddenLayers;

	Hout = arrayptr->LayerOutputs;


	Outputs = arrayptr->Outputs;

	OutLayer = Lhid+1;

	for(i =0; i < Nout; i++)
	{
		error = Outputs[i] - Hout[OutLayer][i];
		arrayptr->MSError[i] += error*error;
		if (userptr->SaveOutput)
		{
			fprintf(userptr->fpout,"%lf ",Hout[OutLayer][i]);
		}
	}
	if (userptr->SaveOutput)
	{
		fprintf(userptr->fpout,"\n ");
	}

	for(i = 0; i < Nout; i++)
	{
		TotalMse += arrayptr->MSError[i];
	}

	return(TotalMse);
}
/***************************************************************************/
void	 WriteErrorLog(struct UserEntry *userptr,struct TopologyInfo *tptr,
				struct MiscArrays *arrayptr)
{
	int i, Nout, Nf, Nv, Layers, *LUnits;
	double Error = 0.0, TotalError = 0.0, *MSError;
	FILE *fplog;
	Layers = tptr->HiddenLayers;
	LUnits = tptr->HiddenLayerUnits;
	Nout = userptr->DesiredOutputs;
	Nf = userptr->Features;
	Nv = userptr->MaxPatterns;
	fplog = userptr->fplog;
	MSError = arrayptr->MSError;

	fprintf(fplog,"\nInput Data Filename	: %s", userptr->FeatFile);
	printf("\nInput Data Filename	: %s", userptr->FeatFile);
	fprintf(fplog,"\n Number of vectors		: %d",Nv);
	printf("\n Number of vectors		: %d",Nv);
	fprintf(fplog,"\n Number of inputs		: %d", Nf);
	printf("\n Number of inputs		: %d", Nf);
	fprintf(fplog,"\n Number of outputs		: %d", Nout);
	printf("\n Number of outputs		: %d", Nout);
	fprintf(fplog,"\n");
	fprintf(fplog,"\nNumber of Hidden Layers in the MLP :%d\n",Layers);
	printf("\nNumber of Hidden Layers in the MLP :%d\n",Layers);

	for (i = 0; i < Layers; i++)
	{
		fprintf(fplog,"\nNumber of Units in Hidden Layer # %d  : %d\n",
		i+1,LUnits[i]);
		printf("\nNumber of Units in Hidden Layer # %d  : %d\n",
		i+1,LUnits[i]);
	}

	for(i=0; i < Nout; i++)
	{
		printf("MSE for output %d = : %13.9f \t %13.9f\n", i+1,
			MSError[i],MSError[i]/Nv);

		fprintf(fplog,	"MSE for output %d = : %13.9f \t %13.9f\n", i+1,
			MSError[i],MSError[i]/Nv);
		TotalError += (MSError[i]);
		Error += (MSError[i]/Nv);
	}
	printf("\nTotal MSE            : %13.9f \t %13.9f \n", TotalError, Error);
	fprintf(fplog,"\nTotal MSE : %13.9f \t %13.9f \n", TotalError,Error);

}
/******************************************************************************/
void WriteLogandGetUserInput(struct UserEntry *userptr)
{
	int ch;
	char ans;

	
	printf("\n Image Processing and Neural Networks laboratory\n");
	printf("          University of Texas at Arlington \n");
	printf("\n       MLP Testing \n");

	do
	{
		userptr->FeatFile =  get_string("\nEnter input  data filename : ");
		userptr->fpfeat =  OpenFile(userptr->FeatFile,"rt");
	}
	while (userptr->fpfeat == NULL);

	do
	{
		userptr->WtFile =  get_string("\nEnter weights filename : ");
		userptr->fpwt =  OpenFile(userptr->WtFile,"rt");
	}
	while (userptr->fpwt == NULL);

	do
	{
		printf("\n Save the Outputs to a File? (y/n) : ");
		ch = getchar();
		ans = (char)ch;
		ch = getchar();
	}
	while(!(ans=='y'||ans!='Y'||ans=='n'||ans=='N'));
	if(ans=='y'||ans=='Y')
	{
		userptr->SaveOutput = TRUE;	
	}
	else
	{
		userptr->SaveOutput = FALSE;	
	}

}
/*****************************************************************************/
void InitializeArrays(struct UserEntry *userptr,struct TopologyInfo *tptr,
								struct MiscArrays *arrayptr)
{
	int i,j,N, TotalLayers, *LUnits, Lhid, Nf, Nout,  OutLayer;

	LUnits = tptr->LayerUnits;
	Lhid = tptr->HiddenLayers;
	OutLayer = Lhid+1;
	TotalLayers = Lhid+2;

	Nout = userptr->DesiredOutputs;
	Nf = userptr->Features;
	N = Nf+1;

	arrayptr->Inputs = AllocateMemory(N+1);

	arrayptr->Outputs = AllocateMemory(Nout);


	/* dimensions are now Nf+Nout+1(threshold)+
	Number of hidden units in all hidden layers	*/

	arrayptr->Wts  =
	(double ****) calloc(TotalLayers, sizeof( double ***));

	for ( i= 1; i < TotalLayers; i++)
	{
		arrayptr->Wts[i]  = (double ***) calloc(i, sizeof( double **));

		if ( arrayptr->Wts[i] == NULL )
		{
			printf("Unable to allocate memory for Wts[i]! \n");
		}
	}
	for ( i= 1; i < TotalLayers; i++)
	{
		for ( j= 0; j < i; j++)
		{
			arrayptr->Wts[i][j] = AllocateMatrixMemory(LUnits[i],LUnits[j]);
		}
	}

		arrayptr->Thresholds =
			(double **) calloc(TotalLayers, sizeof( double *));
		arrayptr->Nets =
			(double **) calloc(TotalLayers, sizeof( double *));
		arrayptr->LayerOutputs =
			(double **) calloc(TotalLayers, sizeof( double *));


	arrayptr->LayerOutputs[0] = AllocateMemory(LUnits[0]);

	for ( i = 1; i < TotalLayers; i++)
	{
		arrayptr->Thresholds[i] = AllocateMemory(LUnits[i]);
		arrayptr->Nets[i] = AllocateMemory(LUnits[i]);
		arrayptr->LayerOutputs[i] = AllocateMemory(LUnits[i]);
	}
}
/*****************************************************************************/

char *get_string(char *input_string)
{
	 char *str;                          

	 str = (char *) malloc(80);
	 if(!str) {
		  printf("\nString allocation error in get_string\n");
		  exit(1);
	 }
	 printf(" %s ",input_string);
	 gets(str);

	 return(str);
}
/***************************************************************************/
FILE *OpenFile(char *filename, char *mode)
{
	FILE *fp;
	/* no user interaction here ! if unable to open , return NULL pointer
		announcing error; if successful, seek beginning of opened file
		indicating success	*/
	if((fp = fopen(filename,mode))== NULL)
	{
		printf("%s file could not be opened \n", filename);
		return(NULL);
	}
	/* if successful, position the file pointer at the beginning	*/
	fseek(fp, (long)0, 0);
	return(fp);
}
/***************************************************************************/
double *AllocateMemory(int size)
{
	double  *dptr;

	dptr = (double *)calloc(size, sizeof(double));
	if(dptr == NULL)
	{
		printf("Not enough memory to allocate buffer. \n");
		printf("Aborting......\n");
		/* terminate program if out of memory 	*/
		exit(1);
	}
	return(dptr);
}
/*****************************************************************************/
double **AllocateMatrixMemory(int m, int n)
{
	double **dptr;
	int i;

	dptr = (double **)calloc(m, sizeof(double*));
	if(dptr == NULL)
	{
		printf("Not enough memory to allocate buffer. \n");
		printf("Aborting......\n");
		/* terminate program if out of memory 	*/
		exit(1);
	}
	for(i = 0; i < m ; i++)
	{
		dptr[i] = (double *)calloc(n, sizeof(double));
		if(dptr[i] == NULL)
		{
			printf("Not enough memory to allocate buffer. \n");
			printf("Aborting......\n");
			/* terminate program if out of memory 	*/
			exit(1);
		}
	}
	return(dptr);
}
/*****************************************************************************/