+static void fast_convolute2(FIREqualizerContext *av_restrict s, const float *av_restrict kernel_buf, FFTComplex *av_restrict conv_buf,
+ OverlapIndex *av_restrict idx, float *av_restrict data0, float *av_restrict data1, int nsamples)
+{
+ if (nsamples <= s->nsamples_max) {
+ FFTComplex *buf = conv_buf + idx->buf_idx * s->rdft_len;
+ FFTComplex *obuf = conv_buf + !idx->buf_idx * s->rdft_len + idx->overlap_idx;
+ int center = s->fir_len/2;
+ int k;
+ float tmp;
+
+ memset(buf, 0, center * sizeof(*buf));
+ for (k = 0; k < nsamples; k++) {
+ buf[center+k].re = data0[k];
+ buf[center+k].im = data1[k];
+ }
+ memset(buf + center + nsamples, 0, (s->rdft_len - nsamples - center) * sizeof(*buf));
+ av_fft_permute(s->fft_ctx, buf);
+ av_fft_calc(s->fft_ctx, buf);
+
+ /* swap re <-> im, do backward fft using forward fft_ctx */
+ /* normalize with 0.5f */
+ tmp = buf[0].re;
+ buf[0].re = 0.5f * kernel_buf[0] * buf[0].im;
+ buf[0].im = 0.5f * kernel_buf[0] * tmp;
+ for (k = 1; k < s->rdft_len/2; k++) {
+ int m = s->rdft_len - k;
+ tmp = buf[k].re;
+ buf[k].re = 0.5f * kernel_buf[k] * buf[k].im;
+ buf[k].im = 0.5f * kernel_buf[k] * tmp;
+ tmp = buf[m].re;
+ buf[m].re = 0.5f * kernel_buf[k] * buf[m].im;
+ buf[m].im = 0.5f * kernel_buf[k] * tmp;
+ }
+ tmp = buf[k].re;
+ buf[k].re = 0.5f * kernel_buf[k] * buf[k].im;
+ buf[k].im = 0.5f * kernel_buf[k] * tmp;
+
+ av_fft_permute(s->fft_ctx, buf);
+ av_fft_calc(s->fft_ctx, buf);
+
+ for (k = 0; k < s->rdft_len - idx->overlap_idx; k++) {
+ buf[k].re += obuf[k].re;
+ buf[k].im += obuf[k].im;
+ }
+
+ /* swapped re <-> im */
+ for (k = 0; k < nsamples; k++) {
+ data0[k] = buf[k].im;
+ data1[k] = buf[k].re;
+ }
+ idx->buf_idx = !idx->buf_idx;
+ idx->overlap_idx = nsamples;
+ } else {
+ while (nsamples > s->nsamples_max * 2) {
+ fast_convolute2(s, kernel_buf, conv_buf, idx, data0, data1, s->nsamples_max);
+ data0 += s->nsamples_max;
+ data1 += s->nsamples_max;
+ nsamples -= s->nsamples_max;
+ }
+ fast_convolute2(s, kernel_buf, conv_buf, idx, data0, data1, nsamples/2);
+ fast_convolute2(s, kernel_buf, conv_buf, idx, data0 + nsamples/2, data1 + nsamples/2, nsamples - nsamples/2);
+ }
+}
+
+static void dump_fir(AVFilterContext *ctx, FILE *fp, int ch)
+{
+ FIREqualizerContext *s = ctx->priv;
+ int rate = ctx->inputs[0]->sample_rate;
+ int xlog = s->dumpscale == SCALE_LOGLIN || s->dumpscale == SCALE_LOGLOG;
+ int ylog = s->dumpscale == SCALE_LINLOG || s->dumpscale == SCALE_LOGLOG;
+ int x;
+ int center = s->fir_len / 2;
+ double delay = s->zero_phase ? 0.0 : (double) center / rate;
+ double vx, ya, yb;
+
+ s->analysis_buf[0] *= s->rdft_len/2;
+ for (x = 1; x <= center; x++) {
+ s->analysis_buf[x] *= s->rdft_len/2;
+ s->analysis_buf[s->analysis_rdft_len - x] *= s->rdft_len/2;
+ }
+
+ if (ch)
+ fprintf(fp, "\n\n");
+
+ fprintf(fp, "# time[%d] (time amplitude)\n", ch);
+
+ for (x = center; x > 0; x--)
+ fprintf(fp, "%15.10f %15.10f\n", delay - (double) x / rate, (double) s->analysis_buf[s->analysis_rdft_len - x]);
+
+ for (x = 0; x <= center; x++)
+ fprintf(fp, "%15.10f %15.10f\n", delay + (double)x / rate , (double) s->analysis_buf[x]);
+
+ av_rdft_calc(s->analysis_rdft, s->analysis_buf);
+
+ fprintf(fp, "\n\n# freq[%d] (frequency desired_gain actual_gain)\n", ch);
+
+ for (x = 0; x <= s->analysis_rdft_len/2; x++) {
+ int i = (x == s->analysis_rdft_len/2) ? 1 : 2 * x;
+ vx = (double)x * rate / s->analysis_rdft_len;
+ if (xlog)
+ vx = log2(0.05*vx);
+ ya = s->dump_buf[i];
+ yb = s->analysis_buf[i];
+ if (ylog) {
+ ya = 20.0 * log10(fabs(ya));
+ yb = 20.0 * log10(fabs(yb));
+ }
+ fprintf(fp, "%17.10f %17.10f %17.10f\n", vx, ya, yb);
+ }
+}
+