git.sesse.net Git - audiosync/blob - estimate-skew.c

   1 #include <stdio.h>
   2 #include <string.h>
   3 #include <math.h>
   4 #include <stdlib.h>
   5
   6 double frac(double x)
   7 {
   8         return x - (int)x;
   9 }
  10
  11 double filter(double x)
  12 {
  13         static double l = 0.0;
  14         l = 0.001 * x + 0.999 * l;
  15         return l;
  16 }
  17
  18 int main(int argc, char **argv)
  19 {
  20         FILE *in1, *in2, *skew;
  21
  22         // open input1 (reference)
  23         if (strcmp(argv[1], "-") == 0) {
  24                 in1 = stdin;
  25         } else {
  26                 in1 = fopen(argv[1], "rb");
  27                 if (in1 == NULL) {
  28                         perror(argv[1]);
  29                         exit(1);
  30                 }
  31         }
  32
  33         // open input2
  34         if (strcmp(argv[2], "-") == 0) {
  35                 in2 = stdin;
  36         } else {
  37                 in2 = fopen(argv[2], "rb");
  38                 if (in2 == NULL) {
  39                         perror(argv[2]);
  40                         exit(1);
  41                 }
  42         }
  43
  44         // open (estimated) skew
  45         if (strcmp(argv[3], "-") == 0) {
  46                 skew = stderr;
  47         } else {
  48                 skew = fopen(argv[3], "wb");
  49                 if (skew == NULL) {
  50                         perror(argv[3]);
  51                         exit(1);
  52                 }
  53         }
  54
  55         short prev_sample, sample = 0;
  56         int in_pos = -1;
  57         double p = 0.0;
  58         double speed = 1.001;
  59
  60         while (!feof(in1) && !feof(in2)) {
  61                 short refs;
  62
  63                 // read sample from reference
  64                 if (fread(&refs, sizeof(short), 1, in1) != 1) {
  65                         exit(0);
  66                 }
  67
  68                 // read sample from in2 (at current speed)
  69                 p += speed;
  70                 while ((int)(ceil(p)) > in_pos) {
  71                         prev_sample = sample;
  72                         if (fread(&sample, sizeof(short), 1, in2) != 1) {
  73                                 exit(0);
  74                         }
  75                         ++in_pos;
  76                 }
  77
  78                 // linear interpolation (works well since delta_p varies so slowly)
  79                 double t = frac(p);
  80                 double intp_sample = prev_sample * (1.0 - t) + sample * t;
  81
  82                 // subtract the two samples and pull it through a filter. we assume
  83                 // the sines are normalized for now (and that there's no dc skew)
  84                 double offset = filter(intp_sample - refs);
  85                 speed += 1e-8 * offset;
  86                 fwrite(&speed, sizeof(double), 1, skew);
  87
  88                 printf("%f (offset=%f / filt=%f)\n", speed, intp_sample - refs, offset);
  89         }
  90 }