acdcrati.c

Auditory Simulation Development Computing System version 1.5.2, is based upon a unified re-interp

/**********************
 *
 * ACDCRatio.c
 *
 * This program investigates the hair cell receptor potential model.
 * 
 *
 **********************/

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

#include "CRLHeaders.h"

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

#define	PARAMETERS_FILE		"ACDCRatio.par"/* Name of paramters file.*/
#define NUM_CHANNELS		1				/* No. of filter channels. */
#define CHANNEL				0				/* Work filter channel.  */
#define OFFSET_DURATION_SCALE	0.5	/* Point of stimulus for calculations. */

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

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

char	stModuleName[MAXLINE], pEModuleName[MAXLINE], bMModuleName[MAXLINE];

int		numberOfFreqs;
double	lowestCF, highestCF, 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)
{
	char	line[MAXLINE];
	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", rPParFile);
	GetPars_ParFile(fp, "%lf", &lowestCF);
	GetPars_ParFile(fp, "%lf", &highestCF);
	GetPars_ParFile(fp, "%d", &numberOfFreqs);
	GetPars_ParFile(fp, "%lf", &rampInterval);
	fclose(fp);
	Free_ParFile();
	
}

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

int main()
{
	int			i;
	double		acComponent, dcComponent, centreFrequency[NUM_CHANNELS];
	double		offsetTime, peakOpen, troughOpen, restingKProb;
	double		intensity[NUM_CHANNELS];
	FILE		*fp;
	ChanLen		troughs[NUM_CHANNELS], peaks[NUM_CHANNELS];
	CFListPtr	frequencies;
	EarObjectPtr	stimulus = NULL, pEFilter = NULL, bMFilter = NULL;
	EarObjectPtr	recpPotn = NULL;
	
	printf("Starting test program...\n\n");
	
	ReadParsFromFile(PARAMETERS_FILE);
	printf("This test routine investigates the hair cell receptor potential\n");
	printf("AC/DC ratio pure tone response.\n");
	printf("N.B. This program uses the receptor potential Q value from the\n");
	printf("global structure.  Because of this the receptor potential module\n");
	printf("is hard coded into the program.\n");
	printf("The stimulus is ramped with a %lg ms rise time.\n",
	  MSEC(rampInterval));
	printf("The calculation of the ac and dc components begin from\n");
	printf("time offset of %lg the stimulus duration\n", OFFSET_DURATION_SCALE);
	printf("%d frequencies are tested, they range from %lg - %lg "\
	  "Hz.\n", numberOfFreqs, lowestCF, highestCF);
	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("\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);
	recpPotn = Init_EarObject("null");

	/* Set up EarObject connections. */
	
	ConnectOutSignalToIn_EarObject( stimulus, pEFilter );
	ConnectOutSignalToIn_EarObject( pEFilter, bMFilter );
	ConnectOutSignalToIn_EarObject( bMFilter, recpPotn );

	/* Initialising Modules */
	
	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 (!ReadPars_IHC_RecptPotn(rPParFile))
		exit(1);
	PrintPars_IHC_RecptPotn();
	
	/* Start main process. */

	DoProcess(GenerateSignal, stimulus);
	if (!stimulus->outSignal->rampFlag )
		RampUpOutSignal_Ramp(stimulus, Sine_Ramp, rampInterval );
	offsetTime = stimulus->outSignal->length * stimulus->outSignal->dt *
	  OFFSET_DURATION_SCALE;
	if ((fp = fopen(outputFile, "w")) == NULL) {
		fprintf(stderr, "ACDCRatio: Cannot open file '%s'\n", outputFile);
		exit(1);
	}
	fprintf(fp, "Frequency (Hz)\tAC/DC Ratio\n");
	restingKProb = 1.0 / (1.0 + exp(recptPotn.param_q));
	frequencies = GenerateLog_CFList(lowestCF, highestCF, numberOfFreqs);
	for (i = 0; i < frequencies->numChannels; i++) {
		printf("Reading no. %d of %d...\n", i + 1, frequencies->numChannels);
		DoFun1(SetFrequency, stimulus, frequencies->channel[i]);
		DoProcess(GenerateSignal, stimulus);
		if (!stimulus->outSignal->rampFlag )
			RampUpOutSignal_Ramp(stimulus, Sine_Ramp, rampInterval );
		DoProcess(RunModel, pEFilter);
		DoProcess(RunModel, bMFilter);
		CalcIntensities_GenAnalysis(intensity, bMFilter->outSignal,
		  rampInterval);
		RunModel_IHC_RecptPotn(recpPotn);
		FindNextMinIndexes_GenAnalysis(troughs, recpPotn->outSignal,
		  offsetTime);
		FindNextPeakIndexes_GenAnalysis(peaks, recpPotn->outSignal,
		  offsetTime);
		peakOpen = recpPotn->outSignal->channel[CHANNEL][peaks[CHANNEL]];
		troughOpen = recpPotn->outSignal->channel[CHANNEL][troughs[CHANNEL]];
		acComponent = peakOpen - troughOpen;
		dcComponent = (peakOpen + troughOpen) / 2.0 - restingKProb;
		fprintf(fp, "%10.5lf\t%6.2lf\t%lg\n", frequencies->channel[i],
		  acComponent / dcComponent, intensity[CHANNEL]);
	}
	fclose(fp);
	
	FreeAll_EarObject();
	
	return(0);
	
}