avfilter/af_aiir: do not leak memory on failure in convert_zp2tf()
[ffmpeg.git] / libavfilter / af_aiir.c
1 /*
2  * Copyright (c) 2018 Paul B Mahol
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 <float.h>
22
23 #include "libavutil/avassert.h"
24 #include "libavutil/avstring.h"
25 #include "libavutil/opt.h"
26 #include "audio.h"
27 #include "avfilter.h"
28 #include "internal.h"
29
30 typedef struct ThreadData {
31     AVFrame *in, *out;
32 } ThreadData;
33
34 typedef struct Pair {
35     int a, b;
36 } Pair;
37
38 typedef struct BiquadContext {
39     double a0, a1, a2;
40     double b0, b1, b2;
41     double i1, i2;
42     double o1, o2;
43 } BiquadContext;
44
45 typedef struct IIRChannel {
46     int nb_ab[2];
47     double *ab[2];
48     double g;
49     double *cache[2];
50     BiquadContext *biquads;
51     int clippings;
52 } IIRChannel;
53
54 typedef struct AudioIIRContext {
55     const AVClass *class;
56     char *a_str, *b_str, *g_str;
57     double dry_gain, wet_gain;
58     int format;
59     int process;
60     int precision;
61
62     IIRChannel *iir;
63     int channels;
64     enum AVSampleFormat sample_format;
65
66     int (*iir_channel)(AVFilterContext *ctx, void *arg, int ch, int nb_jobs);
67 } AudioIIRContext;
68
69 static int query_formats(AVFilterContext *ctx)
70 {
71     AudioIIRContext *s = ctx->priv;
72     AVFilterFormats *formats;
73     AVFilterChannelLayouts *layouts;
74     enum AVSampleFormat sample_fmts[] = {
75         AV_SAMPLE_FMT_DBLP,
76         AV_SAMPLE_FMT_NONE
77     };
78     int ret;
79
80     layouts = ff_all_channel_counts();
81     if (!layouts)
82         return AVERROR(ENOMEM);
83     ret = ff_set_common_channel_layouts(ctx, layouts);
84     if (ret < 0)
85         return ret;
86
87     sample_fmts[0] = s->sample_format;
88     formats = ff_make_format_list(sample_fmts);
89     if (!formats)
90         return AVERROR(ENOMEM);
91     ret = ff_set_common_formats(ctx, formats);
92     if (ret < 0)
93         return ret;
94
95     formats = ff_all_samplerates();
96     if (!formats)
97         return AVERROR(ENOMEM);
98     return ff_set_common_samplerates(ctx, formats);
99 }
100
101 #define IIR_CH(name, type, min, max, need_clipping)                     \
102 static int iir_ch_## name(AVFilterContext *ctx, void *arg, int ch, int nb_jobs)  \
103 {                                                                       \
104     AudioIIRContext *s = ctx->priv;                                     \
105     const double ig = s->dry_gain;                                      \
106     const double og = s->wet_gain;                                      \
107     ThreadData *td = arg;                                               \
108     AVFrame *in = td->in, *out = td->out;                               \
109     const type *src = (const type *)in->extended_data[ch];              \
110     double *ic = (double *)s->iir[ch].cache[0];                         \
111     double *oc = (double *)s->iir[ch].cache[1];                         \
112     const int nb_a = s->iir[ch].nb_ab[0];                               \
113     const int nb_b = s->iir[ch].nb_ab[1];                               \
114     const double *a = s->iir[ch].ab[0];                                 \
115     const double *b = s->iir[ch].ab[1];                                 \
116     int *clippings = &s->iir[ch].clippings;                             \
117     type *dst = (type *)out->extended_data[ch];                         \
118     int n;                                                              \
119                                                                         \
120     for (n = 0; n < in->nb_samples; n++) {                              \
121         double sample = 0.;                                             \
122         int x;                                                          \
123                                                                         \
124         memmove(&ic[1], &ic[0], (nb_b - 1) * sizeof(*ic));              \
125         memmove(&oc[1], &oc[0], (nb_a - 1) * sizeof(*oc));              \
126         ic[0] = src[n] * ig;                                            \
127         for (x = 0; x < nb_b; x++)                                      \
128             sample += b[x] * ic[x];                                     \
129                                                                         \
130         for (x = 1; x < nb_a; x++)                                      \
131             sample -= a[x] * oc[x];                                     \
132                                                                         \
133         oc[0] = sample;                                                 \
134         sample *= og;                                                   \
135         if (need_clipping && sample < min) {                            \
136             *clippings++;                                               \
137             dst[n] = min;                                               \
138         } else if (need_clipping && sample > max) {                     \
139             *clippings++;                                               \
140             dst[n] = max;                                               \
141         } else {                                                        \
142             dst[n] = sample;                                            \
143         }                                                               \
144     }                                                                   \
145                                                                         \
146     return 0;                                                           \
147 }
148
149 IIR_CH(s16p, int16_t, INT16_MIN, INT16_MAX, 1)
150 IIR_CH(s32p, int32_t, INT32_MIN, INT32_MAX, 1)
151 IIR_CH(fltp, float,         -1.,        1., 0)
152 IIR_CH(dblp, double,        -1.,        1., 0)
153
154 #define SERIAL_IIR_CH(name, type, min, max, need_clipping)                  \
155 static int iir_ch_serial_## name(AVFilterContext *ctx, void *arg, int ch, int nb_jobs)  \
156 {                                                                       \
157     AudioIIRContext *s = ctx->priv;                                     \
158     const double ig = s->dry_gain;                                      \
159     const double og = s->wet_gain;                                      \
160     ThreadData *td = arg;                                               \
161     AVFrame *in = td->in, *out = td->out;                               \
162     const type *src = (const type *)in->extended_data[ch];              \
163     type *dst = (type *)out->extended_data[ch];                         \
164     IIRChannel *iir = &s->iir[ch];                                      \
165     int *clippings = &iir->clippings;                                   \
166     int nb_biquads = (FFMAX(iir->nb_ab[0], iir->nb_ab[1]) + 1) / 2;     \
167     int n, i;                                                           \
168                                                                         \
169     for (i = 0; i < nb_biquads; i++) {                                  \
170         const double a1 = -iir->biquads[i].a1;                          \
171         const double a2 = -iir->biquads[i].a2;                          \
172         const double b0 = iir->biquads[i].b0;                           \
173         const double b1 = iir->biquads[i].b1;                           \
174         const double b2 = iir->biquads[i].b2;                           \
175         double i1 = iir->biquads[i].i1;                                 \
176         double i2 = iir->biquads[i].i2;                                 \
177         double o1 = iir->biquads[i].o1;                                 \
178         double o2 = iir->biquads[i].o2;                                 \
179                                                                         \
180         for (n = 0; n < in->nb_samples; n++) {                          \
181             double sample = ig * (i ? dst[n] : src[n]);                 \
182             double o0 = sample * b0 + i1 * b1 + i2 * b2 + o1 * a1 + o2 * a2; \
183                                                                         \
184             i2 = i1;                                                    \
185             i1 = src[n];                                                \
186             o2 = o1;                                                    \
187             o1 = o0;                                                    \
188             o0 *= og;                                                   \
189                                                                         \
190             if (need_clipping && o0 < min) {                            \
191                 *clippings++;                                           \
192                 dst[n] = min;                                           \
193             } else if (need_clipping && o0 > max) {                     \
194                 *clippings++;                                           \
195                 dst[n] = max;                                           \
196             } else {                                                    \
197                 dst[n] = o0;                                            \
198             }                                                           \
199         }                                                               \
200         iir->biquads[i].i1 = i1;                                        \
201         iir->biquads[i].i2 = i2;                                        \
202         iir->biquads[i].o1 = o1;                                        \
203         iir->biquads[i].o2 = o2;                                        \
204     }                                                                   \
205                                                                         \
206     return 0;                                                           \
207 }
208
209 SERIAL_IIR_CH(s16p, int16_t, INT16_MIN, INT16_MAX, 1)
210 SERIAL_IIR_CH(s32p, int32_t, INT32_MIN, INT32_MAX, 1)
211 SERIAL_IIR_CH(fltp, float,         -1.,        1., 0)
212 SERIAL_IIR_CH(dblp, double,        -1.,        1., 0)
213
214 static void count_coefficients(char *item_str, int *nb_items)
215 {
216     char *p;
217
218     if (!item_str)
219         return;
220
221     *nb_items = 1;
222     for (p = item_str; *p && *p != '|'; p++) {
223         if (*p == ' ')
224             (*nb_items)++;
225     }
226 }
227
228 static int read_gains(AVFilterContext *ctx, char *item_str, int nb_items)
229 {
230     AudioIIRContext *s = ctx->priv;
231     char *p, *arg, *old_str, *prev_arg = NULL, *saveptr = NULL;
232     int i;
233
234     p = old_str = av_strdup(item_str);
235     if (!p)
236         return AVERROR(ENOMEM);
237     for (i = 0; i < nb_items; i++) {
238         if (!(arg = av_strtok(p, "|", &saveptr)))
239             arg = prev_arg;
240
241         if (!arg) {
242             av_freep(&old_str);
243             return AVERROR(EINVAL);
244         }
245
246         p = NULL;
247         if (sscanf(arg, "%lf", &s->iir[i].g) != 1) {
248             av_log(ctx, AV_LOG_ERROR, "Invalid gains supplied: %s\n", arg);
249             av_freep(&old_str);
250             return AVERROR(EINVAL);
251         }
252
253         prev_arg = arg;
254     }
255
256     av_freep(&old_str);
257
258     return 0;
259 }
260
261 static int read_tf_coefficients(AVFilterContext *ctx, char *item_str, int nb_items, double *dst)
262 {
263     char *p, *arg, *old_str, *saveptr = NULL;
264     int i;
265
266     p = old_str = av_strdup(item_str);
267     if (!p)
268         return AVERROR(ENOMEM);
269     for (i = 0; i < nb_items; i++) {
270         if (!(arg = av_strtok(p, " ", &saveptr)))
271             break;
272
273         p = NULL;
274         if (sscanf(arg, "%lf", &dst[i]) != 1) {
275             av_log(ctx, AV_LOG_ERROR, "Invalid coefficients supplied: %s\n", arg);
276             av_freep(&old_str);
277             return AVERROR(EINVAL);
278         }
279     }
280
281     av_freep(&old_str);
282
283     return 0;
284 }
285
286 static int read_zp_coefficients(AVFilterContext *ctx, char *item_str, int nb_items, double *dst)
287 {
288     char *p, *arg, *old_str, *saveptr = NULL;
289     int i;
290
291     p = old_str = av_strdup(item_str);
292     if (!p)
293         return AVERROR(ENOMEM);
294     for (i = 0; i < nb_items; i++) {
295         if (!(arg = av_strtok(p, " ", &saveptr)))
296             break;
297
298         p = NULL;
299         if (sscanf(arg, "%lf %lfi", &dst[i*2], &dst[i*2+1]) != 2) {
300             av_log(ctx, AV_LOG_ERROR, "Invalid coefficients supplied: %s\n", arg);
301             av_freep(&old_str);
302             return AVERROR(EINVAL);
303         }
304     }
305
306     av_freep(&old_str);
307
308     return 0;
309 }
310
311 static int read_channels(AVFilterContext *ctx, int channels, uint8_t *item_str, int ab)
312 {
313     AudioIIRContext *s = ctx->priv;
314     char *p, *arg, *old_str, *prev_arg = NULL, *saveptr = NULL;
315     int i, ret;
316
317     p = old_str = av_strdup(item_str);
318     if (!p)
319         return AVERROR(ENOMEM);
320     for (i = 0; i < channels; i++) {
321         IIRChannel *iir = &s->iir[i];
322
323         if (!(arg = av_strtok(p, "|", &saveptr)))
324             arg = prev_arg;
325
326         if (!arg) {
327             av_freep(&old_str);
328             return AVERROR(EINVAL);
329         }
330
331         count_coefficients(arg, &iir->nb_ab[ab]);
332
333         p = NULL;
334         iir->cache[ab] = av_calloc(iir->nb_ab[ab] + 1, sizeof(double));
335         iir->ab[ab] = av_calloc(iir->nb_ab[ab] * (s->format + 1), sizeof(double));
336         if (!iir->ab[ab] || !iir->cache[ab]) {
337             av_freep(&old_str);
338             return AVERROR(ENOMEM);
339         }
340
341         if (s->format) {
342             ret = read_zp_coefficients(ctx, arg, iir->nb_ab[ab], iir->ab[ab]);
343         } else {
344             ret = read_tf_coefficients(ctx, arg, iir->nb_ab[ab], iir->ab[ab]);
345         }
346         if (ret < 0) {
347             av_freep(&old_str);
348             return ret;
349         }
350         prev_arg = arg;
351     }
352
353     av_freep(&old_str);
354
355     return 0;
356 }
357
358 static void multiply(double wre, double wim, int npz, double *coeffs)
359 {
360     double nwre = -wre, nwim = -wim;
361     double cre, cim;
362     int i;
363
364     for (i = npz; i >= 1; i--) {
365         cre = coeffs[2 * i + 0];
366         cim = coeffs[2 * i + 1];
367
368         coeffs[2 * i + 0] = (nwre * cre - nwim * cim) + coeffs[2 * (i - 1) + 0];
369         coeffs[2 * i + 1] = (nwre * cim + nwim * cre) + coeffs[2 * (i - 1) + 1];
370     }
371
372     cre = coeffs[0];
373     cim = coeffs[1];
374     coeffs[0] = nwre * cre - nwim * cim;
375     coeffs[1] = nwre * cim + nwim * cre;
376 }
377
378 static int expand(AVFilterContext *ctx, double *pz, int nb, double *coeffs)
379 {
380     int i;
381
382     coeffs[0] = 1.0;
383     coeffs[1] = 0.0;
384
385     for (i = 0; i < nb; i++) {
386         coeffs[2 * (i + 1)    ] = 0.0;
387         coeffs[2 * (i + 1) + 1] = 0.0;
388     }
389
390     for (i = 0; i < nb; i++)
391         multiply(pz[2 * i], pz[2 * i + 1], nb, coeffs);
392
393     for (i = 0; i < nb + 1; i++) {
394         if (fabs(coeffs[2 * i + 1]) > FLT_EPSILON) {
395             av_log(ctx, AV_LOG_ERROR, "coeff: %lf of z^%d is not real; poles/zeros are not complex conjugates.\n",
396                    coeffs[2 * i + 1], i);
397             return AVERROR(EINVAL);
398         }
399     }
400
401     return 0;
402 }
403
404 static int convert_zp2tf(AVFilterContext *ctx, int channels)
405 {
406     AudioIIRContext *s = ctx->priv;
407     int ch, i, j, ret = 0;
408
409     for (ch = 0; ch < channels; ch++) {
410         IIRChannel *iir = &s->iir[ch];
411         double *topc, *botc;
412
413         topc = av_calloc((iir->nb_ab[0] + 1) * 2, sizeof(*topc));
414         botc = av_calloc((iir->nb_ab[1] + 1) * 2, sizeof(*botc));
415         if (!topc || !botc) {
416             ret = AVERROR(ENOMEM);
417             goto fail;
418         }
419
420         ret = expand(ctx, iir->ab[0], iir->nb_ab[0], botc);
421         if (ret < 0) {
422             goto fail;
423         }
424
425         ret = expand(ctx, iir->ab[1], iir->nb_ab[1], topc);
426         if (ret < 0) {
427             goto fail;
428         }
429
430         for (j = 0, i = iir->nb_ab[1]; i >= 0; j++, i--) {
431             iir->ab[1][j] = topc[2 * i];
432         }
433         iir->nb_ab[1]++;
434
435         for (j = 0, i = iir->nb_ab[0]; i >= 0; j++, i--) {
436             iir->ab[0][j] = botc[2 * i];
437         }
438         iir->nb_ab[0]++;
439
440 fail:
441         av_free(topc);
442         av_free(botc);
443         if (ret < 0)
444             break;
445     }
446
447     return ret;
448 }
449
450 static int decompose_zp2biquads(AVFilterContext *ctx, int channels)
451 {
452     AudioIIRContext *s = ctx->priv;
453     int ch, ret;
454
455     for (ch = 0; ch < channels; ch++) {
456         IIRChannel *iir = &s->iir[ch];
457         int nb_biquads = (FFMAX(iir->nb_ab[0], iir->nb_ab[1]) + 1) / 2;
458         int current_biquad = 0;
459
460         iir->biquads = av_calloc(nb_biquads, sizeof(BiquadContext));
461         if (!iir->biquads)
462             return AVERROR(ENOMEM);
463
464         while (nb_biquads--) {
465             Pair outmost_pole = { -1, -1 };
466             Pair nearest_zero = { -1, -1 };
467             double zeros[4] = { 0 };
468             double poles[4] = { 0 };
469             double b[6] = { 0 };
470             double a[6] = { 0 };
471             double min_distance = DBL_MAX;
472             double max_mag = 0;
473             int i;
474
475             for (i = 0; i < iir->nb_ab[0]; i++) {
476                 double mag;
477
478                 if (isnan(iir->ab[0][2 * i]) || isnan(iir->ab[0][2 * i + 1]))
479                     continue;
480                 mag = hypot(iir->ab[0][2 * i], iir->ab[0][2 * i + 1]);
481
482                 if (mag > max_mag) {
483                     max_mag = mag;
484                     outmost_pole.a = i;
485                 }
486             }
487
488             for (i = 0; i < iir->nb_ab[1]; i++) {
489                 if (isnan(iir->ab[0][2 * i]) || isnan(iir->ab[0][2 * i + 1]))
490                     continue;
491
492                 if (iir->ab[0][2 * i    ] ==  iir->ab[0][2 * outmost_pole.a    ] &&
493                     iir->ab[0][2 * i + 1] == -iir->ab[0][2 * outmost_pole.a + 1]) {
494                     outmost_pole.b = i;
495                     break;
496                 }
497             }
498
499             av_log(ctx, AV_LOG_VERBOSE, "outmost_pole is %d.%d\n", outmost_pole.a, outmost_pole.b);
500
501             if (outmost_pole.a < 0 || outmost_pole.b < 0)
502                 return AVERROR(EINVAL);
503
504             for (i = 0; i < iir->nb_ab[1]; i++) {
505                 double distance;
506
507                 if (isnan(iir->ab[1][2 * i]) || isnan(iir->ab[1][2 * i + 1]))
508                     continue;
509                 distance = hypot(iir->ab[0][2 * outmost_pole.a    ] - iir->ab[1][2 * i    ],
510                                  iir->ab[0][2 * outmost_pole.a + 1] - iir->ab[1][2 * i + 1]);
511
512                 if (distance < min_distance) {
513                     min_distance = distance;
514                     nearest_zero.a = i;
515                 }
516             }
517
518             for (i = 0; i < iir->nb_ab[1]; i++) {
519                 if (isnan(iir->ab[1][2 * i]) || isnan(iir->ab[1][2 * i + 1]))
520                     continue;
521
522                 if (iir->ab[1][2 * i    ] ==  iir->ab[1][2 * nearest_zero.a    ] &&
523                     iir->ab[1][2 * i + 1] == -iir->ab[1][2 * nearest_zero.a + 1]) {
524                     nearest_zero.b = i;
525                     break;
526                 }
527             }
528
529             av_log(ctx, AV_LOG_VERBOSE, "nearest_zero is %d.%d\n", nearest_zero.a, nearest_zero.b);
530
531             if (nearest_zero.a < 0 || nearest_zero.b < 0)
532                 return AVERROR(EINVAL);
533
534             poles[0] = iir->ab[0][2 * outmost_pole.a    ];
535             poles[1] = iir->ab[0][2 * outmost_pole.a + 1];
536
537             zeros[0] = iir->ab[1][2 * nearest_zero.a    ];
538             zeros[1] = iir->ab[1][2 * nearest_zero.a + 1];
539
540             if (nearest_zero.a == nearest_zero.b && outmost_pole.a == outmost_pole.b) {
541                 zeros[2] = 0;
542                 zeros[3] = 0;
543
544                 poles[2] = 0;
545                 poles[3] = 0;
546             } else {
547                 poles[2] = iir->ab[0][2 * outmost_pole.b    ];
548                 poles[3] = iir->ab[0][2 * outmost_pole.b + 1];
549
550                 zeros[2] = iir->ab[1][2 * nearest_zero.b    ];
551                 zeros[3] = iir->ab[1][2 * nearest_zero.b + 1];
552             }
553
554             ret = expand(ctx, zeros, 2, b);
555             if (ret < 0)
556                 return ret;
557
558             ret = expand(ctx, poles, 2, a);
559             if (ret < 0)
560                 return ret;
561
562             iir->ab[0][2 * outmost_pole.a] = iir->ab[0][2 * outmost_pole.a + 1] = NAN;
563             iir->ab[0][2 * outmost_pole.b] = iir->ab[0][2 * outmost_pole.b + 1] = NAN;
564             iir->ab[1][2 * nearest_zero.a] = iir->ab[1][2 * nearest_zero.a + 1] = NAN;
565             iir->ab[1][2 * nearest_zero.b] = iir->ab[1][2 * nearest_zero.b + 1] = NAN;
566
567             iir->biquads[current_biquad].a0 = 1.0;
568             iir->biquads[current_biquad].a1 = a[2] / a[4];
569             iir->biquads[current_biquad].a2 = a[0] / a[4];
570             iir->biquads[current_biquad].b0 = b[4] / a[4] * (current_biquad ? 1.0 : iir->g);
571             iir->biquads[current_biquad].b1 = b[2] / a[4] * (current_biquad ? 1.0 : iir->g);
572             iir->biquads[current_biquad].b2 = b[0] / a[4] * (current_biquad ? 1.0 : iir->g);
573
574             av_log(ctx, AV_LOG_VERBOSE, "a=%lf %lf %lf:b=%lf %lf %lf\n",
575                    iir->biquads[current_biquad].a0,
576                    iir->biquads[current_biquad].a1,
577                    iir->biquads[current_biquad].a2,
578                    iir->biquads[current_biquad].b0,
579                    iir->biquads[current_biquad].b1,
580                    iir->biquads[current_biquad].b2);
581
582             current_biquad++;
583         }
584     }
585
586     return 0;
587 }
588
589 static int config_output(AVFilterLink *outlink)
590 {
591     AVFilterContext *ctx = outlink->src;
592     AudioIIRContext *s = ctx->priv;
593     AVFilterLink *inlink = ctx->inputs[0];
594     int ch, ret, i;
595
596     s->channels = inlink->channels;
597     s->iir = av_calloc(s->channels, sizeof(*s->iir));
598     if (!s->iir)
599         return AVERROR(ENOMEM);
600
601     ret = read_gains(ctx, s->g_str, inlink->channels);
602     if (ret < 0)
603         return ret;
604
605     ret = read_channels(ctx, inlink->channels, s->a_str, 0);
606     if (ret < 0)
607         return ret;
608
609     ret = read_channels(ctx, inlink->channels, s->b_str, 1);
610     if (ret < 0)
611         return ret;
612
613     if (s->format == 0)
614         av_log(ctx, AV_LOG_WARNING, "tf coefficients format is not recommended for too high number of zeros/poles.\n");
615
616     if (s->format == 1 && s->process == 0) {
617         av_log(ctx, AV_LOG_WARNING, "Direct processsing is not recommended for zp coefficients format.\n");
618
619         ret = convert_zp2tf(ctx, inlink->channels);
620         if (ret < 0)
621             return ret;
622     } else if (s->format == 0 && s->process == 1) {
623         av_log(ctx, AV_LOG_ERROR, "Serial cascading is not implemented for transfer function.\n");
624         return AVERROR_PATCHWELCOME;
625     } else if (s->format == 1 && s->process == 1) {
626         if (inlink->format == AV_SAMPLE_FMT_S16P)
627             av_log(ctx, AV_LOG_WARNING, "Serial cascading is not recommended for i16 precision.\n");
628
629         ret = decompose_zp2biquads(ctx, inlink->channels);
630         if (ret < 0)
631             return ret;
632     }
633
634     for (ch = 0; ch < inlink->channels; ch++) {
635         IIRChannel *iir = &s->iir[ch];
636
637         for (i = 1; i < iir->nb_ab[0]; i++) {
638             iir->ab[0][i] /= iir->ab[0][0];
639         }
640
641         for (i = 0; i < iir->nb_ab[1]; i++) {
642             iir->ab[1][i] *= iir->g / iir->ab[0][0];
643         }
644     }
645
646     switch (inlink->format) {
647     case AV_SAMPLE_FMT_DBLP: s->iir_channel = s->process == 1 ? iir_ch_serial_dblp : iir_ch_dblp; break;
648     case AV_SAMPLE_FMT_FLTP: s->iir_channel = s->process == 1 ? iir_ch_serial_fltp : iir_ch_fltp; break;
649     case AV_SAMPLE_FMT_S32P: s->iir_channel = s->process == 1 ? iir_ch_serial_s32p : iir_ch_s32p; break;
650     case AV_SAMPLE_FMT_S16P: s->iir_channel = s->process == 1 ? iir_ch_serial_s16p : iir_ch_s16p; break;
651     }
652
653     return 0;
654 }
655
656 static int filter_frame(AVFilterLink *inlink, AVFrame *in)
657 {
658     AVFilterContext *ctx = inlink->dst;
659     AudioIIRContext *s = ctx->priv;
660     AVFilterLink *outlink = ctx->outputs[0];
661     ThreadData td;
662     AVFrame *out;
663     int ch;
664
665     if (av_frame_is_writable(in)) {
666         out = in;
667     } else {
668         out = ff_get_audio_buffer(outlink, in->nb_samples);
669         if (!out) {
670             av_frame_free(&in);
671             return AVERROR(ENOMEM);
672         }
673         av_frame_copy_props(out, in);
674     }
675
676     td.in  = in;
677     td.out = out;
678     ctx->internal->execute(ctx, s->iir_channel, &td, NULL, outlink->channels);
679
680     for (ch = 0; ch < outlink->channels; ch++) {
681         if (s->iir[ch].clippings > 0)
682             av_log(ctx, AV_LOG_WARNING, "Channel %d clipping %d times. Please reduce gain.\n",
683                    ch, s->iir[ch].clippings);
684         s->iir[ch].clippings = 0;
685     }
686
687     if (in != out)
688         av_frame_free(&in);
689
690     return ff_filter_frame(outlink, out);
691 }
692
693 static av_cold int init(AVFilterContext *ctx)
694 {
695     AudioIIRContext *s = ctx->priv;
696
697     if (!s->a_str || !s->b_str || !s->g_str) {
698         av_log(ctx, AV_LOG_ERROR, "Valid coefficients are mandatory.\n");
699         return AVERROR(EINVAL);
700     }
701
702     switch (s->precision) {
703     case 0: s->sample_format = AV_SAMPLE_FMT_DBLP; break;
704     case 1: s->sample_format = AV_SAMPLE_FMT_FLTP; break;
705     case 2: s->sample_format = AV_SAMPLE_FMT_S32P; break;
706     case 3: s->sample_format = AV_SAMPLE_FMT_S16P; break;
707     default: return AVERROR_BUG;
708     }
709
710     return 0;
711 }
712
713 static av_cold void uninit(AVFilterContext *ctx)
714 {
715     AudioIIRContext *s = ctx->priv;
716     int ch;
717
718     if (s->iir) {
719         for (ch = 0; ch < s->channels; ch++) {
720             IIRChannel *iir = &s->iir[ch];
721             av_freep(&iir->ab[0]);
722             av_freep(&iir->ab[1]);
723             av_freep(&iir->cache[0]);
724             av_freep(&iir->cache[1]);
725             av_freep(&iir->biquads);
726         }
727     }
728     av_freep(&s->iir);
729 }
730
731 static const AVFilterPad inputs[] = {
732     {
733         .name         = "default",
734         .type         = AVMEDIA_TYPE_AUDIO,
735         .filter_frame = filter_frame,
736     },
737     { NULL }
738 };
739
740 static const AVFilterPad outputs[] = {
741     {
742         .name         = "default",
743         .type         = AVMEDIA_TYPE_AUDIO,
744         .config_props = config_output,
745     },
746     { NULL }
747 };
748
749 #define OFFSET(x) offsetof(AudioIIRContext, x)
750 #define AF AV_OPT_FLAG_AUDIO_PARAM|AV_OPT_FLAG_FILTERING_PARAM
751
752 static const AVOption aiir_options[] = {
753     { "z", "set B/numerator/zeros coefficients",   OFFSET(b_str),    AV_OPT_TYPE_STRING, {.str="1+0i 1-0i"}, 0, 0, AF },
754     { "p", "set A/denominator/poles coefficients", OFFSET(a_str),    AV_OPT_TYPE_STRING, {.str="1+0i 1-0i"}, 0, 0, AF },
755     { "k", "set channels gains",                   OFFSET(g_str),    AV_OPT_TYPE_STRING, {.str="1|1"}, 0, 0, AF },
756     { "dry", "set dry gain",                       OFFSET(dry_gain), AV_OPT_TYPE_DOUBLE, {.dbl=1},     0, 1, AF },
757     { "wet", "set wet gain",                       OFFSET(wet_gain), AV_OPT_TYPE_DOUBLE, {.dbl=1},     0, 1, AF },
758     { "f", "set coefficients format",              OFFSET(format),   AV_OPT_TYPE_INT,    {.i64=1},     0, 1, AF, "format" },
759     { "tf", "transfer function",                   0,                AV_OPT_TYPE_CONST,  {.i64=0},     0, 0, AF, "format" },
760     { "zp", "Z-plane zeros/poles",                 0,                AV_OPT_TYPE_CONST,  {.i64=1},     0, 0, AF, "format" },
761     { "r", "set kind of processing",               OFFSET(process),  AV_OPT_TYPE_INT,    {.i64=1},     0, 1, AF, "process" },
762     { "d", "direct",                               0,                AV_OPT_TYPE_CONST,  {.i64=0},     0, 0, AF, "process" },
763     { "s", "serial cascading",                     0,                AV_OPT_TYPE_CONST,  {.i64=1},     0, 0, AF, "process" },
764     { "e", "set precision",                        OFFSET(precision),AV_OPT_TYPE_INT,    {.i64=0},     0, 3, AF, "precision" },
765     { "dbl", "double-precision floating-point",    0,                AV_OPT_TYPE_CONST,  {.i64=0},     0, 0, AF, "precision" },
766     { "flt", "single-precision floating-point",    0,                AV_OPT_TYPE_CONST,  {.i64=1},     0, 0, AF, "precision" },
767     { "i32", "32-bit integers",                    0,                AV_OPT_TYPE_CONST,  {.i64=2},     0, 0, AF, "precision" },
768     { "i16", "16-bit integers",                    0,                AV_OPT_TYPE_CONST,  {.i64=3},     0, 0, AF, "precision" },
769     { NULL },
770 };
771
772 AVFILTER_DEFINE_CLASS(aiir);
773
774 AVFilter ff_af_aiir = {
775     .name          = "aiir",
776     .description   = NULL_IF_CONFIG_SMALL("Apply Infinite Impulse Response filter with supplied coefficients."),
777     .priv_size     = sizeof(AudioIIRContext),
778     .priv_class    = &aiir_class,
779     .init          = init,
780     .uninit        = uninit,
781     .query_formats = query_formats,
782     .inputs        = inputs,
783     .outputs       = outputs,
784     .flags         = AVFILTER_FLAG_SLICE_THREADS,
785 };