--- a
+++ b/HRV.src/pwr.c
@@ -0,0 +1,127 @@
+/* pwr.c		Joe Mietus		Oct 7 2008 */
+
+/*
+pwr :
+Usage: pwr [-r] [-L] lo1 hi1 ... [lo10 hi10]
+  Calculate total (and relative) power in fft between lo and hi.
+  options :
+  [-r] : print ratio of powers to total
+  [-L] : print powers on one line
+*/
+
+
+#include <stdio.h>
+#include <stdlib.h>
+
+#define MAXBANDS 10
+
+main(argc, argv)
+int argc;
+char *argv[];
+{
+
+    int i, n, nbands, rflag, lflag;
+    double lo[MAXBANDS], hi[MAXBANDS], pr[MAXBANDS], tot;
+    double freq[2], mag[2], hbin[2], pwr;
+
+
+    if (argc < 3) {
+        usage(argv[0]);
+        exit(1);
+    }
+
+    rflag = lflag = 0;
+    for (i=1, argc--; argv[i][0] == '-'; i++, argc--) {
+        switch(argv[i][1]) {
+        case 'r': rflag = 1;
+                  break;
+        case 'L': lflag = 1;
+                  break;
+        default: usage(argv[0]);
+                 exit(1);
+        }
+    }
+        
+    if (argc % 2 != 0) {
+        usage(argv[0]);
+        exit(1);
+    }
+
+    if ((nbands = argc/2) > 10) {
+        fprintf(stderr, "%s : max power bands = MAXBANDS\n");
+        exit(1);
+    }
+
+    for (n=0; n<nbands; n++) {
+        lo[n] = atof(argv[i++]);
+        hi[n] = atof(argv[i++]);
+        pr[n] = 0.0;
+    }
+
+    tot = 0.0;
+
+    if (scanf("%lf%lf", &freq[0], &mag[0]) != 2)
+        exit(2);
+    if (scanf("%lf%lf", &freq[1], &mag[1]) != 2)
+        exit(2);
+
+    pwr = mag[0]*mag[0];
+    tot += pwr;
+    for (n=0; n<nbands; n++) {
+        if (freq[0] >= lo[n] && freq[0] <= hi[n])
+            pr[n] += pwr;
+    }
+    freq[0] = freq[1];
+    mag[0] = mag[1];
+
+    while (scanf("%lf%lf", &freq[1], &mag[1]) == 2) {
+        pwr = mag[0]*mag[0];
+        tot += pwr;
+        for (n=0; n<nbands; n++) {
+            if (freq[0] >= lo[n] && freq[0] <= hi[n])
+                pr[n] += pwr;
+        }
+        freq[0] = freq[1];
+	mag[0] = mag[1];
+    }
+
+    pwr = mag[0]*mag[0];
+    tot += pwr;
+    for (n=0; n<nbands; n++) {
+        if (freq[0] >= lo[n] && freq[0] <= hi[n])
+        pr[n] += pwr;
+    }
+
+    if (lflag) {
+        printf("%g ", tot);
+        for (n=0; n<nbands; n++) {
+            printf("%g ", pr[n]);
+            if (rflag)
+                printf("(%.3f) ", pr[n]/tot);
+        }
+        printf("\n");
+    }
+    else {
+        printf("Total = %g\n", tot);
+        for (n=0; n<nbands; n++) {
+            printf("%g - %g = %g", lo[n], hi[n], pr[n]);
+            if (rflag)
+                printf(" (%.3f)\n", pr[n]/tot);
+            else
+                printf("\n");
+        }
+    }
+}
+
+
+usage(prog)
+char *prog;
+{
+    fprintf(stderr, "Usage: %s [-r] [-L] lo1 hi1 ... [lo10 hi10]\n", prog);
+    fprintf(stderr, " Calculate total (and relative) power");
+    fprintf(stderr, " in fft between lo and hi.\n");
+    fprintf(stderr, " options :\n");
+    fprintf(stderr, " [-r] : print ratio of powers to total\n");
+    fprintf(stderr, " [-L] : print powers on one line\n");
+
+}