#include #include #include #include #include #include #include #include #define FREQUENCY 20 int compare(const void *a, const void *b) { return *(double*)a < *(double*)b ? -1 : *(double*)a > *(double*)b; } void sortabs(double *output, const double *input, unsigned int length) { int i; for (i = 0; i < length; i++) output[i] = input[i] > 0 ? input[i] : -input[i]; qsort(output, length, sizeof(double), &compare); } #define FILTERSIZE 2048 static double x[FILTERSIZE+1]; static double y[FILTERSIZE+1]; static double xh[FILTERSIZE+1]; static double yh[FILTERSIZE+1]; static double st[FILTERSIZE+1]; /* Designed with http://www-users.cs.york.ac.uk/~fisher/mkfilter/trad.html */ double lowpassfilter(double input) { double max = 0, res; double *xv = x + FILTERSIZE - 6; double *yv = y + FILTERSIZE - 6; memmove(x, x+1, sizeof(double)*FILTERSIZE); memmove(y, y+1, sizeof(double)*FILTERSIZE); #if 1 sortabs(st, y, FILTERSIZE); max = st[FILTERSIZE/2] * 4; #endif if (max == 0) max = 32767; /* 0.125 Hz */ xv[6] = input / 1.881202115e+10; yv[6] = (xv[0] + xv[6]) + 6 * (xv[1] + xv[5]) + 15 * (xv[2] + xv[4]) + 20 * xv[3] + ( -0.8592119586 * yv[0]) + ( 5.2856359974 * yv[1]) + (-13.5499271760 * yv[2]) + ( 18.5280691620 * yv[3]) + (-14.2528406660 * yv[4]) + ( 5.8482746376 * yv[5]); res = yv[6]; if (res > max) res = max; if (res < -max) res = -max; return res; } /* Designed with http://www-users.cs.york.ac.uk/~fisher/mkfilter/trad.html */ double bandpassfilter(double input) { double max = 0, res; double *xv = x + FILTERSIZE - 14; double *yv = y + FILTERSIZE - 14; memmove(x, x+1, sizeof(double)*FILTERSIZE); memmove(y, y+1, sizeof(double)*FILTERSIZE); #if 1 sortabs(st, y, FILTERSIZE); max = st[FILTERSIZE/2] * 4; #endif if (max == 0) max = 32767; /* 0.005 - 0.125 Hz */ xv[14] = input / 7.955649637e+02; yv[14] = (xv[14] - xv[0]) + 7 * (xv[2] - xv[12]) + 21 * (xv[10] - xv[4]) + 35 * (xv[6] - xv[8]) + ( -0.8441448498 * yv[0]) + ( 11.9606713110 * yv[1]) + (-78.6840116620 * yv[2]) + (318.5452188700 * yv[3]) + (-886.6156672800 * yv[4]) + (1794.7502032000 * yv[5]) + (-2724.8391359000 * yv[6]) + (3151.9940083000 * yv[7]) + (-2791.5985237000 * yv[8]) + (1883.7712556000 * yv[9]) + (-953.3921265300 * yv[10]) + (350.9289275400 * yv[11]) + (-88.8068314720 * yv[12]) + ( 13.8301565290 * yv[13]); res = yv[14]; if (res > max) res = max; if (res < -max) res = -max; return res; } /* Designed with http://www-users.cs.york.ac.uk/~fisher/mkfilter/trad.html */ double highpassfilter(double input) { double *xv = xh + FILTERSIZE - 6; double *yv = yh + FILTERSIZE - 6; memmove(xh, xh+1, sizeof(double)*FILTERSIZE); memmove(yh, yh+1, sizeof(double)*FILTERSIZE); /* 0.032 Hz */ xv[6] = input / 1.019611013; yv[6] = (xv[0] + xv[6]) - 6 * (xv[1] + xv[5]) + 15 * (xv[2] + xv[4]) - 20 * xv[3] + ( -0.9619023031 * yv[0]) + ( 5.8087761720 * yv[1]) + (-14.6160726220 * yv[2]) + ( 19.6145837910 * yv[3]) + (-14.8065428880 * yv[4]) + ( 5.9611578499 * yv[5]); return yv[6]; } int16_t filter(int16_t *input, int16_t *output, int n) { int i; for (i = 0; i < n; i++) { double f; #if 0 f = highpassfilter(lowpassfilter((double)input[i])); #else f = lowpassfilter((double)input[i]); #endif if (f > 32767) f = 32767; if (f < -32767) f = -32767; output[i] = (int16_t)f; } } void usage(const char **argv) { printf("Usage: %s [-d] [log interval] [output file] [input device] [raw output file]\n", argv[0]); printf(" -d = deamon mode\n"); printf(" Default values:\n"); printf(" log interval: 60 (seconds)\n"); printf(" output file: stdout\n"); printf(" input device: stdin\n"); } void printlog(FILE *output, time_t ts, int16_t max, int16_t min) { char timestamp[15]; struct tm* scratch_time = gmtime(&ts); strftime(timestamp, sizeof(timestamp), "%Y%m%d%H%M%S", scratch_time); fprintf(output, "%s\t%d\t%d\n", timestamp, (int)max, (int)min); } /* Log raw data: 16 bit big endian. 0x8000 marks that a 32 bit big endian timestamp (seconds since 1970-01-01 00:00:00 UTC) follows */ void lograw(FILE *output, int16_t *log, unsigned int n) { int i = 0; uint32_t t = (uint32_t)time(0) - n / FREQUENCY; fputc(0x80, output); fputc(0x00, output); fputc((t >> 24) & 0xff, output); fputc((t >> 16) & 0xff, output); fputc((t >> 8) & 0xff, output); fputc((t >> 0) & 0xff, output); while (i < n) { int16_t v = log[i] == -32768 ? -32767 : log[i]; fputc((v >> 8) & 0xff, output); fputc((v >> 0) & 0xff, output); i++; } } void process(const int16_t *log, unsigned int i, int16_t *max, int16_t *min) { *max = -32768; *min = 32767; while (i--) { if (*max < log[i]) *max = log[i]; if (*min > log[i]) *min = log[i]; } } int main(int argc, const char** argv) { FILE *input; FILE *output; FILE *rawoutput = 0; unsigned int interval; time_t logtime; int16_t *log, *filtered; unsigned int logindex; const char **cmdline = argv; int forked = 0; int i = 0; int level; for (i = 0; i <= FILTERSIZE; i++) x[i] = y[i] = xh[i] = yh[i] = 0.0; while (argc > 1 && argv[1] && *argv[1] == '-') { switch (argv[1][1]) { case 'd': if (!forked) { forked = 1; umask(0); setpgrp(); if (fork()) return 0; } break; default: printf("%s: Unknown option -%c\n", cmdline[0], argv[1][1]); usage(cmdline); return -1; } argc--; argv++; } input = argc > 3 && argv[3] ? fopen(argv[3], "r") : stdin; interval = argc > 1 && argv[1] ? atoi(argv[1]) : 60; logtime = time(NULL); logindex = 0; log = malloc(FREQUENCY * interval * sizeof(int16_t) + sizeof(int16_t)); filtered = malloc(FREQUENCY * interval * sizeof(int16_t) + sizeof(int16_t)); if (!log) { printf("%s: Could allocate memory.\n", cmdline[0]); return -1; } if (!input) { printf("%s: Could not open input.\n", cmdline[0]); return -1; } if (!(argc > 2 && argv[2])) output = stdout; fscanf(input, "%d", &level); while (1) { time_t now; fscanf(input, "%d", &level); log[logindex] = level; logindex += logindex <= FREQUENCY * interval; if (logindex > FREQUENCY * interval) memmove(log, log + 1, FREQUENCY * interval * sizeof(int16_t)); now = time(NULL); if (now - logtime >= interval) { int i; double avg; int16_t max, min; rawoutput = argc > 4 && argv[4] ? fopen(argv[4], "a") : 0; if (rawoutput) { lograw(rawoutput, log, logindex); fclose(rawoutput); } filter(log, filtered, logindex); process(filtered, logindex, &max, &min); if (output != stdout) output = fopen(argv[2], "a"); if (output) { printlog(output, now/2 + logtime/2, max, min); if (output != stdout) fclose(output); } logindex = 0; logtime = now; } } fclose(input); return 0; }