/* Foetal Monitoring System   AJF   October 1998
   Example code, non-functional, for educational purposes only, may contain nuts, etc. */

#include <stdio.h>

#define NFGAIN	1.003141582e+00f    /* notch filter gain, reported by mkfilter */
#define BPFGAIN 1.008760041e+02f    /* bandpass filter gain, reported by mkfilter */

#define ALARM	1		    /* bit in flags byte which controls alarm */

#define schar signed char
#define uchar unsigned char

extern schar adc;	/* input port; a read from here blocks until a new sample arrives  */
extern uchar flags;	/* alarm port */

main()
  { float xv[7], yv[7], zv[7];
    for (int i=1; i<7; i++) xv[i] = yv[i] = zv[i] = 0.0f;
    float pow = 1.0;
    flags = 0;
    for (;;)
      { /* shift */
	for (int i=0; i<6; i++)
	  { xv[i] = xv[i+1];
	    yv[i] = yv[i+1];
	    zv[i] = zv[i+1];
	  }
	/* read new input */
	xv[6] = (float) adc / (128.0f * NFGAIN * BPFGAIN);  /* scale so output is in range +-1 */
	/* notch filter step */
	yv[6] = (xv[4] + xv[6]) - 1.6180339887f * xv[5]
	      + ( -0.9937365126f * yv[4]) + ( 1.6129667210f * yv[5]);
	/* bandpass filter step */
	yv[6] = (xv[6] - xv[0]) + 3.0f * (xv[2] - xv[4])
	      + ( -0.3665332263f * yv[0]) + ( 2.5568853665f * yv[1])
	      + ( -7.4804429297f * yv[2]) + (11.7643366740f * yv[3])
	      + (-10.4822135550f * yv[4]) + ( 5.0079676498f * yv[5]);
	/* output stage */
	printf("%g\n", yv[6]);              /* highly unrealistic; a real prog wd e.g. plot on a chart recorder */
	float y2 = yv[6] * yv[6];	    /* instantaneous signal power */
	pow = (0.95f*pow) + (0.05f*y2);	    /* simple 1-pole lowpass filter; see Lecture 1 */
	if (pow < 0.2) flags |= ALARM;	    /* switch alarm on */
	else flags &= ~ALARM;		    /* switch alarm off */
      }
  }