Browse code

avfilter: add audio emphasis filter

Signed-off-by: Paul B Mahol <onemda@gmail.com>

Paul B Mahol authored on 2015/11/30 21:36:58
Showing 7 changed files
... ...
@@ -40,6 +40,7 @@ version <next>:
40 40
 - apulsator filter
41 41
 - sidechaingate audio filter
42 42
 - mipsdspr1 option has been renamed to mipsdsp
43
+- aemphasis filter
43 44
 
44 45
 
45 46
 version 2.8:
... ...
@@ -1051,6 +1051,21 @@ int main(void){ $func(); }
1051 1051
 EOF
1052 1052
 }
1053 1053
 
1054
+check_complexfunc(){
1055
+    log check_complexfunc "$@"
1056
+    func=$1
1057
+    narg=$2
1058
+    shift 2
1059
+    test $narg = 2 && args="f, g" || args="f * I"
1060
+    disable $func
1061
+    check_ld "cc" "$@" <<EOF && enable $func
1062
+#include <complex.h>
1063
+#include <math.h>
1064
+float foo(complex float f, complex float g) { return $func($args); }
1065
+int main(void){ return (int) foo; }
1066
+EOF
1067
+}
1068
+
1054 1069
 check_mathfunc(){
1055 1070
     log check_mathfunc "$@"
1056 1071
     func=$1
... ...
@@ -1768,6 +1783,11 @@ INTRINSICS_LIST="
1768 1768
     intrinsics_neon
1769 1769
 "
1770 1770
 
1771
+COMPLEX_FUNCS="
1772
+    cabs
1773
+    cexp
1774
+"
1775
+
1771 1776
 MATH_FUNCS="
1772 1777
     atanf
1773 1778
     atan2f
... ...
@@ -1903,6 +1923,7 @@ HAVE_LIST="
1903 1903
     $ARCH_FEATURES
1904 1904
     $ATOMICS_LIST
1905 1905
     $BUILTIN_LIST
1906
+    $COMPLEX_FUNCS
1906 1907
     $HAVE_LIST_CMDLINE
1907 1908
     $HAVE_LIST_PUB
1908 1909
     $HEADERS_LIST
... ...
@@ -2785,6 +2806,7 @@ unix_protocol_deps="sys_un_h"
2785 2785
 unix_protocol_select="network"
2786 2786
 
2787 2787
 # filters
2788
+aemphasis_filter_deps="cabs cexp"
2788 2789
 amovie_filter_deps="avcodec avformat"
2789 2790
 aresample_filter_deps="swresample"
2790 2791
 ass_filter_deps="libass"
... ...
@@ -5324,6 +5346,10 @@ for func in $MATH_FUNCS; do
5324 5324
     eval check_mathfunc $func \${${func}_args:-1}
5325 5325
 done
5326 5326
 
5327
+for func in $COMPLEX_FUNCS; do
5328
+    eval check_complexfunc $func \${${func}_args:-1}
5329
+done
5330
+
5327 5331
 # these are off by default, so fail if requested and not available
5328 5332
 enabled avfoundation_indev && { check_header_oc AVFoundation/AVFoundation.h || disable avfoundation_indev; }
5329 5333
 enabled avfoundation_indev && { check_lib2 CoreGraphics/CoreGraphics.h CGGetActiveDisplayList -framework CoreGraphics ||
... ...
@@ -528,6 +528,52 @@ aecho=0.8:0.9:1000|1800:0.3|0.25
528 528
 @end example
529 529
 @end itemize
530 530
 
531
+@section aemphasis
532
+Audio emphasis filter creates or restores material directly taken from LPs or
533
+emphased CDs with different filter curves. E.g. to store music on vinyl the
534
+signal has to be altered by a filter first to even out the disadvantages of
535
+this recording medium.
536
+Once the material is played back the inverse filter has to be applied to
537
+restore the distortion of the frequency response.
538
+
539
+The filter accepts the following options:
540
+
541
+@table @option
542
+@item level_in
543
+Set input gain.
544
+
545
+@item level_out
546
+Set output gain.
547
+
548
+@item mode
549
+Set filter mode. For restoring material use @code{reproduction} mode, otherwise
550
+use @code{production} mode. Default is @code{reproduction} mode.
551
+
552
+@item type
553
+Set filter type. Selects medium. Can be one of the following:
554
+
555
+@table @option
556
+@item col
557
+select Columbia.
558
+@item emi
559
+select EMI.
560
+@item bsi
561
+select BSI (78RPM).
562
+@item riaa
563
+select RIAA.
564
+@item cd
565
+select Compact Disc (CD).
566
+@item 50fm
567
+select 50µs (FM).
568
+@item 75fm
569
+select 75µs (FM).
570
+@item 50kf
571
+select 50µs (FM-KF).
572
+@item 75kf
573
+select 75µs (FM-KF).
574
+@end table
575
+@end table
576
+
531 577
 @section aeval
532 578
 
533 579
 Modify an audio signal according to the specified expressions.
... ...
@@ -27,6 +27,7 @@ OBJS-$(CONFIG_ACOMPRESSOR_FILTER)            += af_sidechaincompress.o
27 27
 OBJS-$(CONFIG_ACROSSFADE_FILTER)             += af_afade.o
28 28
 OBJS-$(CONFIG_ADELAY_FILTER)                 += af_adelay.o
29 29
 OBJS-$(CONFIG_AECHO_FILTER)                  += af_aecho.o
30
+OBJS-$(CONFIG_AEMPHASIS_FILTER)              += af_aemphasis.o
30 31
 OBJS-$(CONFIG_AEVAL_FILTER)                  += aeval.o
31 32
 OBJS-$(CONFIG_AFADE_FILTER)                  += af_afade.o
32 33
 OBJS-$(CONFIG_AFORMAT_FILTER)                += af_aformat.o
33 34
new file mode 100644
... ...
@@ -0,0 +1,370 @@
0
+/*
1
+ * Copyright (c) 2001-2010 Krzysztof Foltman, Markus Schmidt, Thor Harald Johansen, Damien Zammit and others
2
+ *
3
+ * This file is part of FFmpeg.
4
+ *
5
+ * FFmpeg is free software; you can redistribute it and/or
6
+ * modify it under the terms of the GNU Lesser General Public
7
+ * License as published by the Free Software Foundation; either
8
+ * version 2.1 of the License, or (at your option) any later version.
9
+ *
10
+ * FFmpeg is distributed in the hope that it will be useful,
11
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
12
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
13
+ * Lesser General Public License for more details.
14
+ *
15
+ * You should have received a copy of the GNU Lesser General Public
16
+ * License along with FFmpeg; if not, write to the Free Software
17
+ * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18
+ */
19
+
20
+#include <complex.h>
21
+
22
+#include "libavutil/opt.h"
23
+#include "avfilter.h"
24
+#include "internal.h"
25
+#include "audio.h"
26
+
27
+typedef struct BiquadCoeffs {
28
+    double a0, a1, a2, b1, b2;
29
+} BiquadCoeffs;
30
+
31
+typedef struct BiquadD2 {
32
+    double a0, a1, a2, b1, b2, w1, w2;
33
+} BiquadD2;
34
+
35
+typedef struct RIAACurve {
36
+    BiquadD2 r1;
37
+    BiquadD2 brickw;
38
+    int use_brickw;
39
+} RIAACurve;
40
+
41
+typedef struct AudioEmphasisContext {
42
+    const AVClass *class;
43
+    int mode, type;
44
+    double level_in, level_out;
45
+
46
+    RIAACurve *rc;
47
+} AudioEmphasisContext;
48
+
49
+#define OFFSET(x) offsetof(AudioEmphasisContext, x)
50
+#define FLAGS AV_OPT_FLAG_AUDIO_PARAM|AV_OPT_FLAG_FILTERING_PARAM
51
+
52
+static const AVOption aemphasis_options[] = {
53
+    { "level_in",      "set input gain", OFFSET(level_in),  AV_OPT_TYPE_DOUBLE, {.dbl=1}, 0, 64, FLAGS },
54
+    { "level_out",    "set output gain", OFFSET(level_out), AV_OPT_TYPE_DOUBLE, {.dbl=1}, 0, 64, FLAGS },
55
+    { "mode",         "set filter mode", OFFSET(mode), AV_OPT_TYPE_INT,   {.i64=0}, 0, 1, FLAGS, "mode" },
56
+    { "reproduction",              NULL,            0, AV_OPT_TYPE_CONST, {.i64=0}, 0, 0, FLAGS, "mode" },
57
+    { "production",                NULL,            0, AV_OPT_TYPE_CONST, {.i64=1}, 0, 0, FLAGS, "mode" },
58
+    { "type",         "set filter type", OFFSET(type), AV_OPT_TYPE_INT,   {.i64=4}, 0, 8, FLAGS, "type" },
59
+    { "col",                 "Columbia",            0, AV_OPT_TYPE_CONST, {.i64=0}, 0, 0, FLAGS, "type" },
60
+    { "emi",                      "EMI",            0, AV_OPT_TYPE_CONST, {.i64=1}, 0, 0, FLAGS, "type" },
61
+    { "bsi",              "BSI (78RPM)",            0, AV_OPT_TYPE_CONST, {.i64=2}, 0, 0, FLAGS, "type" },
62
+    { "riaa",                    "RIAA",            0, AV_OPT_TYPE_CONST, {.i64=3}, 0, 0, FLAGS, "type" },
63
+    { "cd",         "Compact Disc (CD)",            0, AV_OPT_TYPE_CONST, {.i64=4}, 0, 0, FLAGS, "type" },
64
+    { "50fm",               "50µs (FM)",            0, AV_OPT_TYPE_CONST, {.i64=5}, 0, 0, FLAGS, "type" },
65
+    { "75fm",               "75µs (FM)",            0, AV_OPT_TYPE_CONST, {.i64=6}, 0, 0, FLAGS, "type" },
66
+    { "50kf",            "50µs (FM-KF)",            0, AV_OPT_TYPE_CONST, {.i64=7}, 0, 0, FLAGS, "type" },
67
+    { "75kf",            "75µs (FM-KF)",            0, AV_OPT_TYPE_CONST, {.i64=8}, 0, 0, FLAGS, "type" },
68
+    { NULL }
69
+};
70
+
71
+AVFILTER_DEFINE_CLASS(aemphasis);
72
+
73
+static inline double biquad(BiquadD2 *bq, double in)
74
+{
75
+    double n = in;
76
+    double tmp = n - bq->w1 * bq->b1 - bq->w2 * bq->b2;
77
+    double out = tmp * bq->a0 + bq->w1 * bq->a1 + bq->w2 * bq->a2;
78
+
79
+    bq->w2 = bq->w1;
80
+    bq->w1 = tmp;
81
+
82
+    return out;
83
+}
84
+
85
+static int filter_frame(AVFilterLink *inlink, AVFrame *in)
86
+{
87
+    AVFilterContext *ctx = inlink->dst;
88
+    AVFilterLink *outlink = ctx->outputs[0];
89
+    AudioEmphasisContext *s = ctx->priv;
90
+    const double *src = (const double *)in->data[0];
91
+    const double level_out = s->level_out;
92
+    const double level_in = s->level_in;
93
+    AVFrame *out;
94
+    double *dst;
95
+    int n, c;
96
+
97
+    if (av_frame_is_writable(in)) {
98
+        out = in;
99
+    } else {
100
+        out = ff_get_audio_buffer(inlink, in->nb_samples);
101
+        if (!out) {
102
+            av_frame_free(&in);
103
+            return AVERROR(ENOMEM);
104
+        }
105
+        av_frame_copy_props(out, in);
106
+    }
107
+    dst = (double *)out->data[0];
108
+
109
+    for (n = 0; n < in->nb_samples; n++) {
110
+        for (c = 0; c < inlink->channels; c++)
111
+            dst[c] = level_out * biquad(&s->rc[c].r1, s->rc[c].use_brickw ? biquad(&s->rc[c].brickw, src[c] * level_in) : src[c] * level_in);
112
+        dst += inlink->channels;
113
+        src += inlink->channels;
114
+    }
115
+
116
+    if (in != out)
117
+        av_frame_free(&in);
118
+    return ff_filter_frame(outlink, out);
119
+}
120
+
121
+static int query_formats(AVFilterContext *ctx)
122
+{
123
+    AVFilterChannelLayouts *layouts;
124
+    AVFilterFormats *formats;
125
+    static const enum AVSampleFormat sample_fmts[] = {
126
+        AV_SAMPLE_FMT_DBL,
127
+        AV_SAMPLE_FMT_NONE
128
+    };
129
+    int ret;
130
+
131
+    layouts = ff_all_channel_counts();
132
+    if (!layouts)
133
+        return AVERROR(ENOMEM);
134
+    ret = ff_set_common_channel_layouts(ctx, layouts);
135
+    if (ret < 0)
136
+        return ret;
137
+
138
+    formats = ff_make_format_list(sample_fmts);
139
+    if (!formats)
140
+        return AVERROR(ENOMEM);
141
+    ret = ff_set_common_formats(ctx, formats);
142
+    if (ret < 0)
143
+        return ret;
144
+
145
+    formats = ff_all_samplerates();
146
+    if (!formats)
147
+        return AVERROR(ENOMEM);
148
+    return ff_set_common_samplerates(ctx, formats);
149
+}
150
+
151
+static inline void set_highshelf_rbj(BiquadD2 *bq, double freq, double q, double peak, double sr)
152
+{
153
+    double A = sqrt(peak);
154
+    double w0 = freq * 2 * M_PI / sr;
155
+    double alpha = sin(w0) / (2 * q);
156
+    double cw0 = cos(w0);
157
+    double tmp = 2 * sqrt(A) * alpha;
158
+    double b0 = 0, ib0 = 0;
159
+
160
+    bq->a0 =    A*( (A+1) + (A-1)*cw0 + tmp);
161
+    bq->a1 = -2*A*( (A-1) + (A+1)*cw0);
162
+    bq->a2 =    A*( (A+1) + (A-1)*cw0 - tmp);
163
+        b0 =        (A+1) - (A-1)*cw0 + tmp;
164
+    bq->b1 =    2*( (A-1) - (A+1)*cw0);
165
+    bq->b2 =        (A+1) - (A-1)*cw0 - tmp;
166
+
167
+    ib0     = 1 / b0;
168
+    bq->b1 *= ib0;
169
+    bq->b2 *= ib0;
170
+    bq->a0 *= ib0;
171
+    bq->a1 *= ib0;
172
+    bq->a2 *= ib0;
173
+}
174
+
175
+static inline void set_lp_rbj(BiquadD2 *bq, double fc, double q, double sr, double gain)
176
+{
177
+    double omega = 2.0 * M_PI * fc / sr;
178
+    double sn = sin(omega);
179
+    double cs = cos(omega);
180
+    double alpha = sn/(2 * q);
181
+    double inv = 1.0/(1.0 + alpha);
182
+
183
+    bq->a2 = bq->a0 = gain * inv * (1.0 - cs) * 0.5;
184
+    bq->a1 = bq->a0 + bq->a0;
185
+    bq->b1 = (-2.0 * cs * inv);
186
+    bq->b2 = ((1.0 - alpha) * inv);
187
+}
188
+
189
+static double freq_gain(BiquadCoeffs *c, double freq, double sr)
190
+{
191
+    double complex z, w;
192
+
193
+    freq *= 2.0 * M_PI / sr;
194
+    w = 0 + I * freq;
195
+    z = 1.0 / cexp(w);
196
+
197
+    return cabs(((double complex)c->a0 + c->a1 * z + c->a2 * z*z) /
198
+                ((double complex)1.0 + c->b1 * z + c->b2 * z*z));
199
+}
200
+
201
+static int config_input(AVFilterLink *inlink)
202
+{
203
+    double i, j, k, g, t, a0, a1, a2, b1, b2, tau1, tau2, tau3;
204
+    double cutfreq, gain1kHz, gc, sr = inlink->sample_rate;
205
+    AVFilterContext *ctx = inlink->dst;
206
+    AudioEmphasisContext *s = ctx->priv;
207
+    BiquadCoeffs coeffs;
208
+    int ch;
209
+
210
+    s->rc = av_calloc(inlink->channels, sizeof(*s->rc));
211
+    if (!s->rc)
212
+        return AVERROR(ENOMEM);
213
+
214
+    switch (s->type) {
215
+    case 0: //"Columbia"
216
+        i = 100.;
217
+        j = 500.;
218
+        k = 1590.;
219
+        break;
220
+    case 1: //"EMI"
221
+        i = 70.;
222
+        j = 500.;
223
+        k = 2500.;
224
+        break;
225
+    case 2: //"BSI(78rpm)"
226
+        i = 50.;
227
+        j = 353.;
228
+        k = 3180.;
229
+        break;
230
+    case 3: //"RIAA"
231
+    default:
232
+        tau1 = 0.003180;
233
+        tau2 = 0.000318;
234
+        tau3 = 0.000075;
235
+        i = 1. / (2. * M_PI * tau1);
236
+        j = 1. / (2. * M_PI * tau2);
237
+        k = 1. / (2. * M_PI * tau3);
238
+        break;
239
+    case 4: //"CD Mastering"
240
+        tau1 = 0.000050;
241
+        tau2 = 0.000015;
242
+        tau3 = 0.0000001;// 1.6MHz out of audible range for null impact
243
+        i = 1. / (2. * M_PI * tau1);
244
+        j = 1. / (2. * M_PI * tau2);
245
+        k = 1. / (2. * M_PI * tau3);
246
+        break;
247
+    case 5: //"50µs FM (Europe)"
248
+        tau1 = 0.000050;
249
+        tau2 = tau1 / 20;// not used
250
+        tau3 = tau1 / 50;//
251
+        i = 1. / (2. * M_PI * tau1);
252
+        j = 1. / (2. * M_PI * tau2);
253
+        k = 1. / (2. * M_PI * tau3);
254
+        break;
255
+    case 6: //"75µs FM (US)"
256
+        tau1 = 0.000075;
257
+        tau2 = tau1 / 20;// not used
258
+        tau3 = tau1 / 50;//
259
+        i = 1. / (2. * M_PI * tau1);
260
+        j = 1. / (2. * M_PI * tau2);
261
+        k = 1. / (2. * M_PI * tau3);
262
+        break;
263
+    }
264
+
265
+    i *= 2 * M_PI;
266
+    j *= 2 * M_PI;
267
+    k *= 2 * M_PI;
268
+
269
+    t = 1. / sr;
270
+
271
+    //swap a1 b1, a2 b2
272
+    if (s->type == 7 || s->type == 8) {
273
+        s->rc[0].use_brickw = 0;
274
+        double tau = (s->type == 7 ? 0.000050 : 0.000075);
275
+        double f = 1.0 / (2 * M_PI * tau);
276
+        double nyq = sr * 0.5;
277
+        double gain = sqrt(1.0 + nyq * nyq / (f * f)); // gain at Nyquist
278
+        double cfreq = sqrt((gain - 1.0) * f * f); // frequency
279
+        double q = 1.0;
280
+
281
+        if (s->type == 8)
282
+            q = pow((sr / 3269.0) + 19.5, -0.25); // somewhat poor curve-fit
283
+        if (s->type == 7)
284
+            q = pow((sr / 4750.0) + 19.5, -0.25);
285
+        if (s->mode == 0)
286
+            set_highshelf_rbj(&s->rc[0].r1, cfreq, q, 1. / gain, sr);
287
+        else
288
+            set_highshelf_rbj(&s->rc[0].r1, cfreq, q, gain, sr);
289
+    } else {
290
+        s->rc[0].use_brickw = 1;
291
+        if (s->mode == 0) { // Reproduction
292
+            g  = 1. / (4.+2.*i*t+2.*k*t+i*k*t*t);
293
+            a0 = (2.*t+j*t*t)*g;
294
+            a1 = (2.*j*t*t)*g;
295
+            a2 = (-2.*t+j*t*t)*g;
296
+            b1 = (-8.+2.*i*k*t*t)*g;
297
+            b2 = (4.-2.*i*t-2.*k*t+i*k*t*t)*g;
298
+        } else {  // Production
299
+            g  = 1. / (2.*t+j*t*t);
300
+            a0 = (4.+2.*i*t+2.*k*t+i*k*t*t)*g;
301
+            a1 = (-8.+2.*i*k*t*t)*g;
302
+            a2 = (4.-2.*i*t-2.*k*t+i*k*t*t)*g;
303
+            b1 = (2.*j*t*t)*g;
304
+            b2 = (-2.*t+j*t*t)*g;
305
+        }
306
+
307
+        coeffs.a0 = a0;
308
+        coeffs.a1 = a1;
309
+        coeffs.a2 = a2;
310
+        coeffs.b1 = b1;
311
+        coeffs.b2 = b2;
312
+
313
+        // the coeffs above give non-normalized value, so it should be normalized to produce 0dB at 1 kHz
314
+        // find actual gain
315
+        // Note: for FM emphasis, use 100 Hz for normalization instead
316
+        gain1kHz = freq_gain(&coeffs, 1000.0, sr);
317
+        // divide one filter's x[n-m] coefficients by that value
318
+        gc = 1.0 / gain1kHz;
319
+        s->rc[0].r1.a0 = coeffs.a0 * gc;
320
+        s->rc[0].r1.a1 = coeffs.a1 * gc;
321
+        s->rc[0].r1.a2 = coeffs.a2 * gc;
322
+        s->rc[0].r1.b1 = coeffs.b1;
323
+        s->rc[0].r1.b2 = coeffs.b2;
324
+    }
325
+
326
+    cutfreq = FFMIN(0.45 * sr, 21000.);
327
+    set_lp_rbj(&s->rc[0].brickw, cutfreq, 0.707, sr, 1.);
328
+
329
+    for (ch = 1; ch < inlink->channels; ch++) {
330
+        memcpy(&s->rc[ch], &s->rc[0], sizeof(RIAACurve));
331
+    }
332
+
333
+    return 0;
334
+}
335
+
336
+static av_cold void uninit(AVFilterContext *ctx)
337
+{
338
+    AudioEmphasisContext *s = ctx->priv;
339
+    av_freep(&s->rc);
340
+}
341
+
342
+static const AVFilterPad avfilter_af_aemphasis_inputs[] = {
343
+    {
344
+        .name         = "default",
345
+        .type         = AVMEDIA_TYPE_AUDIO,
346
+        .config_props = config_input,
347
+        .filter_frame = filter_frame,
348
+    },
349
+    { NULL }
350
+};
351
+
352
+static const AVFilterPad avfilter_af_aemphasis_outputs[] = {
353
+    {
354
+        .name = "default",
355
+        .type = AVMEDIA_TYPE_AUDIO,
356
+    },
357
+    { NULL }
358
+};
359
+
360
+AVFilter ff_af_aemphasis = {
361
+    .name          = "aemphasis",
362
+    .description   = NULL_IF_CONFIG_SMALL("Audio emphasis."),
363
+    .priv_size     = sizeof(AudioEmphasisContext),
364
+    .priv_class    = &aemphasis_class,
365
+    .uninit        = uninit,
366
+    .query_formats = query_formats,
367
+    .inputs        = avfilter_af_aemphasis_inputs,
368
+    .outputs       = avfilter_af_aemphasis_outputs,
369
+};
... ...
@@ -49,6 +49,7 @@ void avfilter_register_all(void)
49 49
     REGISTER_FILTER(ACROSSFADE,     acrossfade,     af);
50 50
     REGISTER_FILTER(ADELAY,         adelay,         af);
51 51
     REGISTER_FILTER(AECHO,          aecho,          af);
52
+    REGISTER_FILTER(AEMPHASIS,      aemphasis,      af);
52 53
     REGISTER_FILTER(AEVAL,          aeval,          af);
53 54
     REGISTER_FILTER(AFADE,          afade,          af);
54 55
     REGISTER_FILTER(AFORMAT,        aformat,        af);
... ...
@@ -30,7 +30,7 @@
30 30
 #include "libavutil/version.h"
31 31
 
32 32
 #define LIBAVFILTER_VERSION_MAJOR   6
33
-#define LIBAVFILTER_VERSION_MINOR  19
33
+#define LIBAVFILTER_VERSION_MINOR  20
34 34
 #define LIBAVFILTER_VERSION_MICRO 100
35 35
 
36 36
 #define LIBAVFILTER_VERSION_INT AV_VERSION_INT(LIBAVFILTER_VERSION_MAJOR, \