timefilter: internally compute feedback factors.

The feedback factors for the timefilter are directly computed from
the expected period. This commit changes the init function to accept
the period itself and compute the feedback factors internally,
rather than having all client code duplicate the formulas.

This commit also actually fixes the formulas: the current code had
sqrt(2*o), but the correct formula, both theoretically and according
to experimental testing, is sqrt(2)*o.

Furthermore, it adds an exponential to feedback factors larger than
1 with large periods.

libavdevice/alsa-audio-dec.c

@@ -59,7 +59,6 @@ static av_cold int audio_read_header(AVFormatContext *s1)
     AVStream *st;
     int ret;
     enum CodecID codec_id;
-    double o;
 
     st = avformat_new_stream(s1, NULL);
     if (!st) {
@@ -81,9 +80,9 @@ static av_cold int audio_read_header(AVFormatContext *s1)
     st->codec->sample_rate = s->sample_rate;
     st->codec->channels    = s->channels;
     avpriv_set_pts_info(st, 64, 1, 1000000);  /* 64 bits pts in us */
-    o = 2 * M_PI * s->period_size / s->sample_rate * 1.5; // bandwidth: 1.5Hz
+    /* microseconds instead of seconds, MHz instead of Hz */
     s->timefilter = ff_timefilter_new(1000000.0 / s->sample_rate,
-                                      sqrt(2 * o), o * o);
+                                      s->period_size, 1.5E-6);
     if (!s->timefilter)
         goto fail;
 

libavdevice/jack_audio.c

@@ -145,7 +145,6 @@ static int start_jack(AVFormatContext *context)
     JackData *self = context->priv_data;
     jack_status_t status;
     int i, test;
-    double o, period;
 
     /* Register as a JACK client, using the context filename as client name. */
     self->client = jack_client_open(context->filename, JackNullOption, &status);
@@ -181,9 +180,7 @@ static int start_jack(AVFormatContext *context)
     jack_set_xrun_callback(self->client, xrun_callback, self);
 
     /* Create time filter */
-    period            = (double) self->buffer_size / self->sample_rate;
+    self->timefilter  = ff_timefilter_new (1.0 / self->sample_rate, self->buffer_size, 1.5);
-    o                 = 2 * M_PI * 1.5 * period; /// bandwidth: 1.5Hz
-    self->timefilter  = ff_timefilter_new (1.0 / self->sample_rate, sqrt(2 * o), o * o);
 
     /* Create FIFO buffers */
     self->filled_pkts = av_fifo_alloc(FIFO_PACKETS_NUM * sizeof(AVPacket));

libavdevice/timefilter.c

@@ -36,14 +36,21 @@ struct TimeFilter {
     int count;
 };
 
-TimeFilter *ff_timefilter_new(double clock_period,
+/* 1 - exp(-x) using a 3-order power series */
-                              double feedback2_factor,
+static double qexpneg(double x)
-                              double feedback3_factor)
+{
+    return 1 - 1 / (1 + x * (1 + x / 2 * (1 + x / 3)));
+}
+
+TimeFilter *ff_timefilter_new(double time_base,
+                              double period,
+                              double bandwidth)
 {
     TimeFilter *self       = av_mallocz(sizeof(TimeFilter));
-    self->clock_period     = clock_period;
+    double o               = 2 * M_PI * bandwidth * period * time_base;
-    self->feedback2_factor = feedback2_factor;
+    self->clock_period     = time_base;
-    self->feedback3_factor = feedback3_factor;
+    self->feedback2_factor = qexpneg(M_SQRT2 * o);
+    self->feedback3_factor = qexpneg(o * o);
     return self;
 }
 

libavdevice/timefilter.h

@@ -45,16 +45,18 @@ typedef struct TimeFilter TimeFilter;
  *
  * Unless you know what you are doing, you should set these as follow:
  *
- * o = 2 * M_PI * bandwidth * period
+ * o = 2 * M_PI * bandwidth * period_in_seconds
- * feedback2_factor = sqrt(2 * o)
+ * feedback2_factor = sqrt(2) * o
  * feedback3_factor = o * o
  *
  * Where bandwidth is up to you to choose. Smaller values will filter out more
  * of the jitter, but also take a longer time for the loop to settle. A good
  * starting point is something between 0.3 and 3 Hz.
  *
- * @param clock_period period of the hardware clock in seconds
+ * @param time_base   period of the hardware clock in seconds
- *        (for example 1.0/44100)
+ *                    (for example 1.0/44100)
+ * @param period      expected update interval, in input units
+ * @param brandwidth  filtering bandwidth, in Hz
  *
  * For more details about these parameters and background concepts please see:
  * http://www.kokkinizita.net/papers/usingdll.pdf