Example 9: A Signal Averager

The following program is considerably more complex than the previous examples in this chapter. It reads an annotation file (for which the annotator name is specified in its first argument, and the record name in the second argument) and selects beats of a specified type to be averaged. The program selects segments of the signals that are within 50 milliseconds of the time of the specified beat annotations, subtracts a baseline estimate from each sample, and calculates an average waveform (by default, the average normal QRS complex).

 1  #include <stdio.h>
 2  #include <ecg/db.h>
 3  #include <ecg/ecgmap.h>
 4
 5  main(argc, argv)
 6  int argc;
 7  char *argv[];
 8  {
 9      int btype, i, j, nbeats = 0, nsig, hwindow, window;
10      DB_Sample v[DB_MAXSIG], vb[DB_MAXSIG];
11      long stoptime = 0L, *sum[DB_MAXSIG];
12      DB_Anninfo a;
13      DB_Annotation annot;
14      static DB_Siginfo s[DB_MAXSIG];
15      void *calloc();
16
17      if (argc < 3) {
18          fprintf(stderr,
19                  "usage: %s annotator record [beat-type from to]\n",
20                  argv[0]);
21          exit(1);
22      }
23      a.name = argv[1]; a.stat = READ;
24      if ((nsig = dbinit(argv[2], &a, 1, s, DB_MAXSIG)) < 1) exit(2);
25      hwindow = strtim(".05"); window = 2*hwindow + 1;
26      btype = (argc > 3) ? strann(argv[3]) : NORMAL;
27      if (argc > 4) iannsettime(strtim(argv[4]));
28      if (argc > 5) {
29          if ((stoptime = strtim(argv[5])) < 0L)
30              stoptime = -stoptime;
31          if (s[0].nsamp > 0L && stoptime > s[0].nsamp)
32              stoptime = s[0].nsamp;
33      }
34      else stoptime = s[0].nsamp;
35      if (stoptime > 0L) stoptime -= hwindow;
36      for (i = 0; i < nsig; i++)
37          if ((sum[i]=(long *)calloc(window,sizeof(long))) == NULL) {
38              fprintf(stderr, "%s: insufficient memory\n", argv[0]);
39              exit(3);
40          }
41      while (getann(0, &annot) == 0 && annot.time < hwindow)
42          ;
43      do {
44          if (annot.anntyp != btype) continue;
45          isigsettime(annot.time - hwindow - 1);
46          (void)getvec(vb);
47          for (j = 0; j < window && getvec(v) > 0; j++)
48              for (i = 0; i < nsig; i++)
49                  sum[i][j] += v[i] - vb[i];
50          nbeats++;
51      } while (getann(0, &annot) == 0 &&
52               (stoptime == 0L || annot.time < stoptime));
53      if (nbeats < 1) {
54          fprintf(stderr, "%s: no `%s' beats found\n",
55                  argv[0], annstr(btype));
56          exit(4);
57      }
58      printf("Average of %d `%s' beats:\n", nbeats, annstr(btype));
59      for (j = 0; j < window; j++)
60          for (i = 0; i < nsig; i++)
61              printf("%g%c", (double)sum[i][j]/nbeats, 
62                     (i == nsig-1) ? '\n' : '\t');
63      exit(0);
64  }

Notes:

Line 25:

The "half-window" is 50 milliseconds wide, and the "window" (the duration of a segment to be entered into the average) is one sample more than twice that amount (i.e., 50 milliseconds to either side of the fiducial point defined by the annotation).

Line 26:

If a third argument is present on the command line, it is taken as an annotation code mnemonic for the desired beat type; otherwise, the program will average NORMAL QRS complexes.

Line 27:

If a fourth argument is present on the command line, it is taken as the start time; we arrange for the first annotation to be read by getann to be the first annotation that occurs after the chosen start time.

Lines 28--35:

This code similarly determines when the averaging should stop. Unless no stop time was specified on the command line and the signal length is not defined in the `header' file for the record, stoptime will have a positive value in line 35, which makes a tiny adjustment so that if a beat annotation occurs within 50 milliseconds of the end of the averaging period, the beat will not be included in the average.

Lines 36--40:

Here we allocate memory for the sum vectors that will be used to store the running totals. See K&R, page 167, for a description of calloc.

Lines 41-42:

This code addresses the (admittedly unlikely) prospect that the first annotation(s) might occur within the first 50 milliseconds of the record; any such annotations will be excluded from the average.

Lines 43--52:

Here we read annotations (the first is already in annot when we enter the loop, and subsequent annotations are read in line 51); select the desired ones (line 44); skip to the correct spot in the signals (line 45; the sample selected there is the one just before the beginning of the window); read a sample from each signal (line 46) into the vb vector, which will be used as a crude baseline estimate; read window samples from each signal (line 47), subtracting the baseline from each and adding the result into the running totals; update a beat counter (line 50); and check for loop termination conditions (line 52).

Lines 53--62:

This is the output section. If no beats of type btype were found, obviously no average can be printed; note that the message goes to the standard error output, so the user will notice it even if the standard output has been redirected to a file. In the usual case, the averages are printed out as a table, with a column allocated to each signal. Note the cast in line 61 (necessary to preserve precision), and the trick used in line 62 to print a tab after each column but the last in each line.

Go to the first, previous, next, last section, table of contents.