lavu/tx: support in-place FFT transforms

This commit adds support for in-place FFT transforms. Since our
internal transforms were all in-place anyway, this only changes
the permutation on the input.

Unfortunately, research papers were of no help here. All focused
on dry hardware implementations, where permutes are free, or on
software implementations where binary bloat is of no concern so
storing dozen times the transforms for each permutation and version
is not considered bad practice.
Still, for a pure C implementation, it's only around 28% slower
than the multi-megabyte FFTW3 in unaligned mode.

Unlike a closed permutation like with PFA, split-radix FFT bit-reversals
contain multiple NOPs, multiple simple swaps, and a few chained swaps,
so regular single-loop single-state permute loops were not possible.
Instead, we filter out parts of the input indices which are redundant.
This allows for a single branch, and with some clever AVX512 asm,
could possibly be SIMD'd without refactoring.

The inplace_idx array is guaranteed to never be larger than the
revtab array, and in practice only requires around log2(len) entries.

The power-of-two MDCTs can be done in-place as well. And it's
possible to eliminate a copy in the compound MDCTs too, however
it'll be slower than doing them out of place, and we'd need to dirty
the input array.

1 files changed, 33 insertions, 3 deletions

diff --git a/libavutil/tx_template.c b/libavutil/tx_template.c
index 155e879f8e..f43173e920 100644
--- a/libavutil/tx_template.c
+++ b/libavutil/tx_template.c
@@ -392,8 +392,28 @@ static void monolithic_fft(AVTXContext *s, void *_out, void *_in,
     FFTComplex *in = _in;
     FFTComplex *out = _out;
     int m = s->m, mb = av_log2(m);
-    for (int i = 0; i < m; i++)
-        out[s->revtab[i]] = in[i];
+
+    if (s->flags & AV_TX_INPLACE) {
+        FFTComplex tmp;
+        int src, dst, *inplace_idx = s->inplace_idx;
+
+        out = in;
+        src = *inplace_idx++;
+
+        do {
+            tmp = out[src];
+            dst = s->revtab[src];
+            do {
+                FFSWAP(FFTComplex, tmp, out[dst]);
+                dst = s->revtab[dst];
+            } while (dst != src); /* Can be > as well, but is less predictable */
+            out[dst] = tmp;
+        } while ((src = *inplace_idx++));
+    } else {
+        for (int i = 0; i < m; i++)
+            out[s->revtab[i]] = in[i];
+    }
+
     fft_dispatch[mb](out);
 }
 
@@ -678,6 +698,7 @@ int TX_NAME(ff_tx_init_mdct_fft)(AVTXContext *s, av_tx_fn *tx,
     s->m = m;
     s->inv = inv;
     s->type = type;
+    s->flags = flags;
 
     /* If we weren't able to split the length into factors we can handle,
      * resort to using the naive and slow FT. This also filters out
@@ -685,6 +706,8 @@ int TX_NAME(ff_tx_init_mdct_fft)(AVTXContext *s, av_tx_fn *tx,
     if (len > 1 || m == 1) {
         if (is_mdct && (l & 1)) /* Odd (i)MDCTs are not supported yet */
             return AVERROR(ENOSYS);
+        if (flags & AV_TX_INPLACE) /* Neither are in-place naive transforms */
+            return AVERROR(ENOSYS);
         s->n = l;
         s->m = 1;
         *tx = naive_fft;
@@ -716,7 +739,14 @@ int TX_NAME(ff_tx_init_mdct_fft)(AVTXContext *s, av_tx_fn *tx,
     if (n != 1)
         init_cos_tabs(0);
     if (m != 1) {
-        ff_tx_gen_ptwo_revtab(s);
+        if ((err = ff_tx_gen_ptwo_revtab(s)))
+            return err;
+        if (flags & AV_TX_INPLACE) {
+            if (is_mdct) /* In-place MDCTs are not supported yet */
+                return AVERROR(ENOSYS);
+            if ((err = ff_tx_gen_ptwo_inplace_revtab_idx(s)))
+                return err;
+        }
         for (int i = 4; i <= av_log2(m); i++)
             init_cos_tabs(i);
     }