Mercurial > hg > aim92

comparison saitools/saisummary.c @ 0:5242703e91d3 tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Initial checkin for AIM92 aimR8.2 (last updated May 1997).
author tomwalters
date Fri, 20 May 2011 15:19:45 +0100
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:5242703e91d3
1 /*
2 Copyright (c) Applied Psychology Unit, Medical Research Council. 1993
3 ===========================================================================
4
5 Permission to use, copy, modify, and distribute this software without fee
6 is hereby granted for research purposes, provided that this copyright
7 notice appears in all copies and in all supporting documentation, and that
8 the software is not redistributed for any fee (except for a nominal
9 shipping charge). Anyone wanting to incorporate all or part of this
10 software in a commercial product must obtain a license from the Medical
11 Research Council.
12
13 The MRC makes no representations about the suitability of this
14 software for any purpose. It is provided "as is" without express or
15 implied warranty.
16
17 THE MRC DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE, INCLUDING
18 ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS, IN NO EVENT SHALL
19 THE A.P.U. BE LIABLE FOR ANY SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES
20 OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS,
21 WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION,
22 ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS
23 SOFTWARE.
24 */
25 /*---------------------------------------------------------------------------*/
26 /*
27 * saisummary
28 *
29 * builds a summary .sai
30 * Output is one channel.
31 *
32 * M. Akeroyd. May 1993. Version 2.0
33 * Revised Winter 1994.
34 * revised Autumn 1993.
35 * variuos fiddles added Summer 1994.
36 * (and some new options)
37 * and some more Spring 1995.
38 * and more Summer 1995.
39 */
40
41
42 #include <stdio.h>
43 #include <string.h>
44 #include <stdlib.h>
45
46 #include "tip.h"
47
48
49
50 /* Function declarations */
51 /*---------------------------------------------------------------------------*/
52
53
54 void parsecommandline(int argc, char *argv[], char inputfn[], char outputfn[]);
55
56 int readheader(FILE *inputfp);
57 void checkheader();
58 void writeheader(FILE *outputfp);
59
60 void writedata_output (FILE *outputfp, int samples_to_write);
61 FILE *open_file (char filefn[], FILE *dir_default, int streamtype);
62 void copytooutput_single(int no_peaks, int frame, int channel, int trigger_peak);
63
64 int findpeaksreverse();
65 int findbigpeaks(int total_peaks);
66 void integratebigpeaks(int total_bigpeaks);
67
68 void changeheader(int no_frames, int new_pwidth, int new_nwidth, int framewidth_samples_output);
69 void changeheader_B(int new_mincf, int new_maxcf, int new_channels);
70
71 void write_matrix_hori_time(int framewidth, FILE *outputfigurefp);
72
73
74
75 /* Data arrays: global */
76 /*--------------------------------------------------------------------------*/
77
78
79 short inputdata[MAX_DATA];
80 short outputdata[MAX_DATA];
81 short outputfiguredata[MAX_DATA];
82 short outputgrounddata[MAX_DATA];
83 float summary[MAX_DATA];
84 float framesummary[MAX_DATA];
85
86 short interval[MAX_CHANNELS][MAX_DATA]; /* required for matrix.c */
87 short summation[MAX_CHANNELS]; /* required for matrix.c */
88
89 int clip[MAX_DATA];
90
91 /* other variables */
92 /*-------------------------------------------------------------------------*/
93
94
95 /* variables read from header */
96 /*-------------------------------------------------------------------------*/
97
98 char header[MAX_LINES_HEADER][MAX_LINE_LENGTH];
99 int header_lines;
100
101
102
103 /* .sai parameters */
104 /*-------------------------------------------------------------------------*/
105
106 int no_frames;
107 int frameheight; /* number of channels */
108 int framewidth_samples; /* pwidth + nwidth * samplerate */
109 int frameshift_samples; /* frstep_aid * samplerate */
110 int pwidth; /* in msecs */
111 int nwidth; /* in msecs: NEGATIVE */
112 int width_win; /* pixels */
113 int height_win; /* pixels */
114 long samplerate; /* samples per sec */
115 int mincf; /* Hz */
116 int maxcf; /* Hz */
117
118
119
120
121
122 /* misc */
123 /*---------------- ---------------------------------------------------------*/
124
125 char progname[MAX_STRING_LENGTH];
126 char outputfigurefn[MAX_STRING_LENGTH];
127
128
129 int verboseflag = OFF; /* -v */
130 int summationflag = OFF;
131 int frameflag = OFF; /* -f */
132 int asciiflag = OFF;
133 int greyscaleflag = OFF;
134 int divideflag = OFF;
135
136 int oppositearchflag = OFF; /* -oparch : see saigraph.c for info */
137
138 int nstopflag = ON;
139 int nwarnflag = OFF;
140 int n32760flag = OFF;
141 int n0flag = OFF;
142 int headeroutputflag = ON;
143 int singleframeflag = OFF;
144 int frameaverageflag = ON;
145 int specificchannelflag = OFF;
146 int cutchannels[MAX_CHANNELS];
147 int nchannels=0;
148
149 int framestart;
150 int frameend;
151 int cut_framestart;
152 int cut_frameend;
153
154
155 /* .............. Main .............................*/
156 /* .....................................................................*/
157 /* .....................................................................*/
158
159
160
161
162 void main (int argc, char *argv[])
163 {
164 int sample;
165 int n;
166
167 int frame, channel;
168 int trigger_peak = 0;
169 int no_peak = 0, total_peaks =0;
170 int total_bigpeaks = 0;
171 int trueframeheight = 0;
172 int header_bytes = 0;
173
174 double temp;
175
176 char inputfn[MAX_STRING_LENGTH],
177 outputbasefn[MAX_STRING_LENGTH];
178
179 FILE *inputfp, *outputfigurefp;
180
181 strcpy(progname, argv[0]);
182 strcpy(inputfn, "");
183 strcpy(outputbasefn, "");
184 strcpy(outputfigurefn, "");
185
186 for (n=0; n<MAX_CHANNELS; n++)
187 cutchannels[n] = OFF;
188
189 /* parse command line */
190 parsecommandline(argc, argv, inputfn, outputbasefn);
191
192 /* for (n=0; n<MAX_CHANNELS; n++)
193 if (cutchannels[n] == ON)
194 fprintf(stderr, "%d\n", n);
195 */
196 if (specificchannelflag == OFF)
197 for (n=0; n<MAX_CHANNELS; n++)
198 cutchannels[n] = ON;
199
200 if (strcmp(outputbasefn, "") == 0)
201 strcpy(outputfigurefn, "");
202 else {
203 strcpy(outputfigurefn, outputbasefn);
204 /* strcat(outputfigurefn, OUTPUT_EXT);*/
205 }
206
207 /* open files */
208 /* default directions are:
209 * input = stdin
210 * figure = stdout
211 */
212
213 inputfp = open_file(inputfn, stdin, READ);
214 outputfigurefp = open_file(outputfigurefn, stdout, WRITE);
215
216 /* read Header */
217 header_bytes = readheader(inputfp);
218
219 /* do some error-checking on the header, and set the fread
220 * pointer to the right place */
221 checkheader();
222
223 if (frameflag == ON) {
224 cut_frameend = frameend;
225 cut_framestart = framestart; }
226 else {
227 cut_frameend = no_frames;
228 cut_framestart = 1; }
229
230 if (frameflag == ON)
231 changeheader((cut_frameend - cut_framestart + 1), pwidth, nwidth, framewidth_samples);
232
233 if (singleframeflag == ON)
234 changeheader(1, pwidth, nwidth, framewidth_samples);
235
236 /* reset the variables */
237 /* arguments: mincf maxcf channels */
238 changeheader_B(0, 0, 1);
239
240 /* write header */
241 if (singleframeflag == OFF) {
242 if ((asciiflag == OFF) && (headeroutputflag == ON))
243 writeheader(outputfigurefp);
244 else ;}
245 else {
246 if ((asciiflag == OFF) && (headeroutputflag == ON))
247 writeheader(outputfigurefp);
248 else if (asciiflag == ON){
249 /* begining of ascii header. The rest depends upon the output format. */
250 fprintf(outputfigurefp, "# saisummary output. %s\n", inputfn);
251 fprintf(outputfigurefp, "# \n");
252 fprintf(outputfigurefp, "# samplerate=%i pwidth=%d nwidth=%d frameheight=%i\n", samplerate, pwidth, nwidth, frameheight);
253 fprintf(outputfigurefp, "# \n");
254 }}
255
256
257 /*---------------------------------------------------------------------------*/
258
259 /* Main loop */
260
261
262 if (verboseflag==ON) {
263 fprintf(stderr, " frames ");
264 fflush(stderr);}
265
266 for(sample=0; sample<framewidth_samples; sample++)
267 framesummary[sample] = 0.0;
268
269 for (frame=1; frame<=no_frames; frame++) {
270 if (verboseflag==ON) {
271 fprintf(stderr, " %i ", frame);
272 fflush(stderr);}
273
274 for(sample=0; sample<framewidth_samples; sample++) {
275 interval[1][sample] = outputfiguredata[sample] = 0;
276 summary[sample] = 0.0;
277 clip[sample] = OFF;}
278
279 /*---------------------------------------------------------------------------*/
280
281 nchannels=0;
282 for (channel=1; channel <=frameheight; channel ++) {
283
284 /* clear arrays */
285 for(sample=0; sample<framewidth_samples; sample++)
286 inputdata[sample] = 0;
287
288 /* load data */
289 fread (inputdata, 2, (size_t) framewidth_samples, inputfp);
290
291 /* This next is a simple input check: */
292 /* revised for various warnings */
293 for(sample=0; sample<framewidth_samples; sample++) {
294 /* fprintf(stderr, "data: %i sample: %i channel: %i \n", inputdata[sample], sample, channel);*/
295 if (inputdata[sample] < 0 ) {
296 if (nstopflag == ON) {
297 if (verboseflag == ON) fprintf(stderr, "\n");
298 fprintf(stderr, "%s: negative data: ", progname);
299 fprintf(stderr, "%i sample: %i channel: %i \n", inputdata[sample], sample, channel);
300 fprintf(stderr, "bye.\n");
301 exit(-1);}
302 if (nwarnflag == ON){
303 if (verboseflag == ON) fprintf(stderr, "\n");
304 fprintf(stderr, "warning: negative data: %i (sample: %i channel: %i)\n", inputdata[sample], sample, channel);}
305 if (n0flag == ON)
306 inputdata[sample] = 0;
307 else if (n32760flag == ON)
308 inputdata[sample] = 32760;}}
309
310 /* add data to the summary autoco: outputfigurefp */
311 if ((frame >= cut_framestart) && (frame <= cut_frameend)){
312 if (cutchannels[channel] == ON){
313 nchannels ++;
314 /*fprintf(stderr, "c%d ", channel);*/
315 for(sample=0; sample<framewidth_samples; sample++)
316 summary[sample] += (double) inputdata[sample];}
317 /*fprintf(stderr, "%d hello\n", channel);*/
318 }
319 } /* channel */
320
321 /* fprintf(stderr, "c%d ", nchannels);*/
322 /* divide everything by no. of channels */
323 if (divideflag == ON) {
324 for(sample=0; sample<framewidth_samples; sample++){
325 temp = (double) summary[sample] / nchannels;
326 if ( ( temp > 32760.0)||( temp < 0)) {
327 if (clip[sample] == OFF) {
328 fprintf(stderr, "\nclipping: time interval %i\n", sample);
329 clip[sample] = ON; }
330 temp = 32761.0; }
331 interval[1][sample] = outputfiguredata[sample] = (short) temp; }}
332 else {
333 for(sample=0; sample<framewidth_samples; sample++){
334 temp = (double) summary[sample];
335 if (( temp > 32760.0) || (temp < 0.0)){
336 if (clip[sample] == OFF) {
337 fprintf(stderr, "\nclipping: time interval %i\n", sample);
338 clip[sample] = ON; }
339 temp = 32761.0; }
340 interval[1][sample] = outputfiguredata[sample] = (short) temp; }}
341
342 /* This next is a simple output check: */
343 for(sample=0; sample<framewidth_samples; sample++)
344 if (outputfiguredata[sample] < 0 ) {
345 fprintf(stderr, " %s : something's gone wrong: the output data is negative. \n", progname);
346 exit(-1); }
347
348 /* if required, sum the values */
349 if (summationflag == ON) {
350 /* clear ... */
351 summation[1] = 0;
352 clip[0] = OFF;
353 /* set ... */
354 for (sample=0; sample < framewidth_samples; sample++) {
355 summation[1] += interval[1][sample];
356 /* clip check */
357 if (( summation[1] > 32760 )||(summation[1] < 0 )) {
358 if (clip[0] == OFF) {
359 fprintf(stderr, "\nclipping: summation\n");
360 clip[0] = ON; }
361 summation[1] = 32761;}}
362
363 }
364
365 /* add to frame average */
366 if ((frame >= cut_framestart) && (frame <= cut_frameend)){
367 for (sample=0; sample < framewidth_samples; sample++)
368 framesummary[sample] += outputfiguredata[sample];
369 }
370
371 /* write output: outputfiguredata[] for the .sai
372 * interval[] & summation[] for the ascii
373 * Note hack for frameheight
374 * BUT: not if only 1 frame is to be output
375 */
376
377 if ((frame >= cut_framestart) && (frame <= cut_frameend)){
378 if (singleframeflag == OFF) {
379 if (asciiflag == OFF)
380 writedata_output(outputfigurefp, framewidth_samples);
381 else {
382 trueframeheight=frameheight;
383 frameheight=1;
384 write_matrix_hori_time(framewidth_samples, outputfigurefp);
385 frameheight=trueframeheight; }}
386 else ;
387 }
388
389 /* -------------------------------------------------------------------------*/
390
391 } /* frame */
392
393 if (singleframeflag == ON) {
394 if (verboseflag == ON) {
395 fprintf(stderr, " averaging ");
396 fflush(stderr);}
397 clip[0] = OFF;
398 for (sample=0; sample < framewidth_samples; sample++) {
399 if (frameaverageflag == ON)
400 temp = framesummary[sample]/(cut_frameend - cut_framestart + 1);
401 else
402 temp = framesummary[sample];
403 /* clip check */
404 if (temp > 32760) {
405 if (clip[0] == OFF) {
406 fprintf(stderr, "\nclipping: frameaverage at sample %i\n", sample);
407 clip[0] = ON; }
408 temp = 32760;}
409 outputfiguredata[sample]= (short)temp;
410 interval[1][sample] = outputfiguredata[sample];}
411 if (asciiflag == OFF)
412 writedata_output(outputfigurefp, framewidth_samples);
413 else {
414 trueframeheight=frameheight;
415 frameheight=1;
416 write_matrix_hori_time(framewidth_samples, outputfigurefp);
417 frameheight=trueframeheight; }}
418 else ;
419
420
421 /* Tidy up and exit */
422 fprintf(stderr, "\n");
423 fclose(inputfp);
424 fclose(outputfigurefp);
425
426 if ( ferror(inputfp) != 0) {
427 fprintf(stderr, " %s : error closing input file.\n", progname);
428 exit(-1);}
429
430 if ( ferror(outputfigurefp) != 0) {
431 fprintf(stderr, " %s : error closing figure file.\n", progname);
432 exit(-1);}
433
434 exit(0);
435
436 } /* Main */
437
438
439
440
441
442
443 /* .....................................................................*/
444 /* .....................................................................*/
445 /* .....................................................................*/
446 /* .....................................................................*/
447
448
449
450
451
452
453 void parsecommandline(int argc, char *argv[], char inputfn[], char outputfn[])
454 {
455 int x=1, helpflag = OFF;
456 int group;
457
458 while (x < argc){
459 if (!strcmp(argv[x], "-o")) {strcpy(outputfn, argv[x+1]); x+=2;}
460 else if (!strcmp(argv[x], "-output")) {strcpy(outputfn, argv[x+1]); x+=2;}
461 else if (!strcmp(argv[x], "-i")) {strcpy(inputfn, argv[x+1]); x+=2;}
462 else if (!strcmp(argv[x], "-input")) {strcpy(inputfn, argv[x+1]); x+=2;}
463 else if (!strcmp(argv[x], "-help")) {helpflag = ON; x+=1;}
464 else if (!strcmp(argv[x], "-h")) {helpflag = ON; x+=1;}
465 else if (!strcmp(argv[x], "-verbose")) {verboseflag = ON; x+=1;}
466 else if (!strcmp(argv[x], "-v")) {verboseflag = ON; x+=1;}
467 else if (!strcmp(argv[x], "-s")) {summationflag = ON; x+=1;}
468 else if (!strcmp(argv[x], "-g")) {greyscaleflag = ON; x+=1;}
469 else if (!strcmp(argv[x], "-gray")) {greyscaleflag = ON; x+=1;}
470 else if (!strcmp(argv[x], "-grey")) {greyscaleflag = ON; x+=1;}
471 else if (!strcmp(argv[x], "-ar")) {asciiflag = ON; x+=1;}
472 else if (!strcmp(argv[x], "-d")) {divideflag = ON; x+=1;}
473 else if (!strcmp(argv[x], "-channels")) {specificchannelflag = ON;
474 x+=1;
475 while (x<argc){
476 if ((int) strspn(argv[x], "-") == 0)
477 cutchannels[atoi(argv[x])]=ON;
478 else break;
479 x++;}}
480 else if (!strcmp(argv[x], "-oneframe")) {singleframeflag = ON; frameaverageflag = ON; x+=1;}
481 else if (!strcmp(argv[x], "-sumframe")) {singleframeflag = ON; frameaverageflag = OFF; x+=1;}
482 else if (!strcmp(argv[x], "-oparch")) {oppositearchflag = ON; x+=1;}
483 else if (!strcmp(argv[x], "-nstop")) {nstopflag = ON; x+=1;}
484 else if (!strcmp(argv[x], "-header")) {headeroutputflag = ON; x+=1;}
485 else if (!strcmp(argv[x], "-noheader")) {headeroutputflag = OFF; x+=1;}
486 else if (!strcmp(argv[x], "-nwarn")) {nwarnflag = ON; nstopflag = OFF; x+=1;}
487 else if (!strcmp(argv[x], "-nignore")) {nwarnflag = OFF; nstopflag = OFF; x+=1;}
488 else if (!strcmp(argv[x], "-n0")) {n0flag = ON; n32760flag = OFF; nstopflag = OFF; x+=1;}
489 else if (!strcmp(argv[x], "-n32760")) {n0flag = OFF; n32760flag = ON; nstopflag = OFF; x+=1;}
490 else if (!strcmp(argv[x], "-frames")) {frameflag=ON; framestart=atoi(argv[x+1]); frameend=atoi(argv[x+2]); x+=3;}
491
492 else {fprintf(stderr, "%s: unknown option %s\n", progname, argv[x]);
493 exit(-1);}
494 }
495
496 if (helpflag == ON)
497 {
498 fprintf(stderr, "\n");
499 fprintf(stderr, " saisummary\n");
500 fprintf(stderr, " -----------------------------------------------------------\n");
501 fprintf(stderr, " Works out a summary image (cf a summary autocorrelogram).\n");
502 fprintf(stderr, " The -oparch option is REQUIRED if reading a Sun .sai files \n");
503 fprintf(stderr, " on a DEC cpu or vice-versa.\n");
504 fprintf(stderr, " Add .sai to both input and output filenames.\n");
505 fprintf(stderr, "\n");
506 fprintf(stderr, " -i <.sai file> input file default = stdin\n");
507 fprintf(stderr, " -o <.sai file> output file default = stdout\n");
508 fprintf(stderr, " -oparch assume input was generated on an opposing architecture cpu\n");
509 fprintf(stderr, "\n");
510 fprintf(stderr, " -frames <start> <end> (inclusive) (default=all)\n");
511 fprintf(stderr, " -channels <> <> ..... specify channel numbers to be included (1++)\n");
512 fprintf(stderr, " -d divide by no. of channels (Warning: precision may be lost).\n");
513 fprintf(stderr, " -ar ASCII output: rows=time intervals (requires -oneframe)\n");
514 fprintf(stderr, " -s include summations (ascii only)\n");
515 fprintf(stderr, "\n");
516 fprintf(stderr, " -header include AIM header in output file (default)\n");
517 fprintf(stderr, " -noheader don't ... ... ... \n");
518 fprintf(stderr, " -v verbose output (to stderr)\n");
519 fprintf(stderr, " -g change .sai format to 'grayscale'\n");
520 fprintf(stderr, " -oneframe average frames together; one (1) output frame\n");
521 fprintf(stderr, " -sumframe sum frames together; one (1) output frame\n");
522 fprintf(stderr, "\n");
523 fprintf(stderr, " -nstop if negative inputs found, stop processing (default)\n");
524 fprintf(stderr, " -nwarn if negative inputs found, continue (but warn user) \n");
525 fprintf(stderr, " -nignore if negative inputs found, continue \n");
526 fprintf(stderr, " -n0 if negative inputs found, replace with 0 \n");
527 fprintf(stderr, " -n32760 if negative inputs found, replace with 32760 \n");
528 fprintf(stderr, "\n\n");
529 exit(-1);}
530 }
531
532 /* The End. */
533 /*------------------------------------------------------------------------------------------------------*/