31 #define RDFT_BITS_MIN 4
32 #define RDFT_BITS_MAX 16
56 #define NB_GAIN_ENTRY_MAX 4096
126 #define OFFSET(x) offsetof(FIREqualizerContext, x)
127 #define FLAGS AV_OPT_FLAG_AUDIO_PARAM|AV_OPT_FLAG_FILTERING_PARAM
128 #define TFLAGS AV_OPT_FLAG_AUDIO_PARAM|AV_OPT_FLAG_FILTERING_PARAM|AV_OPT_FLAG_RUNTIME_PARAM
172 s->analysis_rdft =
s->analysis_irdft =
s->rdft =
s->irdft =
NULL;
174 s->cepstrum_rdft =
NULL;
175 s->cepstrum_irdft =
NULL;
202 if (nsamples <= s->nsamples_max) {
203 float *buf = conv_buf + idx->buf_idx *
s->rdft_len;
204 float *obuf = conv_buf + !idx->buf_idx *
s->rdft_len + idx->overlap_idx;
205 float *tbuf =
s->tx_buf;
206 int center =
s->fir_len/2;
209 memset(buf, 0, center *
sizeof(*
data));
210 memcpy(buf + center,
data, nsamples *
sizeof(*
data));
211 memset(buf + center + nsamples, 0, (
s->rdft_len - nsamples - center) *
sizeof(*
data));
212 s->rdft_fn(
s->rdft, tbuf, buf,
sizeof(
float));
214 for (k = 0; k <=
s->rdft_len/2; k++) {
215 tbuf[2*k] *= kernel_buf[k];
216 tbuf[2*k+1] *= kernel_buf[k];
220 for (k = 0; k <
s->rdft_len - idx->overlap_idx; k++)
222 memcpy(
data, buf, nsamples *
sizeof(*
data));
223 idx->buf_idx = !idx->buf_idx;
224 idx->overlap_idx = nsamples;
226 while (nsamples >
s->nsamples_max * 2) {
228 data +=
s->nsamples_max;
229 nsamples -=
s->nsamples_max;
238 float *restrict
data,
int nsamples)
240 if (nsamples <= s->nsamples_max) {
241 float *buf = conv_buf + idx->buf_idx *
s->rdft_len;
242 float *obuf = conv_buf + !idx->buf_idx *
s->rdft_len + idx->overlap_idx;
243 float *tbuf =
s->tx_buf;
246 memcpy(buf,
data, nsamples *
sizeof(*
data));
247 memset(buf + nsamples, 0, (
s->rdft_len - nsamples) *
sizeof(*
data));
248 s->rdft_fn(
s->rdft, tbuf, buf,
sizeof(
float));
250 for (k = 0; k <
s->rdft_len + 2; k += 2) {
252 re = tbuf[k] * kernel_buf[k] - tbuf[k+1] * kernel_buf[k+1];
253 im = tbuf[k] * kernel_buf[k+1] + tbuf[k+1] * kernel_buf[k];
259 for (k = 0; k <
s->rdft_len - idx->overlap_idx; k++)
261 memcpy(
data, buf, nsamples *
sizeof(*
data));
262 idx->buf_idx = !idx->buf_idx;
263 idx->overlap_idx = nsamples;
265 while (nsamples >
s->nsamples_max * 2) {
267 data +=
s->nsamples_max;
268 nsamples -=
s->nsamples_max;
276 OverlapIndex *restrict idx,
float *restrict data0,
float *restrict data1,
int nsamples)
278 if (nsamples <= s->nsamples_max) {
280 AVComplexFloat *obuf = conv_buf + !idx->buf_idx *
s->rdft_len + idx->overlap_idx;
282 int center =
s->fir_len/2;
286 memset(buf, 0, center *
sizeof(*buf));
287 for (k = 0; k < nsamples; k++) {
288 buf[center+k].
re = data0[k];
289 buf[center+k].
im = data1[k];
291 memset(buf + center + nsamples, 0, (
s->rdft_len - nsamples - center) *
sizeof(*buf));
297 tbuf[0].
re = 0.5f * kernel_buf[0] * tbuf[0].
im;
298 tbuf[0].
im = 0.5f * kernel_buf[0] *
tmp;
299 for (k = 1; k <
s->rdft_len/2; k++) {
300 int m =
s->rdft_len - k;
302 tbuf[k].
re = 0.5f * kernel_buf[k] * tbuf[k].
im;
303 tbuf[k].
im = 0.5f * kernel_buf[k] *
tmp;
305 tbuf[m].
re = 0.5f * kernel_buf[k] * tbuf[m].
im;
306 tbuf[m].
im = 0.5f * kernel_buf[k] *
tmp;
309 tbuf[k].
re = 0.5f * kernel_buf[k] * tbuf[k].
im;
310 tbuf[k].
im = 0.5f * kernel_buf[k] *
tmp;
314 for (k = 0; k <
s->rdft_len - idx->overlap_idx; k++) {
315 buf[k].
re += obuf[k].
re;
316 buf[k].
im += obuf[k].
im;
320 for (k = 0; k < nsamples; k++) {
321 data0[k] = buf[k].
im;
322 data1[k] = buf[k].
re;
324 idx->buf_idx = !idx->buf_idx;
325 idx->overlap_idx = nsamples;
327 while (nsamples >
s->nsamples_max * 2) {
329 data0 +=
s->nsamples_max;
330 data1 +=
s->nsamples_max;
331 nsamples -=
s->nsamples_max;
334 fast_convolute2(
s, kernel_buf, conv_buf, idx, data0 + nsamples/2, data1 + nsamples/2, nsamples - nsamples/2);
341 int rate =
ctx->inputs[0]->sample_rate;
345 int center =
s->fir_len / 2;
346 double delay =
s->zero_phase ? 0.0 : (
double) center / rate;
350 s->analysis_buf[0] *=
s->rdft_len/2;
351 for (x = 1; x <= center; x++) {
352 s->analysis_buf[x] *=
s->rdft_len/2;
353 s->analysis_buf[
s->analysis_rdft_len - x] *=
s->rdft_len/2;
356 for (x = 0; x <
s->fir_len; x++)
357 s->analysis_buf[x] *=
s->rdft_len/2;
363 fprintf(
fp,
"# time[%d] (time amplitude)\n", ch);
366 for (x = center; x > 0; x--)
367 fprintf(
fp,
"%15.10f %15.10f\n", delay - (
double) x / rate, (
double)
s->analysis_buf[
s->analysis_rdft_len - x]);
369 for (x = 0; x <= center; x++)
370 fprintf(
fp,
"%15.10f %15.10f\n", delay + (
double)x / rate , (
double)
s->analysis_buf[x]);
372 for (x = 0; x <
s->fir_len; x++)
373 fprintf(
fp,
"%15.10f %15.10f\n", (
double)x / rate, (
double)
s->analysis_buf[x]);
376 s->analysis_rdft_fn(
s->analysis_rdft,
s->analysis_tbuf,
s->analysis_buf,
sizeof(
float));
378 fprintf(
fp,
"\n\n# freq[%d] (frequency desired_gain actual_gain)\n", ch);
380 for (x = 0; x <=
s->analysis_rdft_len/2; x++) {
382 vx = (
double)x * rate /
s->analysis_rdft_len;
386 yb =
s->min_phase ? hypotf(
s->analysis_tbuf[
i],
s->analysis_tbuf[
i+1]) :
s->analysis_tbuf[
i];
390 ya = 20.0 * log10(
fabs(ya));
391 yb = 20.0 * log10(
fabs(yb));
393 fprintf(
fp,
"%17.10f %17.10f %17.10f\n", vx, ya, yb);
404 s->gain_entry_err =
AVERROR(EINVAL);
410 s->gain_entry_err =
AVERROR(EINVAL);
414 if (
s->nb_gain_entry > 0 && freq <= s->gain_entry_tbl[
s->nb_gain_entry - 1].freq) {
416 s->gain_entry_err =
AVERROR(EINVAL);
420 s->gain_entry_tbl[
s->nb_gain_entry].freq = freq;
421 s->gain_entry_tbl[
s->nb_gain_entry].gain = gain;
428 const double *freq =
key;
431 if (*freq <
entry[0].freq)
433 if (*freq >
entry[1].freq)
448 if (!
s->nb_gain_entry)
451 if (freq <= s->gain_entry_tbl[0].freq)
452 return s->gain_entry_tbl[0].gain;
454 if (freq >=
s->gain_entry_tbl[
s->nb_gain_entry-1].freq)
455 return s->gain_entry_tbl[
s->nb_gain_entry-1].gain;
457 res = bsearch(&freq, &
s->gain_entry_tbl,
s->nb_gain_entry - 1,
sizeof(*res),
gain_entry_compare);
461 d0 = freq - res[0].
freq;
462 d1 = res[1].
freq - freq;
465 return (d0 * res[1].gain + d1 * res[0].gain) /
d;
480 double m0, m1, m2, msum, unit;
482 if (!
s->nb_gain_entry)
485 if (freq <= s->gain_entry_tbl[0].freq)
486 return s->gain_entry_tbl[0].gain;
488 if (freq >=
s->gain_entry_tbl[
s->nb_gain_entry-1].freq)
489 return s->gain_entry_tbl[
s->nb_gain_entry-1].gain;
491 res = bsearch(&freq, &
s->gain_entry_tbl,
s->nb_gain_entry - 1,
sizeof(*res),
gain_entry_compare);
495 m0 = res !=
s->gain_entry_tbl ?
496 unit * (res[0].
gain - res[-1].
gain) / (res[0].freq - res[-1].freq) : 0;
498 m2 = res !=
s->gain_entry_tbl +
s->nb_gain_entry - 2 ?
499 unit * (res[2].
gain - res[1].
gain) / (res[2].freq - res[1].freq) : 0;
502 m0 = msum > 0 ? (
fabs(m0) * m1 +
fabs(m1) * m0) / msum : 0;
504 m1 = msum > 0 ? (
fabs(m1) * m2 +
fabs(m2) * m1) / msum : 0;
508 b = 3 * res[1].
gain - m1 - 2 *
c - 3 *
d;
511 x = (freq - res[0].
freq) / unit;
515 return a * x3 +
b * x2 +
c * x +
d;
540 int k, cepstrum_len =
s->cepstrum_len, rdft_len =
s->rdft_len;
541 double norm = 2.0 / cepstrum_len;
542 double minval = 1e-7 / rdft_len;
544 memset(
s->cepstrum_buf, 0, cepstrum_len *
sizeof(*
s->cepstrum_buf));
545 memset(
s->cepstrum_tbuf, 0, (cepstrum_len + 2) *
sizeof(*
s->cepstrum_tbuf));
546 memcpy(
s->cepstrum_buf, rdft_buf, rdft_len/2 *
sizeof(*rdft_buf));
547 memcpy(
s->cepstrum_buf + cepstrum_len - rdft_len/2, rdft_buf + rdft_len/2, rdft_len/2 *
sizeof(*rdft_buf));
549 s->cepstrum_rdft_fn(
s->cepstrum_rdft,
s->cepstrum_tbuf,
s->cepstrum_buf,
sizeof(
float));
551 for (k = 0; k < cepstrum_len + 2; k += 2) {
552 s->cepstrum_tbuf[k] =
log(
FFMAX(
s->cepstrum_tbuf[k], minval));
553 s->cepstrum_tbuf[k+1] = 0;
556 s->cepstrum_irdft_fn(
s->cepstrum_irdft,
s->cepstrum_buf,
s->cepstrum_tbuf,
sizeof(
AVComplexFloat));
558 memset(
s->cepstrum_buf + cepstrum_len/2 + 1, 0, (cepstrum_len/2 - 1) *
sizeof(*
s->cepstrum_buf));
559 for (k = 1; k <= cepstrum_len/2; k++)
560 s->cepstrum_buf[k] *= 2;
562 s->cepstrum_rdft_fn(
s->cepstrum_rdft,
s->cepstrum_tbuf,
s->cepstrum_buf,
sizeof(
float));
564 for (k = 0; k < cepstrum_len + 2; k += 2) {
565 double mag =
exp(
s->cepstrum_tbuf[k] * norm) * norm;
566 double ph =
s->cepstrum_tbuf[k+1] * norm;
567 s->cepstrum_tbuf[k] = mag * cos(
ph);
568 s->cepstrum_tbuf[k+1] = mag * sin(
ph);
571 s->cepstrum_irdft_fn(
s->cepstrum_irdft,
s->cepstrum_buf,
s->cepstrum_tbuf,
sizeof(
AVComplexFloat));
572 memset(rdft_buf, 0,
s->rdft_len *
sizeof(*rdft_buf));
573 memcpy(rdft_buf,
s->cepstrum_buf,
s->fir_len *
sizeof(*rdft_buf));
576 memset(
s->analysis_buf, 0, (
s->analysis_rdft_len + 2) *
sizeof(*
s->analysis_buf));
577 memcpy(
s->analysis_buf,
s->cepstrum_buf,
s->fir_len *
sizeof(*
s->analysis_buf));
585 const char *gain_entry_func_names[] = {
"entry",
NULL };
586 const char *gain_func_names[] = {
"gain_interpolate",
"cubic_interpolate",
NULL };
591 int ret, k, center, ch;
594 FILE *dump_fp =
NULL;
596 s->nb_gain_entry = 0;
597 s->gain_entry_err = 0;
601 gain_entry_func_names, gain_entry_funcs,
ctx, 0,
ctx);
604 if (
s->gain_entry_err < 0)
605 return s->gain_entry_err;
615 if (
s->dumpfile && (!
s->dump_buf || !
s->analysis_rdft || !(dump_fp =
avpriv_fopen_utf8(
s->dumpfile,
"w"))))
620 inlink->ch_layout.u.mask : 0;
622 for (ch = 0; ch <
inlink->ch_layout.nb_channels; ch++) {
623 float *rdft_buf =
s->kernel_tmp_buf + ch * (
s->rdft_len * 2);
624 float *rdft_tbuf =
s->kernel_tmp_tbuf;
629 for (k = 0; k <=
s->analysis_rdft_len/2; k++) {
635 s->analysis_tbuf[2*k+1] = 0.0;
639 memcpy(
s->dump_buf,
s->analysis_tbuf, (
s->analysis_rdft_len + 2) *
sizeof(*
s->analysis_tbuf));
641 s->analysis_irdft_fn(
s->analysis_irdft,
s->analysis_buf,
s->analysis_tbuf,
sizeof(
AVComplexFloat));
642 center =
s->fir_len / 2;
644 for (k = 0; k <= center; k++) {
645 double u = k * (
M_PI/center);
652 win = 0.5 + 0.5 * cos(
u);
655 win = 0.53836 + 0.46164 * cos(
u);
658 win = 0.42 + 0.5 * cos(
u) + 0.08 * cos(2*
u);
661 win = 0.40897 + 0.5 * cos(
u) + 0.09103 * cos(2*
u);
664 win = 0.4243801 + 0.4973406 * cos(
u) + 0.0782793 * cos(2*
u);
667 win = 0.355768 + 0.487396 * cos(
u) + 0.144232 * cos(2*
u) + 0.012604 * cos(3*
u);
670 win = 0.3635819 + 0.4891775 * cos(
u) + 0.1365995 * cos(2*
u) + 0.0106411 * cos(3*
u);
673 win = 0.35875 + 0.48829 * cos(
u) + 0.14128 * cos(2*
u) + 0.01168 * cos(3*
u);
676 win = (
u <= 0.5 *
M_PI) ? 1.0 : (0.5 + 0.5 * cos(2*
u -
M_PI));
681 s->analysis_buf[k] *= (2.0/
s->analysis_rdft_len) * (2.0/
s->rdft_len) *
win;
683 s->analysis_buf[
s->analysis_rdft_len - k] =
s->analysis_buf[k];
686 memset(
s->analysis_buf + center + 1, 0, (
s->analysis_rdft_len -
s->fir_len) *
sizeof(*
s->analysis_buf));
687 memcpy(rdft_tbuf,
s->analysis_buf,
s->rdft_len/2 *
sizeof(*
s->analysis_buf));
688 memcpy(rdft_tbuf +
s->rdft_len/2,
s->analysis_buf +
s->analysis_rdft_len -
s->rdft_len/2,
s->rdft_len/2 *
sizeof(*
s->analysis_buf));
691 s->rdft_fn(
s->rdft, rdft_buf, rdft_tbuf,
sizeof(
float));
693 for (k = 0; k <
s->rdft_len + 2; k++) {
694 if (
isnan(rdft_buf[k]) ||
isinf(rdft_buf[k])) {
704 for (k = 0; k <=
s->rdft_len/2; k++)
705 rdft_buf[k] = rdft_buf[2*k];
715 memcpy(
s->kernel_buf,
s->kernel_tmp_buf, (
s->multi ?
inlink->ch_layout.nb_channels : 1) * (
s->rdft_len * 2) *
sizeof(*
s->kernel_buf));
722 #define SELECT_GAIN(s) (s->gain_cmd ? s->gain_cmd : s->gain)
723 #define SELECT_GAIN_ENTRY(s) (s->gain_entry_cmd ? s->gain_entry_cmd : s->gain_entry)
729 float iscale,
scale = 1.f;
735 s->frame_nsamples_max = 0;
737 s->fir_len =
FFMAX(2 * (
int)(
inlink->sample_rate *
s->delay) + 1, 3);
738 s->remaining =
s->fir_len - 1;
741 s->rdft_len = 1 << rdft_bits;
742 s->nsamples_max =
s->rdft_len -
s->fir_len + 1;
743 if (
s->nsamples_max * 2 >=
s->fir_len)
758 if (
s->fft2 && !
s->multi &&
inlink->ch_layout.nb_channels > 1 &&
763 int cepstrum_bits = rdft_bits + 2;
780 s->cepstrum_len = 1 << cepstrum_bits;
782 if (!
s->cepstrum_buf)
785 if (!
s->cepstrum_tbuf)
790 s->analysis_rdft_len = 1 << rdft_bits;
791 if (
inlink->sample_rate <=
s->accuracy *
s->analysis_rdft_len)
811 s->analysis_buf =
av_malloc_array((
s->analysis_rdft_len + 2),
sizeof(*
s->analysis_buf));
812 s->analysis_tbuf =
av_malloc_array(
s->analysis_rdft_len + 2,
sizeof(*
s->analysis_tbuf));
813 s->kernel_tmp_buf =
av_malloc_array((
s->rdft_len * 2) * (
s->multi ?
inlink->ch_layout.nb_channels : 1),
sizeof(*
s->kernel_tmp_buf));
815 s->kernel_buf =
av_malloc_array((
s->rdft_len * 2) * (
s->multi ?
inlink->ch_layout.nb_channels : 1),
sizeof(*
s->kernel_buf));
817 s->conv_buf =
av_calloc(2 *
s->rdft_len *
inlink->ch_layout.nb_channels,
sizeof(*
s->conv_buf));
819 if (!
s->analysis_buf || !
s->analysis_tbuf || !
s->kernel_tmp_buf || !
s->kernel_buf || !
s->conv_buf || !
s->conv_idx || !
s->kernel_tmp_tbuf || !
s->tx_buf)
822 av_log(
ctx,
AV_LOG_DEBUG,
"sample_rate = %d, channels = %d, analysis_rdft_len = %d, rdft_len = %d, fir_len = %d, nsamples_max = %d.\n",
823 inlink->sample_rate,
inlink->ch_layout.nb_channels,
s->analysis_rdft_len,
s->rdft_len,
s->fir_len,
s->nsamples_max);
838 for (ch = 0; ch + 1 <
inlink->ch_layout.nb_channels &&
s->fft_ctx; ch += 2) {
844 for ( ; ch <
inlink->ch_layout.nb_channels; ch++) {
846 s->conv_buf + 2 * ch *
s->rdft_len,
s->conv_idx + ch,
850 for (ch = 0; ch <
inlink->ch_layout.nb_channels; ch++) {
852 s->conv_buf + 2 * ch *
s->rdft_len,
s->conv_idx + ch,
860 if (
s->zero_phase && !
s->min_phase)
890 char *res,
int res_len,
int flags)
895 if (!strcmp(cmd,
"gain")) {
910 s->gain_cmd = gain_cmd;
914 }
else if (!strcmp(cmd,
"gain_entry")) {
915 char *gain_entry_cmd;
929 s->gain_entry_cmd = gain_entry_cmd;
957 .
name =
"firequalizer",
965 .priv_class = &firequalizer_class,