tomwalters@0: /* tomwalters@0: Copyright (c) Applied Psychology Unit, Medical Research Council. 1993 tomwalters@0: =========================================================================== tomwalters@0: tomwalters@0: Permission to use, copy, modify, and distribute this software without fee tomwalters@0: is hereby granted for research purposes, provided that this copyright tomwalters@0: notice appears in all copies and in all supporting documentation, and that tomwalters@0: the software is not redistributed for any fee (except for a nominal tomwalters@0: shipping charge). Anyone wanting to incorporate all or part of this tomwalters@0: software in a commercial product must obtain a license from the Medical tomwalters@0: Research Council. tomwalters@0: tomwalters@0: The MRC makes no representations about the suitability of this tomwalters@0: software for any purpose. It is provided "as is" without express or tomwalters@0: implied warranty. tomwalters@0: tomwalters@0: THE MRC DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE, INCLUDING tomwalters@0: ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS, IN NO EVENT SHALL tomwalters@0: THE A.P.U. BE LIABLE FOR ANY SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES tomwalters@0: OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, tomwalters@0: WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, tomwalters@0: ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS tomwalters@0: SOFTWARE. tomwalters@0: */ tomwalters@0: tomwalters@0: /*------------------------------------------------------------------------*/ tomwalters@0: tomwalters@0: /* saigraph.c tomwalters@0: * -------------- tomwalters@0: * tomwalters@0: * Works out the (first-order) temporal intervals between peaks in a .sai file, tomwalters@0: * and saves them in a SINGLE-FRAME .sai tomwalters@0: * In essence, a histogram. tomwalters@0: * tomwalters@0: * Uses ALL peaks: if you just want the locallly maxima peaks, then use tomwalters@0: * trbigpeak first. tomwalters@0: * tomwalters@0: * Output is another .sai file: but only a single frame. The header is tomwalters@0: * modified accordingly. tomwalters@0: * tomwalters@0: * tomwalters@0: * tomwalters@0: * M. Akeroyd. May 1993. Version 2.0 tomwalters@0: * tomwalters@0: * Autumn 1993: added grayscaleflag tomwalters@0: * Winter 1994: allowed weighting, and other general revisions. tomwalters@0: */ tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: #include tomwalters@0: #include tomwalters@0: #include tomwalters@0: tomwalters@0: #include "tip.h" tomwalters@0: tomwalters@0: tomwalters@0: /* Function declarations */ tomwalters@0: /*------------------------------------------------------------*/ tomwalters@0: tomwalters@0: void parsecommandline (int argc, char *argv[], char inputfn[], char outputfn[]); tomwalters@0: tomwalters@0: int readheader(FILE *inputfp); tomwalters@0: void checkheader(); tomwalters@0: void writeheader(FILE *outputfp); tomwalters@0: tomwalters@0: void changeheader(int no_frames, int new_pwidth, int new_nwidth, int framewidth_samples_output); tomwalters@0: tomwalters@0: int findpeaksreverse (); tomwalters@0: void find_intervals(int total_peaks, int channel, int max_interval); tomwalters@0: tomwalters@0: void writedata_output (FILE *outputfp, int samples_to_write); tomwalters@0: void writedata_output_fg (FILE *figurefp, FILE *groundfp, int samples_to_write); tomwalters@0: FILE *open_file (char filefn[], FILE *dir_default, int streamtype); tomwalters@0: tomwalters@0: void copytooutput_fg_single (int total_peaks, int frame, int channel, int trigger_peak); tomwalters@0: tomwalters@0: void write_matrix_hori_time(int framewidth, FILE *outputfigurefp); tomwalters@0: void write_matrix_hori_freq(int framewidth, FILE *outputfigurefp); tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: /* Data arrays: global */ tomwalters@0: /*-----------------------------------------------------------------*/ tomwalters@0: tomwalters@0: tomwalters@0: short inputdata[MAX_DATA]; /* input .sai: a single channel of a frame */ tomwalters@0: short outputdata[MAX_DATA]; /* output .sai: a single channel */ tomwalters@0: short outputfiguredata[MAX_DATA]; /* output figure.sai: a single channel */ tomwalters@0: short outputgrounddata[MAX_DATA]; /* output ground.sai: a single channel */ tomwalters@0: short clip[MAX_CHANNELS]; tomwalters@0: tomwalters@0: short interval[MAX_CHANNELS][MAX_DATA]; /* histogram */ tomwalters@0: short summation[MAX_CHANNELS]; tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: /* other variables */ tomwalters@0: /*----------------------------------------------------------------*/ tomwalters@0: tomwalters@0: struct Peak peak[MAX_PEAKS]; /* used in findpeaksreverse(): tomwalters@0: * single channel */ tomwalters@0: struct NextEvent nextevent[MAX_NEXTEVENTS]; /* used in findpeaksreverse(): tomwalters@0: * single channel */ tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: /* variables read from header */ tomwalters@0: /*---------------------------------------------------------------*/ tomwalters@0: tomwalters@0: char header[MAX_LINES_HEADER][MAX_LINE_LENGTH]; tomwalters@0: int header_lines; tomwalters@0: tomwalters@0: int no_frames; tomwalters@0: int frameheight; /* number of channels */ tomwalters@0: int framewidth_samples; /* pwidth + nwidth * samplerate */ tomwalters@0: int frameshift_samples; /* frstep_aid * samplerate */ tomwalters@0: int pwidth; /* in msecs */ tomwalters@0: int nwidth; /* in msecs: NEGATIVE */ tomwalters@0: int width_win; /* pixels */ tomwalters@0: int height_win; /* pixels */ tomwalters@0: long samplerate; /* samples per sec */ tomwalters@0: int mincf; /* Hz */ tomwalters@0: int maxcf; /* Hz */ tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: /* misc */ tomwalters@0: /*-----------------------------------------------------------*/ tomwalters@0: tomwalters@0: char progname[MAX_STRING_LENGTH]; tomwalters@0: char outputfigurefn[MAX_STRING_LENGTH]; tomwalters@0: tomwalters@0: tomwalters@0: int verboseflag = OFF; /* -v */ tomwalters@0: int removeflag = OFF; /* -r */ tomwalters@0: double removevalue = 0.1; tomwalters@0: int widthflag = OFF; /* -width : if OFF, copy pwdith/nwidth tomwalters@0: * off Input .sai */ tomwalters@0: int total_width; /* -width : total width of new .sai */ tomwalters@0: int asciiflag = OFF; /* -a */ tomwalters@0: int asciireverseflag = OFF; /* -ar */ tomwalters@0: int summationflag = OFF; /* -s */ tomwalters@0: int greyscaleflag = OFF; /* -g */ tomwalters@0: int weightflag = OFF; /* -weight */ tomwalters@0: int oppositearchflag = OFF; /* -oparch: if ON, then the header-reading tomwalters@0: * code will load one (1) extra byte. This tomwalters@0: * is so that Sun .sai/.nap can be read tomwalters@0: * on a DEC, and equally the reverse. tomwalters@0: * It is NOT required if reading Sun on tomwalters@0: * Sun, or DEC on DEC. */ tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: /*............... Main ................................*/ tomwalters@0: /*.........................................................................*/ tomwalters@0: /*.........................................................................*/ tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: void main (int argc, char *argv[]) tomwalters@0: { tomwalters@0: int sample; tomwalters@0: int n; tomwalters@0: tomwalters@0: int frame, channel; tomwalters@0: int trigger_peak = 0; tomwalters@0: int no_peak = 0, total_peaks =0, total_peaks_2; tomwalters@0: int framewidth_samples_input, tomwalters@0: framewidth_samples_output; tomwalters@0: int new_pwidth, tomwalters@0: new_nwidth; tomwalters@0: int header_bytes = 0; tomwalters@0: tomwalters@0: char inputfn[MAX_STRING_LENGTH], tomwalters@0: outputbasefn[MAX_STRING_LENGTH]; tomwalters@0: tomwalters@0: FILE *inputfp, *outputfigurefp; tomwalters@0: tomwalters@0: int clipflag = OFF; tomwalters@0: tomwalters@0: /*-----------------------------------------*/ tomwalters@0: tomwalters@0: strcpy(progname, argv[0]); tomwalters@0: strcpy(inputfn, ""); tomwalters@0: strcpy(outputbasefn, ""); tomwalters@0: strcpy(outputfigurefn, ""); tomwalters@0: tomwalters@0: /* parse command line */ tomwalters@0: parsecommandline(argc, argv, inputfn, outputbasefn); tomwalters@0: tomwalters@0: if (strcmp(outputbasefn, "") == 0) tomwalters@0: strcpy(outputfigurefn, ""); tomwalters@0: else tomwalters@0: strcpy(outputfigurefn, outputbasefn); tomwalters@0: tomwalters@0: tomwalters@0: /* open files */ tomwalters@0: /* default directions are: tomwalters@0: * input = stdin tomwalters@0: * figure = stdout tomwalters@0: */ tomwalters@0: tomwalters@0: inputfp = open_file(inputfn, stdin, READ); tomwalters@0: outputfigurefp = open_file(outputfigurefn, stdout, WRITE); tomwalters@0: tomwalters@0: /* read Header */ tomwalters@0: header_bytes = readheader(inputfp); tomwalters@0: checkheader(); tomwalters@0: tomwalters@0: tomwalters@0: /* reset the variables */ tomwalters@0: framewidth_samples_input = framewidth_samples; tomwalters@0: if (widthflag == OFF) { tomwalters@0: new_pwidth = + pwidth; tomwalters@0: new_nwidth = -nwidth; tomwalters@0: framewidth_samples_output = framewidth_samples_input;} tomwalters@0: else { tomwalters@0: new_pwidth = + NEW_PWIDTH; tomwalters@0: new_nwidth = (new_pwidth + total_width); tomwalters@0: framewidth_samples_output = (int) ((int) (abs(new_pwidth) + abs (new_nwidth)) * samplerate ) / 1000;} tomwalters@0: tomwalters@0: if (framewidth_samples_output >= MAX_DATA) { tomwalters@0: fprintf(stderr, "%s: new frame is too wide at %i: only allowed %i samples\n", progname, framewidth_samples_output, MAX_DATA); tomwalters@0: exit(-1); } tomwalters@0: tomwalters@0: tomwalters@0: /* change header */ tomwalters@0: changeheader(1, new_pwidth, new_nwidth, framewidth_samples_output); tomwalters@0: tomwalters@0: /* write header */ tomwalters@0: if (asciiflag == OFF) tomwalters@0: writeheader(outputfigurefp); tomwalters@0: else { tomwalters@0: /* begining of ascii header. The rest depends on the matrix format ... */ tomwalters@0: fprintf(outputfigurefp, "# saigraph output. %s\n", inputfn); tomwalters@0: fprintf(outputfigurefp, "# \n"); tomwalters@0: fprintf(outputfigurefp, "# samplerate=%i pwidth=%d nwidth=%d frames=%i ", samplerate, pwidth, nwidth, no_frames); tomwalters@0: fprintf(outputfigurefp, " mincf=%i maxcf=%i channels=%i \n", mincf, maxcf, frameheight); tomwalters@0: fprintf(outputfigurefp, "# \n"); tomwalters@0: } tomwalters@0: tomwalters@0: /*------------------------------------------*/ tomwalters@0: tomwalters@0: /* Main loop */ tomwalters@0: tomwalters@0: if (verboseflag == ON) { tomwalters@0: fprintf(stderr, " frames "); tomwalters@0: fflush(stderr);} tomwalters@0: tomwalters@0: tomwalters@0: /* clear arrays */ tomwalters@0: tomwalters@0: for (channel=0; channel 32765) || tomwalters@0: (summation[channel] < 0)){ tomwalters@0: if (clipflag == OFF) { tomwalters@0: fprintf(stderr, " \nclipping: channel %i, summations.\n ", channel); tomwalters@0: clipflag = ON;} tomwalters@0: summation[channel] = 32765; tomwalters@0: } tomwalters@0: } tomwalters@0: } tomwalters@0: tomwalters@0: /*-------------------------------------------*/ tomwalters@0: tomwalters@0: /* write output */ tomwalters@0: tomwalters@0: if (asciiflag == OFF) { tomwalters@0: for (channel = 1; channel <= frameheight; channel ++){ tomwalters@0: for (sample = 0; sample < framewidth_samples_output; sample++) tomwalters@0: outputfiguredata[sample] = interval[channel][sample]; tomwalters@0: writedata_output(outputfigurefp, framewidth_samples_output); tomwalters@0: } tomwalters@0: } tomwalters@0: else { tomwalters@0: if (asciireverseflag == OFF) tomwalters@0: write_matrix_hori_freq(framewidth_samples_output, outputfigurefp); tomwalters@0: else tomwalters@0: write_matrix_hori_time(framewidth_samples_output, outputfigurefp); tomwalters@0: } tomwalters@0: tomwalters@0: /*-------------------------------------------*/ tomwalters@0: tomwalters@0: /* Tidy up and exit */ tomwalters@0: tomwalters@0: if (verboseflag == ON) tomwalters@0: fprintf(stderr, "\n"); tomwalters@0: fclose(inputfp); tomwalters@0: fclose(outputfigurefp); tomwalters@0: tomwalters@0: if ( ferror(inputfp) != 0) { tomwalters@0: fprintf(stderr, " %s : error closing input file.\n", progname); tomwalters@0: exit(-1);} tomwalters@0: tomwalters@0: if ( ferror(outputfigurefp) != 0) { tomwalters@0: fprintf(stderr, " %s : error closing figure file.\n", progname); tomwalters@0: exit(-1);} tomwalters@0: tomwalters@0: exit(0); tomwalters@0: tomwalters@0: } /* Main */ tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: /*---------------------------------------------------------------------------*/ tomwalters@0: /*---------------------------------------------------------------------------*/ tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: void parsecommandline(int argc, char *argv[], char inputfn[], char outputbasefn[]) tomwalters@0: { tomwalters@0: int x=1, helpflag = OFF; tomwalters@0: tomwalters@0: while (x < argc){ tomwalters@0: if (!strcmp(argv[x], "-o")) {strcpy(outputbasefn, argv[x+1]); x+=2;} tomwalters@0: else if (!strcmp(argv[x], "-output")) {strcpy(outputbasefn, argv[x+1]); x+=2;} tomwalters@0: else if (!strcmp(argv[x], "-i")) {strcpy(inputfn, argv[x+1]); x+=2;} tomwalters@0: else if (!strcmp(argv[x], "-input")) {strcpy(inputfn, argv[x+1]); x+=2;} tomwalters@0: else if (!strcmp(argv[x], "-help")) {helpflag = ON; x+=1;} tomwalters@0: else if (!strcmp(argv[x], "-h")) {helpflag = ON; x+=1;} tomwalters@0: else if (!strcmp(argv[x], "-verbose")) {verboseflag = ON; x+=1;} tomwalters@0: else if (!strcmp(argv[x], "-v")) {verboseflag = ON; x+=1;} tomwalters@0: else if (!strcmp(argv[x], "-weight")) {weightflag = NORMAL_WEIGHTING; x+=1;} tomwalters@0: else if (!strcmp(argv[x], "-weightlog")) {weightflag = LOG_WEIGHTING; x+=1;} tomwalters@0: else if (!strcmp(argv[x], "-width")) {widthflag = ON; total_width=atoi(argv[x+1]); x+=2;} tomwalters@0: else if (!strcmp(argv[x], "-r")) { removeflag = ON; removevalue = atof(argv[x+1]); x+=2;} tomwalters@0: else if (!strcmp(argv[x], "-a")) { asciiflag = ON;x+=1;} tomwalters@0: else if (!strcmp(argv[x], "-ar")) { asciiflag = ON; asciireverseflag=ON; x+=1;} tomwalters@0: else if (!strcmp(argv[x], "-s")) { summationflag=ON; x+=1;} tomwalters@0: else if (!strcmp(argv[x], "-g")) { greyscaleflag=ON; x+=1;} tomwalters@0: else if (!strcmp(argv[x], "-grey")) { greyscaleflag=ON; x+=1;} tomwalters@0: else if (!strcmp(argv[x], "-gray")) { greyscaleflag=ON; x+=1;} tomwalters@0: else if (!strcmp(argv[x], "-oparch")) { oppositearchflag=ON; x+=1;} tomwalters@0: else {fprintf(stderr, "%s: unknown option %s\n", progname, argv[x]); tomwalters@0: exit(-1);} tomwalters@0: } tomwalters@0: tomwalters@0: if (helpflag == ON) tomwalters@0: { tomwalters@0: fprintf(stderr, "\n"); tomwalters@0: fprintf(stderr, " %s\n", progname); tomwalters@0: fprintf(stderr, " -----------------------------------------------------------------------------\n"); tomwalters@0: fprintf(stderr, " Works out a histogram of peak intervals, from a .sai\n"); tomwalters@0: fprintf(stderr, " Default output is a .sai file: if -a or -ar specified,\n"); tomwalters@0: tomwalters@0: fprintf(stderr, " then output is a ASCII matrix.\n"); tomwalters@0: fprintf(stderr, "\n"); tomwalters@0: fprintf(stderr, " The -oparch option is REQUIRED if reading a Sun .sai files \n"); tomwalters@0: fprintf(stderr, " on a DEC cpu or vice-versa.\n"); tomwalters@0: fprintf(stderr, " Add .sai to both input and output filenames.\n"); tomwalters@0: tomwalters@0: fprintf(stderr, " -----------------------------------------------------------------------------\n"); tomwalters@0: fprintf(stderr, "\n"); tomwalters@0: fprintf(stderr, " -i <.sai file> input file default = stdin\n"); tomwalters@0: fprintf(stderr, " -o <.sai file> output file default = stdout\n"); tomwalters@0: fprintf(stderr, " -oparch assume input was generated on an opposing architecture\n"); tomwalters@0: fprintf(stderr, "\n"); tomwalters@0: fprintf(stderr, " -a ASCII output: rows = channels\n"); tomwalters@0: fprintf(stderr, " -ar ASCII output: rows = time intervals\n"); tomwalters@0: fprintf(stderr, " -s include summations over time intervals (ASCII options only)\n"); tomwalters@0: fprintf(stderr, "\n"); tomwalters@0: fprintf(stderr, " -width width of graph.sai (msecs).\n"); tomwalters@0: fprintf(stderr, " -r remove any peaks F(n) < F(n+1) & F(n-1) * fract \n"); tomwalters@0: fprintf(stderr, " default = %.3f\n", removevalue); tomwalters@0: fprintf(stderr, " -weight weight histogram entries by mean height of peaks\n"); tomwalters@0: fprintf(stderr, " -weightlog as per -weight, but use 10log(mean height)\n"); tomwalters@0: fprintf(stderr, "\n"); tomwalters@0: fprintf(stderr, " -v verbose output (to stderr)\n"); tomwalters@0: fprintf(stderr, " -g change .sai 'view' format to 'grayscale'.\n"); tomwalters@0: fprintf(stderr, "\n\n"); tomwalters@0: exit(-1);} tomwalters@0: tomwalters@0: } tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: tomwalters@0: /* The End */ tomwalters@0: /*---------------------------------------------------------------------------*/