Mercurial > hg > aim92
view model/image.c @ 0:5242703e91d3 tip
Initial checkin for AIM92 aimR8.2 (last updated May 1997).
| author | tomwalters |
|---|---|
| date | Fri, 20 May 2011 15:19:45 +0100 |
| parents | |
| children |
line wrap: on
line source
/* SCCS VERSION 1.29 WITH 1.28 STRIPPED OUT. Copyright (c) Applied Psychology Unit, Medical Research Council. 1988, 1989 =========================================================================== Permission to use, copy, modify, and distribute this software without fee is hereby granted for research purposes, provided that this copyright notice appears in all copies and in all supporting documentation, and that the software is not redistributed for any fee (except for a nominal shipping charge). Anyone wanting to incorporate all or part of this software in a commercial product must obtain a license from the Medical Research Council. The MRC makes no representations about the suitability of this software for any purpose. It is provided "as is" without express or implied warranty. */ /*************************************************************************** * image.c Version with trigger at right end of sai. * ======= * Input data derived from a cochleagram is used to construct a stabilised * auditory image (sai) of size `chans' * `framewidth' points. * * Organisation of data. * Data is 2-dimensional: in time (rows) and frequency (columns). * * Input data is ordered by columns. The input array consists of blocks (columns) * of `chans' points, each column being one time sample of `chans' channels. * The zero'th point in each column is the lowest frequency channel. * * Output (sai) data is ordered by rows. The output array consists of blocks * (rows) of `framewidth' points. * The zero'th point in each row is the "oldest" time, as the time origin is * on the right edge of the sai. * * Each sai is called a `frame', and its dimensions are `chans' rows by * `framewidth' columns. In the code, the `framewidth' is often represented * as imagewidth/chans, where imagewidth is the total amount of data in an * sai frame, (ie. imagewidth=chans*framewidth). * This module has been modified to have 5 levels of strobing. A Jay Datta 6/03/95 ***************************************************************************/ #include <math.h> #include <string.h> #include "options.h" #include "stitch.h" #include "source.h" #include "model.h" #include "units.h" #include "image.h" #include "calc.h" #define MAX_BUFFER (1l<<15) #ifdef FLOAT #define INPUT_SCALE 0.25 #define DECAY_SCALE 0.875 #else #define INPUT_SHIFT 2 #define DECAY_ROUND 7 #endif #define DECAY_SHIFT 3 #define MAXCHAN 1000 #ifndef lint static char *sccs_id = "@(#)image.c 1.27 J. Holdsworth, M. Allerhand (MRC-APU) 6/6/91" ; #endif #ifdef DSP32 static int segsize = 16 ; #else static int segsize = 50 ; #endif typedef enum { False, True } Bool ; typedef struct Node{int index; short val; struct Node *next;} node; static void doColumn(); static Pointer sai_callback(); static Pointer summary_callback(); static void save_callback(); static void addIn(); static void decayImage(); void output_simple_strobe_info(); void doList_strobe_info(); void output_thresh_info(); void decay_strobe_threshold(); void print_trigger_debugging_info(); void initlist(); node *getnode(); node *insertl(); void inss(); FILE *trigger_file; FILE *infotxt_file; FILE *thresho_file; int locmax_searchtime[MAXCHAN]; int locmax_searchstart[MAXCHAN]; int val[MAXCHAN]; int ltime[MAXCHAN]; int initial_strobe_candidate_time[MAXCHAN]; int ttime[MAXCHAN]; static node *start, *endl, *pp, *qq; static char *exposwitch; #ifdef FLOAT float *S; #else short *S; #endif /*************************************************************************** * struct _sai_state: * The source->info->state structure which is initialised by the * routine Sai, and is then the state argument for the callback * function: sai_callback. ***************************************************************************/ struct _sai_state { struct _pullable_source parent ; Source source; void (*trigger)(); /* trigger algorithm (originally doColumn) */ unsigned chans; /* number of filter channels */ unsigned framestep; /* sai display update period */ DataType *image; /* sai array (DataType = short) */ unsigned imagewidth; /* length of sai array (framewidth*chans) */ unsigned imagenwidth; /* portion of sai array for transient info */ unsigned decay_time; /* a factor of sai decay time constant */ unsigned decay_count; ScalarType *framedecay; /* array of factors for cochleagram frame */ int *cps; /* array of 1/centre-freq (in samples) */ int *isPulse; /* array of bools, to see if data > thres */ int Stcrit; /* int restriction */ int SusThresh; /* Level 2 user settable thresh */ ScalarType triggerdecay; /* prop/point decay for strobe threshold */ int *trigtime; /* time of last trigger pulse */ int *trigheight; /* height of last trigger pulse */ int *isStrobe; /* array of booleans, set when pulse peak found */ float *thresh; /* array of thresholds, one per channel */ float *tlim; /* linear decay value */ int time; /* keep track of sai */ int *previnput; /* lock till new pulse found */ int *def_strobe_candidate_time; int *def_strobe_height; char *switch_info; int tlim1; int stlag; /* strobe nwid first attempt to decouple nwids */ }; /*************************************************************************** * Sai: * Sai is called from SaiEntry (in module model.c). * Routines in this module it uses are: * doColumn, (which uses: addIn) * sai_callback, (which uses: decayImage) * * SaiEntry (in model.c) is the entry-point function to the stabilised image * processing module. Its purpose is to prepare the arguments for Sai by * converting parameter-option strings to integers, and then to call Sai. * * Sai is called with the following arguments (in model.c): * Sai( source, Frameheight(), segment, width, * (int) Samples( decaystr )/Framestep(), (int) Freq( samplestr ), * frequencies, Samples(ttstr)/Framestep(), Samples(cgmstr)/Framestep(), * (int) Freq("250Hz"), (int) Freq("20Hz") ) ; * * The arguments are as follows: * (The option name refers to the description given by "help". The values etc * can be found in table.c as described in model.docs). * * Argument Option name Default Comment * -------- ----------- ------- ------- * chans (see below) Number of filter channels * framestep segment_sai 16ms SAI display update period * framewidth duration_sai 32ms SAI duration (ms) * decay decayrate_sai 32ms SAI decay time constant (SU/ms) 15ms NOW Default * samplerate samplerate 20000. Input wave sample rate (Hz) * cfs (see below) An array of filter centre frequencies * cgmdecay napdecay_sai 16ms NAP (ie cochleagram) decay time constant * * Both "chans" and "cfs" are setup from the routine "updateFrequencies" in * model.c, because both arguments are derived from three basic parameters: * * minstr mincf_afb 220Hz Minimum center frequency (Hz) * maxstr maxcf_afb 4400Hz Maximum center frequency (Hz) * denstr dencf_afb 4. Filter density (filters/critical band) * * The "chans" argument is calculated in the routine "NumberCenterFrequencies", * and the result is copied into the "frameheightstr" by a call to the routine * "setFrameheight". * The "cfs" argument is calculated in the routine "GenerateCenterFrequencies". * (Both CenterFrequency routines are in gamma_tone.c). * ***************************************************************************/ Source Sai(source, chans, framestep, pwidth, nwidth, decay, samplerate, cfs, ttdecay, cgmdecay, tlim1, tlim2, suslevel, susthresh, switchinfo, expoSwitch,stlag) Source source ; int chans, framestep, pwidth, nwidth, decay ; int samplerate; double *cfs ; double ttdecay, cgmdecay ; int tlim1, tlim2; int suslevel; char *susthresh; char *switchinfo; char *expoSwitch; int stlag; { DeclareNew( struct _sai_state *, state ) ; int i, chan, startup, segment, framewidth ; double pts_per_ms, lin_napdec, lin_napfac ; /* roy 11-8-92 roy */ char *ThreshoFile, *InfotxtFile, *TriggerFile; exposwitch=(char *)calloc(5, sizeof(char)); exposwitch=expoSwitch; ThreshoFile = (char *)calloc(256, sizeof(char)); InfotxtFile = (char *)calloc(256, sizeof(char)); TriggerFile = (char *)calloc(256, sizeof(char)); state->switch_info = (char *)calloc(256,sizeof(char)); state->switch_info = switchinfo; if (!strcmp(susthresh, "on")) state->SusThresh=0; else if (!strcmp(susthresh, "off")) state->SusThresh=1; else state->SusThresh=atoi(susthresh); if (!strcmp(state->switch_info, "off")) ; else if (!strcmp(state->switch_info, "on")) infotxt_file=stdout; else{ strcpy(InfotxtFile, switchinfo); strcat(InfotxtFile, ".info"); infotxt_file=fopen(InfotxtFile, "w"); strcpy(TriggerFile, switchinfo); strcat(TriggerFile, ".trigger"); trigger_file=fopen(TriggerFile, "w"); strcpy(ThreshoFile, switchinfo); strcat(ThreshoFile, ".thresh"); thresho_file=fopen(ThreshoFile, "w"); } initlist(); state->tlim1=0; framewidth = pwidth+nwidth; pts_per_ms = (double)samplerate/1000 ; state->stlag=stlag; /* NEW line 95April */ /* Initialise members of the sai state structure */ state->source = NewRetainingSource( source, sizeof ( DataType ) * framewidth * chans ) ; state->chans = chans ; state->trigger = doColumn ; state->framestep = framestep * chans ; state->imagewidth = framewidth * chans ; /* size of sai buffer */ state->imagenwidth = nwidth * chans ; /* size of buffer for transient info */ state->decay_time = decay * -log( 1. - 1. / ( 1 << DECAY_SHIFT ) ) ; state->decay_count = state->decay_time ; state->Stcrit = suslevel; if (stlag==0 && state->Stcrit==5) /* was nwid NEW */ state->Stcrit = suslevel-1; /* state->SusThresh = susthresh; */ if (!strcmp(exposwitch, "on")) fprintf(stderr, "Exponential trigger decay SWITCH is ON\n"); /* if (!strcmp(exposwitch, "off")) printf("SWITCH is OFF\n"); */ if (!strcmp(exposwitch, "on")) state->triggerdecay = (double) 1-((ttdecay / pts_per_ms ) /100.0) ; else state->triggerdecay = (double) (ttdecay /pts_per_ms) / 100.0 ; if (!strcmp(switchinfo, "off")) ; else fprintf(infotxt_file, "\ntriggerdecay = %15e\n", state->triggerdecay); state->isStrobe = NewZeroedArray( int, state->chans, "image.c for Strobe" ) ; state->isPulse = NewZeroedArray( int, state->chans, "image.c for trigger pulse" ) ; state->trigtime = NewZeroedArray( int, state->chans, "image.c for trigger time" ) ; state->trigheight = NewZeroedArray( int, state->chans, "image.c for trigger height" ) ; state->thresh = NewZeroedArray( float, state->chans, "image.c for thresh" ) ; state->tlim = NewZeroedArray( float, state->chans, "image.c for decay" ) ; state->time = 0; /* Time initialized to Zero */ state->previnput = NewZeroedArray( int, state->chans, "image.c for pulse lock" ) ; state->def_strobe_candidate_time = NewZeroedArray( int, state->chans, "image.c for prevtrig time" ); state->def_strobe_height = NewZeroedArray( int, state->chans, "image.c for prev trig" ); for (chan=0; chan < MAXCHAN; chan++) locmax_searchtime[chan]=0; for (chan=0; chan < MAXCHAN; chan++) locmax_searchstart[chan]=0; for (chan=0; chan < MAXCHAN; chan++) val[chan]=0; for (chan=0; chan < MAXCHAN; chan++) initial_strobe_candidate_time[chan]=0; for (chan=0; chan < MAXCHAN; chan++) ltime[chan]=0; for (chan=0; chan < MAXCHAN; chan++) ttime[chan]=0; #ifdef FLOAT S=(float *)calloc(state->chans, sizeof(float)); #else S=(short *)calloc(state->chans, sizeof(short)); #endif /* Declare new arrays for members of state structure */ state->framedecay = NewArray( ScalarType, framewidth, "image.c for decay" ) ; state->image = NewZeroedArray( DataType, state->imagewidth, "image.c for image" ) ; state->cps = NewZeroedArray( int, state->chans, "image.c for channel centre periods" ) ; /* Initialise cochleagram frame decay factors */ if (cgmdecay > 0) { /* 11111 roy 11-8-92 These mods change to a linear nap decay. 1111111 roy */ /* It requires pts_per_ms, lin_napdec and lin_napfac declared just after DeclareNew above */ lin_napdec = ((cgmdecay/100)/pts_per_ms)/pts_per_ms ; /* roy 19-8 roy */ /* cgmdecay is div by pts_per_ms to get back to napdecay in the original units */ /* When napdecay is interpreted as the % to decay per ms, and since we are */ /* operating in points, we need napdec/pts_per_ms, as the dec for each pt within the ms */ /* The decay vector is set to scale the NAP by 1.0 at strobe point. So, in the sai, */ /* it falls to the left AND RISES to the right of 0 ms. **** roy **** */ lin_napfac = 1.0 + nwidth*lin_napdec ; for (i=0 ; i < framewidth ; i++) { if ( lin_napfac > 0 ) state->framedecay[i] = SCALAR( lin_napfac ); else state->framedecay[i] = SCALAR( 0 ); lin_napfac -= lin_napdec ; } } else /* disable decay factor on zero argument */ for (i=0 ; i < framewidth ; i++) state->framedecay[i] = SCALAR( 1.0 ) ; /* Initialise centre-period for each channel */ for (i=0 ; i < state->chans ; i++) state->cps[i] = samplerate / cfs[i]; #if defined( PC ) || defined( THINK_C ) if( (long) ( MAX_BUFFER - (long) state->imagewidth * sizeof ( DataType ) ) < (long) ( state->chans * sizeof ( DataType ) ) ) stitch_error( "Sorry, image larger than maximum buffer size\n" ) ; #endif /* for (i=0; i<state->imagewidth;i++) fprintf(stderr, "Image Buffer is %f, count = %d\n", *(state->image+i), (i+1)); */ return ( SetPullableSource( state, sai_callback, "image.c stabilised image" ) ) ; } /**************************************************************************** * sai_callback * The callback function at source->info->callback in the source * returned by Sai. * A "frame" is a window over the cochleagram, of size chans * points. * A "segment" is a stretch of the cochleagram which is segsize (ie 50) time * samples long, ie segsize columns of the cochleagram. * * This callback routine is executed once for each sai frame which is plotted. * "*bytes" is the number of bytes in the sai, ie: * *bytes = state->imagewidth * sizeof(DataType) ****************************************************************************/ static Pointer sai_callback( state, bytes ) struct _sai_state *state ; ByteCount *bytes ; { register long col, point ; register long points = state->framestep; /* total num points under frame */ #if defined( PC ) int segment = ( MAX_BUFFER - state->imagewidth * sizeof ( DataType ) ) / state->chans / sizeof ( DataType ) * state->chans ; #else int segment = segsize * state->chans ; #endif register DataType *input ; /* cochleagram (input) data array */ static int first=1; static int test=0; int chan; int i; /* #if defined( PC ) test=1; #else if (state->SusThresh == 0) segment= 30 * state->chans ; #endif */ /* test++; */ /* If size argument is zero, then by convention free space and return null */ if( *bytes == 0 ) { Pull( state->source, 0 ) ; Delete( state->image ) ; return ( DeletePullableSource( state ) ) ; } /* Initially pull enough data to buffer transient info ahead of the trigger */ /* Note, the pull operation is not valid when the size arg is zero. */ /* The result would be a segmentation fault on the next pull, below. */ if (state->Stcrit==0 || state->Stcrit==1 || state->Stcrit==2 || state->Stcrit==3) if (first) { if (state->imagenwidth > 0) input = PullItems(state->source, state->imagenwidth, DataType); first = 0; } /* Search frame in segment blocks for trigger-points */ for( point=0 ; point < points ; point += segment ) { if( segment > points - point ) /* Finish at end of frame */ segment = points - point ; /* Get one segment of the input data */ /* Pull 50 columns of data from the source, in buffer *input. */ /* Keep one sai of data, from last pull, in buffer behind *input */ input = PullItems(state->source, segment, DataType); for( col=0 ; col < segment ; col += state->chans ) { /* Retard the "current" data point by transient buffer width */ /* (This retards the trigger by the transient time, about 5ms) */ if (state->Stcrit==0 || state->Stcrit==1 || state->Stcrit==2 || state->Stcrit==3) state->trigger(input-state->imagenwidth, state) ; else{ state->trigger( input, state) ; /* -state->imagenwidth, state ) ; */ /*fprintf(stderr, " %d ", points) ; */ } input += state->chans ; /* next column */ /* decay image periodically, once every `decay_time' columns */ /* Lines moved down */ if( --state->decay_count <= 0 ) { decayImage( state ) ; state->decay_count = state->decay_time ; } } /* for each col */ } /* for each seg */ /* whatever the number of bytes requested, the image is always this size */ *bytes = state->imagewidth * sizeof ( *state->image ) ; /* if (state->Stcrit==1 && state->SusThresh==0){ if (test%5==0) decayImage( state ); test++;} else if (state->Stcrit==1 && state->SusThresh>0){ changed 2 to 3 if (test%10==0) decayImage( state ); test++;} else{ */ return ( (Pointer) state->image ) ; } /**************************************************************************** * decayImage * Periodically attenuate image to give computationally efficient image decay. * Attenuate all points in the sai by a constant factor. ****************************************************************************/ static void decayImage( state ) struct _sai_state *state ; { register DataType *image_ptr, *image_end ; #ifdef FLOAT register FLOAT image_decay = DECAY_SCALE ; #endif image_end = state->image + state->imagewidth ; /* end of sai array */ #ifdef FLOAT for( image_ptr=state->image ; image_ptr < image_end ; ) *image_ptr++ *= image_decay ; #else for( image_ptr=state->image ; image_ptr < image_end ; image_ptr++ ) *image_ptr -= *image_ptr + DECAY_ROUND >> DECAY_SHIFT ; #endif return ; } /**************************************************************************** * doColumn * Trigger algorithm for Sai. * Installed as source->info->state->trigger. * ****************************************************************************/ static void doColumn( input, state ) DataType *input ; struct _sai_state *state ; { register int chan ; float strobelag; int s_lag; short *STI; static int check=0; int storev; char buf[2]; short stipts=-2000; short zeroval=0; STI=(short *)calloc(state->chans, sizeof(short)); strobelag = (double) (state->time)-(state->stlag); s_lag=(state->time)-(state->stlag); /* For each channel up column, add a row of data to sai if the */ /* current point is non-zero. Add the row from the current point */ /* plus state->imagenwidth, to display the transient (nwidth) info. */ /* 000000000000000000000000 stcrit 0000000000000000000000000000000000 */ /* Add on every point : Low Pass Filter */ if (state->Stcrit==0) { for( chan=0 ; chan < state->chans ; chan++ ) { addIn(input+(state->imagenwidth), chan, state); output_simple_strobe_info(state, stipts); } } /* if (state->Stcrit==0) */ /* 111111111111111111111111 stcrit 1111111111111111111111111111111111 */ /* Add on every point greater than "stthresh_ai" (default 0) */ if (state->Stcrit==1) { for( chan=0 ; chan < state->chans ; chan++ ) { if (input[chan] > state->SusThresh ) { addIn(input+(state->imagenwidth), chan, state); output_simple_strobe_info(state, stipts); } else if (input[chan] <=state->SusThresh) output_simple_strobe_info(state, zeroval); } } /* if (state->Stcrit==1) */ /* 222222222222222222222222 stcrit 2222222222222222222222222222222222 */ /* Add on every Nap pulse peak */ if (state->Stcrit==2) { for( chan=0 ; chan < state->chans ; chan++ ) { if (state->isStrobe[chan]) /* Check to see if it is time to addIn */ { addIn(input+((state->imagenwidth)-2*(state->chans)), chan, state); state->isStrobe[chan]=0; /* Reset strobe flag */ doList_strobe_info(state, chan); } /* Find a pulse peak routine */ if (!state->isPulse[chan]) /* Not in pulse, looking for pulse */ { state->thresh[chan]=input[chan]; if (input[chan] > 0 /*state->thresh[chan] */ ) /* Start of a Nap pulse */ state->isPulse[chan] =1; /* set in pulse flag */ } else /* In pulse, looking for peak */ { if (input[chan] < state->thresh[chan]) /* peak found */ { state->isPulse[chan]=0; state->isStrobe[chan]=1; /* set strobe pending flag */ } else state->thresh[chan]=input[chan]; } output_thresh_info(state, chan); } } /* if (state->stcrit==2) */ /* 333333333333333333333333 stcrit 3333333333333333333333333333333333 */ /* Add on pulse peaks that exceed strobe threshold (temporal shadow ) */ /* Primitive Local Max algorithm */ if (state->Stcrit==3) { for( chan=0 ; chan < state->chans ; chan++ ) { if (state->isStrobe[chan]) /* Check to see if it is time to addIn */ { addIn(input+((state->imagenwidth)-2*(state->chans)), chan, state); state->isStrobe[chan] = 0; /* Reset strobe flag */ doList_strobe_info(state, chan); } /* Find a pulse peak routine */ if (!state->isPulse[chan]) /* Not in Pulse, looking for Pulse */ { if (input[chan] > state->thresh[chan] && input[chan] > state->previnput[chan]) { /* Previnput check ensures threshold decay does not cut into pulse */ state->isPulse[chan] = 1; /* set flag: in pulse */ state->thresh[chan] = input[chan]; /* start search for peak */ } else { decay_strobe_threshold(state, chan); if (input[chan]==0) /* in the zeroes between nap pulses */ state->previnput[chan] = 0; /* Set it to Zero previnput */ } } if (state->isPulse[chan]) /* Else in ongoing pulse looking for peak */ { if (input[chan] <= state->thresh[chan] && input[chan] <= state->previnput[chan]) /* pulse peak found */ { state->isPulse[chan] = 0; /* clear flag: not in pulse */ state->isStrobe[chan] = 1; /* set Strobe pending flag */ if (!strcmp(exposwitch, "off")) state->tlim[chan]=state->triggerdecay * state->thresh[chan]; state->trigtime[chan]=(state->time)-1; /* state->previnput[chan]=input[cha */ } else /* still in pulse, so continue search for peak */ state->thresh[chan] = input[chan]; } /* if in Pulse */ state->previnput[chan]=input[chan]; output_thresh_info(state, chan); } /* for all chans... */ } /* if (state->Stcrit==3) */ /* 444444444444444444444444 stcrit 4444444444444444444444444444444444 */ /* Add on pulse peaks that exceed strobe threshold and which are not succeeded by a larger pulse in nwid ms. (Better Local Max mechanism) */ if (state->Stcrit==4) { for( chan=0 ; chan < state->chans ; chan++ ) { if (state->isStrobe[chan]) /* Check to see if it is time to addIn */ { if ( strobelag-1 >= state->trigtime[chan]) /* initiates strobe after nwid ms */ { addIn(input-(state->chans), chan, state); doList_strobe_info(state, chan); state->isStrobe[chan] = 0; /* Reset strobe flag */ state->trigheight[chan] = 0; /* Reset local max value */ } } /* Find a pulse peak routine */ if (!state->isPulse[chan]) /* Not in Pulse, looking for pulse */ { if (input[chan] > state->thresh[chan] && input[chan] > state->previnput[chan]) { /* Previnput check ensures threshold decay does not cut into pulse */ state->isPulse[chan] = 1; /* set flag: in pulse */ state->thresh[chan] = input[chan]; /* start search for peak */ /* state->previnput[chan]= 1; in pulse till the next zero value */ } else { decay_strobe_threshold(state, chan); if (input[chan]==0) /* in the zeroes between nap pulses */ state->previnput[chan] = 0; /* Free previnput */ } } if (state->isPulse[chan]) /* Else in ongoing pulse, looking for peak */ { if (input[chan] <= state->thresh[chan] && input[chan] <= state->previnput[chan]) /* pulse peak found */ { state->isPulse[chan] = 0; /* clear flag: not in pulse */ state->isStrobe[chan] = 1; /* set Strobe pending flag */ if ( state->thresh[chan] > state->trigheight[chan] ) /* Check for LOCAL MAX */ { state->trigheight[chan] = state->thresh[chan]; state->trigtime[chan] = (state->time)-1; /* Set strobe time */ if (!strcmp(exposwitch, "off")) state->tlim[chan]=state->triggerdecay * state->trigheight[chan]; } } else /* still in pulse, so continue search for peak */ state->thresh[chan] = input[chan]; } /* if in Pulse */ state->previnput[chan]=input[chan]; output_thresh_info(state, chan); } /* for all chans... */ } /* if (state->Stcrit==4) */ /* 555555555555555555555555 stcrit 5555555555555555555555555555555555 */ /* Add on pulse peaks that exceed storbe threshold and which are not succeeded by a larger pulse in nwid ms, provided total lag is < 2*nwid ms Local Max mechanism with a timeout */ if (state->Stcrit==5) { for (chan = 0; chan < state->chans; chan++) { if(state->isStrobe[chan]) /* Check to see if it is time to addIn */ { /* check for waiting till the first local max found */ if (strobelag-1 >= state->def_strobe_candidate_time[chan]) /* initiate strobe after nwid ms */ { /* print_trigger_debugging_info(state, chan, s_lag, 1); */ addIn(input-state->chans, chan, state); doList_strobe_info(state, chan); state->def_strobe_candidate_time[chan] = 0; /* Reset strobe time */ state->def_strobe_height[chan] = 0; /* Reset local max value */ state->trigheight[chan]=0; state->isStrobe[chan]=0; /* Reset strobe flag */ } } /* Find a pulse peak routine */ if (!state->isPulse[chan]) /* Not in pulse, looking for pulse */ { if (input[chan] > state->thresh[chan] && input[chan] > state->previnput[chan]) { /* Previnput check ensures threshold decay does not cut into pulse */ state->isPulse[chan] = 1; /* set flag: in pulse */ state->thresh[chan] = input[chan]; /* start search for peak */ /* state->previnput[chan] = 1; in pulse till the next zero value */ } else { decay_strobe_threshold(state, chan); if (input[chan]==0) /* in the zeroes between nap pulses */ state->previnput[chan]=0; /* Free previnput */ } } if (state->isPulse[chan]) /* Else in ongoing pulse, looking for peak */ { if (input[chan] <= state->thresh[chan] && input[chan] <= state->previnput[chan]) /* Pulse Peak Found */ { state->isPulse[chan]=0; /* clear flag: not in pulse */ if (locmax_searchstart[chan]==0 ) /* && state->time>=((ttime[chan]-initial_strobe_candidate_time[chan])+ltime[chan])) AJD 1-2-96 */ { /* start local max time (for timeout) and shift search window */ locmax_searchtime[chan]=1; /* local max time started */ locmax_searchstart[chan]=1; /* local max flag set */ print_trigger_debugging_info(state, chan, s_lag, 2); initial_strobe_candidate_time[chan]=(state->time)-1; /* local max start time noted */ } /* if (locmax_searchstart[chan]==1) */ /* 8rd March, 1995 AJayDatta */ if (state->thresh[chan] > state->trigheight[chan]) /* Check for LOCAL MAX */ { state->trigheight[chan]=state->thresh[chan]; state->trigtime[chan]=(state->time)-1; /* Set strobe time */ if (!strcmp(exposwitch, "off")) state->tlim[chan]=state->triggerdecay * state->trigheight[chan]; if (locmax_searchstart[chan]==1) print_trigger_debugging_info(state, chan, s_lag, 3); } /* if greater than thresh */ /* Reached end of time-out period (While in the peak finding stage !) End of the timeout period can occur at two stages of the algorithm, at a Nap pulse peak or at any point of the pulse train. */ if (locmax_searchtime[chan]==(state->stlag)) /* imagenwidth/state->chans)) */ { print_trigger_debugging_info(state, chan, s_lag, 4); state->def_strobe_candidate_time[chan]=state->trigtime[chan]; state->def_strobe_height[chan]=state->trigheight[chan]; locmax_searchtime[chan]=0; locmax_searchstart[chan]=0; state->trigheight[chan]=0; state->trigtime[chan]=0; state->isStrobe[chan]=1; /* set Strobe pending flag */ ltime[chan]=state->time; /* End of search time */ ttime[chan]=state->def_strobe_candidate_time[chan]; /* Local Max time */ } /* reset search period to zero */ /* as y ms search period ends, update the info for strobe and clear current_ values for the next y ms */ } /* Pulse Peak Loop */ else state->thresh[chan]=input[chan]; } /* if (state->isPulse..) */ if (locmax_searchtime[chan] == (state->stlag)) /* imagenwidth/state->chans)) */ /* && state->def_strobe_candidate_time[chan]==0) */ { /* start local max time (for timeout) */ print_trigger_debugging_info(state, chan, s_lag, 5); state->def_strobe_candidate_time[chan]=state->trigtime[chan]; state->def_strobe_height[chan]=state->trigheight[chan]; locmax_searchtime[chan]=0; locmax_searchstart[chan]=0; state->trigheight[chan]=0; state->trigtime[chan]=0; state->isStrobe[chan]=1; /* set Strobe pending flag */ ltime[chan]=state->time; ttime[chan]=state->def_strobe_candidate_time[chan]; } /* reset search period to zero */ /* as y ms search period ends, update the info for strobe and clear current_ values for the next y ms */ if (locmax_searchstart[chan]==1) (locmax_searchtime[chan])++; state->previnput[chan]=input[chan]; output_thresh_info(state, chan); } /* for all chans */ } /* if (state->Stcrit== 5) */ qq=insertl(state->time); /* insert time to list */ (state->time)++; return; } /* doColumn */ /**************************************************************************** * addIn * add row cochleagram into stabilised image ****************************************************************************/ static void addIn( input, chan, state ) DataType *input ; int chan ; struct _sai_state *state ; { register DataType *channel_ptr, *image_ptr, *end; #ifdef FLOAT register FLOAT input_scale = INPUT_SCALE ; #endif register ScalarType *decayfactor = state->framedecay; /* Initialize channel_ptr to the input data to be added in. */ /* `input' points to a particular column in the cochleagram, */ /* and `chan' indexes a particular channel in this column. */ channel_ptr = input + chan ; end = channel_ptr - state->imagewidth ; /* Initialize image_ptr to end of sai row corresponding to `chan'.*/ /* `state->image' points to the start of the sai, */ /* so increment pointer by chan*framewidth to get to row `chan'. */ /* To get to end of row, increment by ((chan+1)*framewidth)-1. */ /* But framewidth = state->imagewidth / state->chans, and hence: */ /* Added fix; as new spiral.c accepts nwid AJD 17-3-95 */ if (state->Stcrit==4 && state->stlag==0) image_ptr = (state->image + state->imagewidth / state->chans * (chan+1)); /* 13-3-95 removed -1; ajd */ else image_ptr = (state->image + state->imagewidth / state->chans * (chan+1))-1; /* Decrement channel_ptr by columns, from the initial column, */ /* until one complete sai row has been added into, which is when */ /* channel_ptr has decreased by a whole imagewidth. */ while( channel_ptr > end ) { #ifdef FLOAT *image_ptr += ( ( *channel_ptr * input_scale ) * (*decayfactor++) ) ; #else *image_ptr += DESCALE( ( *channel_ptr >> INPUT_SHIFT ) * (*decayfactor++) ) ; #endif /* if (state->chans <= 40) */ if (*image_ptr > 32767.0) *image_ptr = 32767.0; /* 32767 */ image_ptr--; channel_ptr -= state->chans ; /* next (ie previous) input time point */ } return ; } /*************************************************************************** * Summary is not called in this module. * (see SummaryEntry in model.c, which is the entry point for gensas). * Summary computes a row summary spectrogram and is called during the * program "gensas". * Routines in this module it uses are: * summary_callback ****************************************************************************/ struct _summary_state { struct _fillable_source parent ; Source source; int framewidth; /* number of time-points or cols in an sai */ int frameheight; /* number of channels or rows in an sai */ double scale; int *llim; /* integration limits (lower & upper) for sai summary data */ int *ulim; }; Source Summary( source, frameheight, scale, cfs, samplerate, llimstr, ulimstr) Source source ; int frameheight ; double scale ; double *cfs, samplerate ; /* array of channel centre frequencies */ char *llimstr, *ulimstr; /* lower and upper limit strings */ { int i, framewidth, max_ulim=0, min_llim=99999999; DeclareNew( struct _summary_state *, state ) ; /* Allocate new arrays for integration limits */ state->llim = NewArray(int, frameheight, "image.c for lower limit array" ); state->ulim = NewArray(int, frameheight, "image.c for upper limit array" ); for (i=0 ; i<frameheight ; i++) { /* Convert strings to sample points, allowing for cycles units */ state->llim[i] = Cycles( llimstr, cfs[i], Samplerate() ); state->ulim[i] = Cycles( ulimstr, cfs[i], Samplerate() ); /* Check that llim < ulim, and quit if not so */ if (state->llim[i] >= state->ulim[i]) stitch_error("Warning: gensas integration limits badly ordered\n"); /* Find limits on required framewidth */ if (state->ulim[i] > max_ulim) max_ulim = state->ulim[i]; if (state->llim[i] < min_llim) min_llim = state->llim[i]; } if (min_llim > 0) min_llim = 0; framewidth = max_ulim - min_llim; /* Adjust the integration limits for an sai with triggering on the */ /* right. (The parameters llim <= ulim, but the resulting sai indices */ /* state->llim[i] > state->ulim[i]). */ /* This is done by subtracting from (framewidth-1), where framewidth is */ /* the maximum required number of points, (ie ulim in points). */ for (i=0 ; i<frameheight ; i++) { state->llim[i] = framewidth - state->llim[i] ; state->ulim[i] = framewidth - state->ulim[i] ; } framewidth++; /* extra point as array starts from zeroth location */ /* a blockings requests up into requests of the size specifed */ state->source = NewBlockingSource( source, sizeof ( DataType ) * framewidth * frameheight ) ; state->framewidth = framewidth ; state->frameheight = frameheight ; state->scale = scale ; return ( SetFillableSource( state, summary_callback, "image.c summarising sai" ) ) ; } /*********************** Routines supporting Summary ***********************/ static Pointer summary_callback( state, bytes, buffer ) struct _summary_state *state ; ByteCount *bytes ; DataType *buffer ; { register int last = *bytes == 0 ; register int i, j, ulim, llim, point, points=ToPoints(DataType,*bytes) ; register DataType *sairow; #ifdef FLOAT register DataType sum ; #else register long sum ; #endif /* Pull an sai frame */ for( point=0 ; point < points ; ) /* For each channel (row) in the sai, sum the row between the given limits */ for(i=0 ; i<state->frameheight ; i++) { sairow = PullItems(state->source,state->framewidth,DataType); ulim = state->ulim[i]; llim = state->llim[i]; sum = 0; for (j=ulim ; j<=llim ; j++) sum += sairow[j]; /* store the row-sum, scaled and normalized for the range of the sum */ buffer[point++] = sum*state->scale / (llim-ulim); } if( !last ) return ( (Pointer) buffer ) ; else { Delete( state->llim ) ; Delete( state->ulim ) ; return ( DeleteFillableSource( state ) ) ; } } /********************************************************************************/ /* */ /* initlist initialises the linked list pointer before the other list functions */ /* call the list. */ /* */ /********************************************************************************/ void initlist() { static int status=0; if (status==0) start = endl = pp = qq = (node *)malloc(sizeof(node)); status=1; } /********************************************************************************/ /* */ /* Function getnode allocates storage for a link list node and returns a */ /* pointer to that node */ /* */ /********************************************************************************/ node *getnode() { node *q; q=(node *)malloc(sizeof(node)); if (q==NULL) exit(66); return q; } /********************************************************************************/ /* */ /* The function insertl inserts the value x at the end of the linked list and */ /* moves the end point by another node */ /* */ /********************************************************************************/ node *insertl(x) int x; { node *q=endl; endl=getnode(); q->next=endl; q->index=x; q->val=0; /* fprintf(testfile, "insertl: index is %d, val is %d\n", q->index, q->val); */ return q; } /********************************************************************************/ /* */ /* The function inss searches the list from the start, stops 3 nodes before the */ /* value x, assigns the value -2000 to the three nodes till the (x+1)th node. */ /* */ /********************************************************************************/ void inss(x) int x; { node *q=start; int y=x-3; int z=x; while (q!=endl) { if (q->index == y) do { q->val=-2000; q=q->next; /* fprintf(stderr, "Inss(x) index=%d\n",q->index); */ } while (q->index < z); q=q->next; } } /********************************************************************************/ /* */ /* output_simple_strobe_info routine creates a file pointed by trigger file */ /* which writes the points in the Nap pulse which initiates strobing */ /* */ /********************************************************************************/ void output_simple_strobe_info(state, value) struct _sai_state *state; short value; { if (!(!strcmp(state->switch_info, "off") || !strcmp(state->switch_info, "on"))) fwrite(&(value), sizeof(short), 1 , trigger_file); } /********************************************************************************/ /* */ /* doList_strobe_info creates a file (pointed by "trigger_file") which writes */ /* in binary shorts the points in the Nap pulse which initiates a strobe */ /* It uses a linked list data structure to keep track of the actual time which */ /* causes a strobe, rather than the time when the strobing event occurs. */ /* At stcrit 4 or 5, strobing occurs nwid ms after a nap pulse peak has been */ /* located. */ /* */ /********************************************************************************/ void doList_strobe_info(state, chan) struct _sai_state *state; int chan; { int storev; if (state->Stcrit==2) storev=(state->time)-2; else if (state->Stcrit==5) storev=state->def_strobe_candidate_time[chan]; /* fprintf(stderr, "Storev value of doList is = %d\n", storev);} */ else storev=state->trigtime[chan]; inss(storev); while (pp->index<=storev) /* (pp!=qq) */ { if (!(!strcmp(state->switch_info, "off") || !strcmp(state->switch_info, "on"))) fwrite(&(pp->val), sizeof(short), 1, trigger_file); /* fprintf(testfile, "doList: index is %d, val is %d\n", pp->index, pp->val); */ pp=pp->next; } /* pp=qq; */ } /********************************************************************************/ /* */ /* output_thresh_info writes threshold values (in binary shorts or floats) to */ /* the file opened by the file pointer "thresho_file" */ /* */ /********************************************************************************/ void output_thresh_info(state, chan) struct _sai_state *state; int chan; { S[chan]=state->thresh[chan]; if (!(!strcmp(state->switch_info, "off") || !strcmp(state->switch_info, "on"))) #ifdef FLOAT fwrite((S+chan), sizeof(float), 1, thresho_file); #else fwrite((S+chan), sizeof(short), 1, thresho_file); #endif } /********************************************************************************/ /* */ /* decay_strobe_threshold decays the strobe threshold at every clock */ /* tick either linearly or exponentially */ /* */ /********************************************************************************/ void decay_strobe_threshold(state, chan) struct _sai_state *state; int chan; { if (!strcmp(exposwitch, "on")) state->thresh[chan] *= state->triggerdecay; else { state->thresh[chan] -= state->tlim[chan]; if (state->thresh[chan] < 0) state->thresh[chan]=0; } } /********************************************************************************/ /* */ /* This routine prints debugging information of the local max search */ /* with timeout (i.e. stcrit_ai=5) */ /* */ /********************************************************************************/ void print_trigger_debugging_info(state, chan, strobe_lag, stage) struct _sai_state *state; int chan; int strobe_lag; int stage; { if (stage==1) { if (strobe_lag > state->def_strobe_candidate_time[chan]) val[chan]=strobe_lag-state->def_strobe_candidate_time[chan]; if (strobe_lag == state->def_strobe_candidate_time[chan]) val[chan]=0; if (!strcmp(state->switch_info,"off")); else fprintf(infotxt_file, "Chan %d: :: Actual time of trigger: s_lag=%d, trigtime=%d\n", chan, strobe_lag, state->def_strobe_candidate_time[chan]); } /* fprintf(infotxt_file, "Chan %d: Actual time of trigger: s_lag=%d, trigtime=%d, shift=%d\n", chan, strobe_lag, state->def_strobe_candidate_time[chan], val[chan]*state->chans); */ if (stage==2) { if (!strcmp(state->switch_info,"off")); else fprintf(infotxt_file, "Chan %d: Loc Max Start Time = %d, trigcurrent %d, Threshold=%f, shift=%d\n", chan, state->time, state->trigheight[chan], state->thresh[chan], (ttime[chan]-initial_strobe_candidate_time[chan])); } if (stage==3) { if (!strcmp(state->switch_info,"off")); else fprintf(infotxt_file, "Chan %d: Loc Max is %d Time is %d\n", chan, state->trigheight[chan], state->trigtime[chan] ); } if (stage==4) { if (!strcmp(state->switch_info,"off")); else fprintf(infotxt_file, "Chan %d: At the End of Search (At peak), trigtime = %d time = %d height = %d\n", chan, state->trigtime[chan], state->time, state->trigheight[chan]); } if (stage==5) { if (!strcmp(state->switch_info,"off")); else fprintf(infotxt_file, "Chan %d: At the End of Search, trigtime = %d time = %d height = %d\n", chan, state->trigtime[chan], state->time, state->trigheight[chan]); } }