LCOV - code coverage report
Current view: top level - src/libavfilter - af_aemphasis.c (source / functions) Hit Total Coverage
Test: coverage.info Lines: 142 193 73.6 %
Date: 2017-09-22 00:04:37 Functions: 8 8 100.0 %

          Line data    Source code
       1             : /*
       2             :  * Copyright (c) 2001-2010 Krzysztof Foltman, Markus Schmidt, Thor Harald Johansen, Damien Zammit and others
       3             :  *
       4             :  * This file is part of FFmpeg.
       5             :  *
       6             :  * FFmpeg is free software; you can redistribute it and/or
       7             :  * modify it under the terms of the GNU Lesser General Public
       8             :  * License as published by the Free Software Foundation; either
       9             :  * version 2.1 of the License, or (at your option) any later version.
      10             :  *
      11             :  * FFmpeg is distributed in the hope that it will be useful,
      12             :  * but WITHOUT ANY WARRANTY; without even the implied warranty of
      13             :  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
      14             :  * Lesser General Public License for more details.
      15             :  *
      16             :  * You should have received a copy of the GNU Lesser General Public
      17             :  * License along with FFmpeg; if not, write to the Free Software
      18             :  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
      19             :  */
      20             : 
      21             : #include "libavutil/opt.h"
      22             : #include "avfilter.h"
      23             : #include "internal.h"
      24             : #include "audio.h"
      25             : 
      26             : typedef struct BiquadCoeffs {
      27             :     double a0, a1, a2, b1, b2;
      28             : } BiquadCoeffs;
      29             : 
      30             : typedef struct BiquadD2 {
      31             :     double a0, a1, a2, b1, b2, w1, w2;
      32             : } BiquadD2;
      33             : 
      34             : typedef struct RIAACurve {
      35             :     BiquadD2 r1;
      36             :     BiquadD2 brickw;
      37             :     int use_brickw;
      38             : } RIAACurve;
      39             : 
      40             : typedef struct AudioEmphasisContext {
      41             :     const AVClass *class;
      42             :     int mode, type;
      43             :     double level_in, level_out;
      44             : 
      45             :     RIAACurve *rc;
      46             : } AudioEmphasisContext;
      47             : 
      48             : #define OFFSET(x) offsetof(AudioEmphasisContext, x)
      49             : #define FLAGS AV_OPT_FLAG_AUDIO_PARAM|AV_OPT_FLAG_FILTERING_PARAM
      50             : 
      51             : static const AVOption aemphasis_options[] = {
      52             :     { "level_in",      "set input gain", OFFSET(level_in),  AV_OPT_TYPE_DOUBLE, {.dbl=1}, 0, 64, FLAGS },
      53             :     { "level_out",    "set output gain", OFFSET(level_out), AV_OPT_TYPE_DOUBLE, {.dbl=1}, 0, 64, FLAGS },
      54             :     { "mode",         "set filter mode", OFFSET(mode), AV_OPT_TYPE_INT,   {.i64=0}, 0, 1, FLAGS, "mode" },
      55             :     { "reproduction",              NULL,            0, AV_OPT_TYPE_CONST, {.i64=0}, 0, 0, FLAGS, "mode" },
      56             :     { "production",                NULL,            0, AV_OPT_TYPE_CONST, {.i64=1}, 0, 0, FLAGS, "mode" },
      57             :     { "type",         "set filter type", OFFSET(type), AV_OPT_TYPE_INT,   {.i64=4}, 0, 8, FLAGS, "type" },
      58             :     { "col",                 "Columbia",            0, AV_OPT_TYPE_CONST, {.i64=0}, 0, 0, FLAGS, "type" },
      59             :     { "emi",                      "EMI",            0, AV_OPT_TYPE_CONST, {.i64=1}, 0, 0, FLAGS, "type" },
      60             :     { "bsi",              "BSI (78RPM)",            0, AV_OPT_TYPE_CONST, {.i64=2}, 0, 0, FLAGS, "type" },
      61             :     { "riaa",                    "RIAA",            0, AV_OPT_TYPE_CONST, {.i64=3}, 0, 0, FLAGS, "type" },
      62             :     { "cd",         "Compact Disc (CD)",            0, AV_OPT_TYPE_CONST, {.i64=4}, 0, 0, FLAGS, "type" },
      63             :     { "50fm",               "50µs (FM)",            0, AV_OPT_TYPE_CONST, {.i64=5}, 0, 0, FLAGS, "type" },
      64             :     { "75fm",               "75µs (FM)",            0, AV_OPT_TYPE_CONST, {.i64=6}, 0, 0, FLAGS, "type" },
      65             :     { "50kf",            "50µs (FM-KF)",            0, AV_OPT_TYPE_CONST, {.i64=7}, 0, 0, FLAGS, "type" },
      66             :     { "75kf",            "75µs (FM-KF)",            0, AV_OPT_TYPE_CONST, {.i64=8}, 0, 0, FLAGS, "type" },
      67             :     { NULL }
      68             : };
      69             : 
      70             : AVFILTER_DEFINE_CLASS(aemphasis);
      71             : 
      72     1587600 : static inline double biquad(BiquadD2 *bq, double in)
      73             : {
      74     1587600 :     double n = in;
      75     1587600 :     double tmp = n - bq->w1 * bq->b1 - bq->w2 * bq->b2;
      76     1587600 :     double out = tmp * bq->a0 + bq->w1 * bq->a1 + bq->w2 * bq->a2;
      77             : 
      78     1587600 :     bq->w2 = bq->w1;
      79     1587600 :     bq->w1 = tmp;
      80             : 
      81     1587600 :     return out;
      82             : }
      83             : 
      84         518 : static int filter_frame(AVFilterLink *inlink, AVFrame *in)
      85             : {
      86         518 :     AVFilterContext *ctx = inlink->dst;
      87         518 :     AVFilterLink *outlink = ctx->outputs[0];
      88         518 :     AudioEmphasisContext *s = ctx->priv;
      89         518 :     const double *src = (const double *)in->data[0];
      90         518 :     const double level_out = s->level_out;
      91         518 :     const double level_in = s->level_in;
      92             :     AVFrame *out;
      93             :     double *dst;
      94             :     int n, c;
      95             : 
      96         518 :     if (av_frame_is_writable(in)) {
      97         518 :         out = in;
      98             :     } else {
      99           0 :         out = ff_get_audio_buffer(inlink, in->nb_samples);
     100           0 :         if (!out) {
     101           0 :             av_frame_free(&in);
     102           0 :             return AVERROR(ENOMEM);
     103             :         }
     104           0 :         av_frame_copy_props(out, in);
     105             :     }
     106         518 :     dst = (double *)out->data[0];
     107             : 
     108      529718 :     for (n = 0; n < in->nb_samples; n++) {
     109     1587600 :         for (c = 0; c < inlink->channels; c++)
     110     1058400 :             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);
     111      529200 :         dst += inlink->channels;
     112      529200 :         src += inlink->channels;
     113             :     }
     114             : 
     115         518 :     if (in != out)
     116           0 :         av_frame_free(&in);
     117         518 :     return ff_filter_frame(outlink, out);
     118             : }
     119             : 
     120           2 : static int query_formats(AVFilterContext *ctx)
     121             : {
     122             :     AVFilterChannelLayouts *layouts;
     123             :     AVFilterFormats *formats;
     124             :     static const enum AVSampleFormat sample_fmts[] = {
     125             :         AV_SAMPLE_FMT_DBL,
     126             :         AV_SAMPLE_FMT_NONE
     127             :     };
     128             :     int ret;
     129             : 
     130           2 :     layouts = ff_all_channel_counts();
     131           2 :     if (!layouts)
     132           0 :         return AVERROR(ENOMEM);
     133           2 :     ret = ff_set_common_channel_layouts(ctx, layouts);
     134           2 :     if (ret < 0)
     135           0 :         return ret;
     136             : 
     137           2 :     formats = ff_make_format_list(sample_fmts);
     138           2 :     if (!formats)
     139           0 :         return AVERROR(ENOMEM);
     140           2 :     ret = ff_set_common_formats(ctx, formats);
     141           2 :     if (ret < 0)
     142           0 :         return ret;
     143             : 
     144           2 :     formats = ff_all_samplerates();
     145           2 :     if (!formats)
     146           0 :         return AVERROR(ENOMEM);
     147           2 :     return ff_set_common_samplerates(ctx, formats);
     148             : }
     149             : 
     150           1 : static inline void set_highshelf_rbj(BiquadD2 *bq, double freq, double q, double peak, double sr)
     151             : {
     152           1 :     double A = sqrt(peak);
     153           1 :     double w0 = freq * 2 * M_PI / sr;
     154           1 :     double alpha = sin(w0) / (2 * q);
     155           1 :     double cw0 = cos(w0);
     156           1 :     double tmp = 2 * sqrt(A) * alpha;
     157           1 :     double b0 = 0, ib0 = 0;
     158             : 
     159           1 :     bq->a0 =    A*( (A+1) + (A-1)*cw0 + tmp);
     160           1 :     bq->a1 = -2*A*( (A-1) + (A+1)*cw0);
     161           1 :     bq->a2 =    A*( (A+1) + (A-1)*cw0 - tmp);
     162           1 :         b0 =        (A+1) - (A-1)*cw0 + tmp;
     163           1 :     bq->b1 =    2*( (A-1) - (A+1)*cw0);
     164           1 :     bq->b2 =        (A+1) - (A-1)*cw0 - tmp;
     165             : 
     166           1 :     ib0     = 1 / b0;
     167           1 :     bq->b1 *= ib0;
     168           1 :     bq->b2 *= ib0;
     169           1 :     bq->a0 *= ib0;
     170           1 :     bq->a1 *= ib0;
     171           1 :     bq->a2 *= ib0;
     172           1 : }
     173             : 
     174           2 : static inline void set_lp_rbj(BiquadD2 *bq, double fc, double q, double sr, double gain)
     175             : {
     176           2 :     double omega = 2.0 * M_PI * fc / sr;
     177           2 :     double sn = sin(omega);
     178           2 :     double cs = cos(omega);
     179           2 :     double alpha = sn/(2 * q);
     180           2 :     double inv = 1.0/(1.0 + alpha);
     181             : 
     182           2 :     bq->a2 = bq->a0 = gain * inv * (1.0 - cs) * 0.5;
     183           2 :     bq->a1 = bq->a0 + bq->a0;
     184           2 :     bq->b1 = (-2.0 * cs * inv);
     185           2 :     bq->b2 = ((1.0 - alpha) * inv);
     186           2 : }
     187             : 
     188           1 : static double freq_gain(BiquadCoeffs *c, double freq, double sr)
     189             : {
     190             :     double zr, zi;
     191             : 
     192           1 :     freq *= 2.0 * M_PI / sr;
     193           1 :     zr = cos(freq);
     194           1 :     zi = -sin(freq);
     195             : 
     196             :     /* |(a0 + a1*z + a2*z^2)/(1 + b1*z + b2*z^2)| */
     197           2 :     return hypot(c->a0 + c->a1*zr + c->a2*(zr*zr-zi*zi), c->a1*zi + 2*c->a2*zr*zi) /
     198           1 :            hypot(1 + c->b1*zr + c->b2*(zr*zr-zi*zi), c->b1*zi + 2*c->b2*zr*zi);
     199             : }
     200             : 
     201           2 : static int config_input(AVFilterLink *inlink)
     202             : {
     203             :     double i, j, k, g, t, a0, a1, a2, b1, b2, tau1, tau2, tau3;
     204           2 :     double cutfreq, gain1kHz, gc, sr = inlink->sample_rate;
     205           2 :     AVFilterContext *ctx = inlink->dst;
     206           2 :     AudioEmphasisContext *s = ctx->priv;
     207             :     BiquadCoeffs coeffs;
     208             :     int ch;
     209             : 
     210           2 :     s->rc = av_calloc(inlink->channels, sizeof(*s->rc));
     211           2 :     if (!s->rc)
     212           0 :         return AVERROR(ENOMEM);
     213             : 
     214           2 :     switch (s->type) {
     215           0 :     case 0: //"Columbia"
     216           0 :         i = 100.;
     217           0 :         j = 500.;
     218           0 :         k = 1590.;
     219           0 :         break;
     220           0 :     case 1: //"EMI"
     221           0 :         i = 70.;
     222           0 :         j = 500.;
     223           0 :         k = 2500.;
     224           0 :         break;
     225           0 :     case 2: //"BSI(78rpm)"
     226           0 :         i = 50.;
     227           0 :         j = 353.;
     228           0 :         k = 3180.;
     229           0 :         break;
     230           1 :     case 3: //"RIAA"
     231             :     default:
     232           1 :         tau1 = 0.003180;
     233           1 :         tau2 = 0.000318;
     234           1 :         tau3 = 0.000075;
     235           1 :         i = 1. / (2. * M_PI * tau1);
     236           1 :         j = 1. / (2. * M_PI * tau2);
     237           1 :         k = 1. / (2. * M_PI * tau3);
     238           1 :         break;
     239           0 :     case 4: //"CD Mastering"
     240           0 :         tau1 = 0.000050;
     241           0 :         tau2 = 0.000015;
     242           0 :         tau3 = 0.0000001;// 1.6MHz out of audible range for null impact
     243           0 :         i = 1. / (2. * M_PI * tau1);
     244           0 :         j = 1. / (2. * M_PI * tau2);
     245           0 :         k = 1. / (2. * M_PI * tau3);
     246           0 :         break;
     247           1 :     case 5: //"50µs FM (Europe)"
     248           1 :         tau1 = 0.000050;
     249           1 :         tau2 = tau1 / 20;// not used
     250           1 :         tau3 = tau1 / 50;//
     251           1 :         i = 1. / (2. * M_PI * tau1);
     252           1 :         j = 1. / (2. * M_PI * tau2);
     253           1 :         k = 1. / (2. * M_PI * tau3);
     254           1 :         break;
     255           0 :     case 6: //"75µs FM (US)"
     256           0 :         tau1 = 0.000075;
     257           0 :         tau2 = tau1 / 20;// not used
     258           0 :         tau3 = tau1 / 50;//
     259           0 :         i = 1. / (2. * M_PI * tau1);
     260           0 :         j = 1. / (2. * M_PI * tau2);
     261           0 :         k = 1. / (2. * M_PI * tau3);
     262           0 :         break;
     263             :     }
     264             : 
     265           2 :     i *= 2 * M_PI;
     266           2 :     j *= 2 * M_PI;
     267           2 :     k *= 2 * M_PI;
     268             : 
     269           2 :     t = 1. / sr;
     270             : 
     271             :     //swap a1 b1, a2 b2
     272           3 :     if (s->type == 7 || s->type == 8) {
     273           1 :         double tau = (s->type == 7 ? 0.000050 : 0.000075);
     274           1 :         double f = 1.0 / (2 * M_PI * tau);
     275           1 :         double nyq = sr * 0.5;
     276           1 :         double gain = sqrt(1.0 + nyq * nyq / (f * f)); // gain at Nyquist
     277           1 :         double cfreq = sqrt((gain - 1.0) * f * f); // frequency
     278           1 :         double q = 1.0;
     279             : 
     280           1 :         if (s->type == 8)
     281           1 :             q = pow((sr / 3269.0) + 19.5, -0.25); // somewhat poor curve-fit
     282           1 :         if (s->type == 7)
     283           0 :             q = pow((sr / 4750.0) + 19.5, -0.25);
     284           1 :         if (s->mode == 0)
     285           1 :             set_highshelf_rbj(&s->rc[0].r1, cfreq, q, 1. / gain, sr);
     286             :         else
     287           0 :             set_highshelf_rbj(&s->rc[0].r1, cfreq, q, gain, sr);
     288           1 :         s->rc[0].use_brickw = 0;
     289             :     } else {
     290           1 :         s->rc[0].use_brickw = 1;
     291           1 :         if (s->mode == 0) { // Reproduction
     292           1 :             g  = 1. / (4.+2.*i*t+2.*k*t+i*k*t*t);
     293           1 :             a0 = (2.*t+j*t*t)*g;
     294           1 :             a1 = (2.*j*t*t)*g;
     295           1 :             a2 = (-2.*t+j*t*t)*g;
     296           1 :             b1 = (-8.+2.*i*k*t*t)*g;
     297           1 :             b2 = (4.-2.*i*t-2.*k*t+i*k*t*t)*g;
     298             :         } else {  // Production
     299           0 :             g  = 1. / (2.*t+j*t*t);
     300           0 :             a0 = (4.+2.*i*t+2.*k*t+i*k*t*t)*g;
     301           0 :             a1 = (-8.+2.*i*k*t*t)*g;
     302           0 :             a2 = (4.-2.*i*t-2.*k*t+i*k*t*t)*g;
     303           0 :             b1 = (2.*j*t*t)*g;
     304           0 :             b2 = (-2.*t+j*t*t)*g;
     305             :         }
     306             : 
     307           1 :         coeffs.a0 = a0;
     308           1 :         coeffs.a1 = a1;
     309           1 :         coeffs.a2 = a2;
     310           1 :         coeffs.b1 = b1;
     311           1 :         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           1 :         gain1kHz = freq_gain(&coeffs, 1000.0, sr);
     317             :         // divide one filter's x[n-m] coefficients by that value
     318           1 :         gc = 1.0 / gain1kHz;
     319           1 :         s->rc[0].r1.a0 = coeffs.a0 * gc;
     320           1 :         s->rc[0].r1.a1 = coeffs.a1 * gc;
     321           1 :         s->rc[0].r1.a2 = coeffs.a2 * gc;
     322           1 :         s->rc[0].r1.b1 = coeffs.b1;
     323           1 :         s->rc[0].r1.b2 = coeffs.b2;
     324             :     }
     325             : 
     326           2 :     cutfreq = FFMIN(0.45 * sr, 21000.);
     327           2 :     set_lp_rbj(&s->rc[0].brickw, cutfreq, 0.707, sr, 1.);
     328             : 
     329           4 :     for (ch = 1; ch < inlink->channels; ch++) {
     330           2 :         memcpy(&s->rc[ch], &s->rc[0], sizeof(RIAACurve));
     331             :     }
     332             : 
     333           2 :     return 0;
     334             : }
     335             : 
     336           2 : static av_cold void uninit(AVFilterContext *ctx)
     337             : {
     338           2 :     AudioEmphasisContext *s = ctx->priv;
     339           2 :     av_freep(&s->rc);
     340           2 : }
     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             : };

Generated by: LCOV version 1.13