/**********************
 *
 * TwoTSuppress.c
 *
 * This program investigates the two tone suppression of a filter.
 * A probe pure tone is held fixed while a test tone is sweeped across
 * a range of frequencies.
 *
 **********************/

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

#include "CRLHeaders.h"

/******************************************************************************/
/****************************** Constant definitions **************************/
/******************************************************************************/

#define	PARAMETERS_FILE		"TwoTSuppress.par" /* Name of paramters file. */
#define NUM_CHANNELS		1		/* No. of filter channels. */
#define CHANNEL				0		/* Work filter channel. */
#define INTENSITY_ACCURACY	0.1		/* Accuracy of intensity points calculated*/
#define NUM_PURE_TONES		2		/* No. of pure tones in test. */
#define PROBE				0		/* Probe tone no. in lists. */
#define	TEST				1		/* Test tone no. in lists. */

/******************************************************************************/
/****************************** Global variables ******************************/
/******************************************************************************/

char	outputFile[MAXLINE], stParFile[MAXLINE], pEParFile[MAXLINE];
char	bMParFile[MAXLINE], rPParFile[MAXLINE], hCParFile[MAXLINE];
char	sGParFile[MAXLINE];

char	stModuleName[MAXLINE], pEModuleName[MAXLINE], bMModuleName[MAXLINE];
char	rPModuleName[MAXLINE], hCModuleName[MAXLINE], sGModuleName[MAXLINE];

double	lowestCFFrequency, highestCFFrequency, eRBDensity;
double	initialTestIntensity, finalTestIntensity, intensityDecrement;
double	rampInterval;

/******************************************************************************/
/****************************** Functions and subroutines *********************/
/******************************************************************************/

/****************************** ReadParsFromFile ******************************/

/*
 * This program reads a specified number of parameters from a file.
 * It expects there to be one parameter per line.
 */
 
void
ReadParsFromFile(char *fileName)
{
	FILE	*fp;
	
	if ((fp = fopen(fileName, "r")) == NULL) {
		NotifyError("ReadTestPars: Cannot open data file '%s'.\n", fileName);
		exit(1);
	}
	printf("Reading parameters from file: %s\n", fileName);
	Init_ParFile();
	GetPars_ParFile(fp, "%s", outputFile);
	GetPars_ParFile(fp, "%s %s", stParFile, stModuleName);
	GetPars_ParFile(fp, "%s %s", pEParFile, pEModuleName);
	GetPars_ParFile(fp, "%s %s", bMParFile, bMModuleName);
	GetPars_ParFile(fp, "%s %s", rPParFile, rPModuleName);
	GetPars_ParFile(fp, "%s %s", hCParFile, hCModuleName);
	GetPars_ParFile(fp, "%s %s", sGParFile, sGModuleName);
	GetPars_ParFile(fp, "%lf", &lowestCFFrequency);
	GetPars_ParFile(fp, "%lf", &highestCFFrequency);
	GetPars_ParFile(fp, "%lf", &eRBDensity);
	GetPars_ParFile(fp, "%lf", &initialTestIntensity);
	GetPars_ParFile(fp, "%lf", &finalTestIntensity);
	GetPars_ParFile(fp, "%lf", &intensityDecrement);
	GetPars_ParFile(fp, "%lf", &rampInterval);
	fclose(fp);
	Free_ParFile();
	
}

/******************************************************************************/
/****************************** Main Body *************************************/
/******************************************************************************/

int main()
{
	int			i;
	clock_t		startTime;
	double		sStateAverage[NUM_CHANNELS], probeIntensity, signalDuration;
	FILE		*fp;
	CFListPtr	frequencies;
	EarObjectPtr	stimulus = NULL, pEFilter = NULL, bMFilter = NULL;
	EarObjectPtr	recpPotn = NULL, hairCell = NULL, spikeGen = NULL;
	
	printf("Starting test program...\n\n");
	
	ReadParsFromFile(PARAMETERS_FILE);
	printf("This test routine investigates two-tone suppression for an\n");
	printf("auditory periphery model.\n");
	printf("The %s stimuli used in this investigation are ramped \n",
	  stModuleName);
	printf("with a %g ms rise time.\n", MSEC(rampInterval));
	printf("The probe tone is the first in the multi-pure tone list shown "\
	  "below.\n");
	printf("The test tone varies over the range %g - %g dB SPL, in\n",
	  initialTestIntensity, finalTestIntensity);
	printf("decrements of %g dB SPL.\n", intensityDecrement);
	
	/* GenerateERB is used to create a frequency list, which is not directly
	 * used by the filter itself. */
	frequencies = GenerateERB_CFList(lowestCFFrequency, highestCFFrequency,
	  eRBDensity);
	printf("The %s filter is tested over the range %g - %g Hz, \n",
	  bMModuleName, lowestCFFrequency, highestCFFrequency);
	printf("using an ERB scale with ERB-density %g ERBs/Hz, producing %d\n",
	  eRBDensity, frequencies->numChannels);
	printf("different frequency values.\n\n");
	printf("The model process contains the following modules:\n\n");
	printf("\tStimulus generation:\t%s\n", stModuleName);
	printf("\tOuter-/middle-ear:\t%s\n", pEModuleName);
	printf("\tBasilar membrane:\t%s\n", bMModuleName);
	printf("\tIHC receptor pot.:\t%s\n", rPModuleName);
	printf("\tInner hair cell (IHC):\t%s\n", hCModuleName);
	printf("\tAuditory nerve spike generation:\t%s\n", sGModuleName);
	printf("\n");
	   
	/* Initialising EarObjects. */
	
	if ((stimulus = Init_EarObject(stModuleName)) == NULL)
		exit(1);
	if ((pEFilter = Init_EarObject(pEModuleName)) == NULL)
		exit(1);
	if ((bMFilter = Init_EarObject(bMModuleName)) == NULL)
		exit(1);
	if ((recpPotn = Init_EarObject(rPModuleName)) == NULL)
		exit(1);
	if ((hairCell = Init_EarObject(hCModuleName)) == NULL)
		exit(1);
	if ((spikeGen = Init_EarObject(sGModuleName)) == NULL)
		exit(1);

	/* Set up EarObject connections. */
	
	ConnectOutSignalToIn_EarObject( stimulus, pEFilter );
	ConnectOutSignalToIn_EarObject( pEFilter, bMFilter );
	ConnectOutSignalToIn_EarObject( bMFilter, recpPotn );
	ConnectOutSignalToIn_EarObject( recpPotn, hairCell );
	ConnectOutSignalToIn_EarObject( hairCell, spikeGen );
	
	/* Initialising Modules */
	
	printf("Module parameters...\n\n" );

	if (!DoFun1( ReadPars, stimulus, stParFile))
		exit(1);
	DoFun( PrintPars, stimulus );

	if (!DoFun1( ReadPars, pEFilter, pEParFile))
		exit(1);
	DoFun( PrintPars, pEFilter );

	if (!DoFun1( ReadPars, bMFilter, bMParFile))
		exit(1);
	DoFun( PrintPars, bMFilter );

	if (!DoFun1( ReadPars, recpPotn, rPParFile))
		exit(1);
	DoFun( PrintPars, recpPotn );
	
	if (!DoFun1( ReadPars, hairCell, hCParFile))
		exit(1);
	DoFun( PrintPars, hairCell );

	if (!DoFun1( ReadPars, spikeGen, sGParFile))
		exit(1);
	DoFun( PrintPars, spikeGen );

	/* Start main process. */

	if ((fp = fopen(outputFile, "w")) == NULL) {
		fprintf(stderr, "TwoTSuppress: Cannot open file '%s'\n", outputFile);
		exit(1);
	}
	fprintf(fp, "dB SPL/Hz");
	for (i = 0; i < frequencies->numChannels; i++)
		fprintf(fp, "\t%6.3f", frequencies->channel[i]);
	fprintf(fp, "\n");
	
	startTime = clock();
	for (probeIntensity = initialTestIntensity; probeIntensity >=
	  finalTestIntensity; probeIntensity -= intensityDecrement) {
		DoFun2(SetIndividualIntensity, stimulus, TEST, probeIntensity);
		fprintf(fp, "%7.2lf", probeIntensity);
		printf("Intensity = %g dB SPL...\n", probeIntensity);
		for (i = 0; i < frequencies->numChannels; i++) {
			DoFun2(SetIndividualFreq, stimulus, TEST, frequencies->channel[i]);
			DoProcess(GenerateSignal, stimulus);
			if (!stimulus->outSignal->rampFlag )
				RampUpOutSignal_Ramp(stimulus, Sine_Ramp, rampInterval );
			DoProcess(RunModel, pEFilter);
			DoProcess(RunModel, bMFilter);
			DoProcess(RunModel, recpPotn);
			DoProcess(RunModel, hairCell);
			DoProcess(RunModel, spikeGen);
			signalDuration = GetDuration_SignalData(spikeGen->outSignal);
			CalcAverages_GenAnalysis(sStateAverage, spikeGen->outSignal,
			  signalDuration / 2.0, signalDuration);
			fprintf(fp, "\t%6.2lf", sStateAverage[CHANNEL] /
			  spikeGen->outSignal->dt);
		}
		fprintf(fp, "\n");
	}
	fclose(fp);
	printf("The process took %lu seconds to run.\n", (clock() -
	  startTime) / CLOCKS_PER_SEC);
	FreeAll_EarObject();
	
	return(0);
	
}