2
$\begingroup$

I've got a sliding DFT implementation that appears to be working (judging from an output plot). I would like to be able to invert this implementation using the standard tricks of swapping the real and imaginary components or conjugating my data. However, when I attempt this my output looks nothing like the original. I also end up with a strong imaginary component where that was not expected. (I'm not entirely sure on the proper range of output for each pass through this filter.) Can you see anything obviously wrong with this code? Thanks for your time.

    [TestMethod]
    public void TestInverseDft()
    {
        var rand = new Random(42);
        var dataSize = 50000;
        var data = new List<Complex>(dataSize);
        for(int i = 0; i < dataSize; i++)
            data.Add(new Complex(rand.NextDouble() * 2.0 - 1.0, 0.0));

        var clone = new List<Complex>(data);

        var n = SlidingDFT(data, 44100, 12, 55, 13);

        Swap(data);
        SlidingDFT(data, 44100, 12, 55, 13);
        Swap(data);

        for(int i = 15000; i < dataSize - 15000; i++)
            Assert.AreEqual(clone[i].Real, data[i].Real / n, 0.00001);
    }

    public static int SlidingDFT(IList<Complex> data, int sampleRate, int binsPerOctave, double minFrequency, int stepsFromMin)
    {
        var Q = 1.0 / (Math.Pow(2.0, 1.0 / binsPerOctave) - 1.0); // denominator * Fk gives size of bin 
        var k = (int)(Q + 0.5);

        // frequency at k:
        var Fk = minFrequency * Math.Pow(2.0, (double)stepsFromMin / binsPerOctave);
        var Nk = (int)(Q * sampleRate / Fk + 0.5);
        Complex exp = Complex.Exp(new Complex(0.0, 2.0 * Math.PI * k / Nk));
        var current = data[0];
        for(int x = 0; x < data.Count - 1; x++)
        {
            var lookahead = x + Nk >= data.Count ? new Complex() : data[x + Nk];
            var fx = (lookahead - current) + data[x]; //  * 0.9999999 if unstable?
            current = data[x + 1];
            data[x + 1] = fx * exp;
        }
        return Nk;
    }

    public static void Conjugate(IList<Complex> data)
    {
        for(int i = 0; i < data.Count; i++)
            data[i] = Complex.Conjugate(data[i]);
    }

    public static void Swap(IList<Complex> data)
    {
        for(int i = 0; i < data.Count; i++)
            data[i] = new Complex(data[i].Imaginary, data[i].Real);
    }
$\endgroup$
  • $\begingroup$ What do you mean by "invert your implementation"? $\endgroup$ – Jim Clay Apr 20 '12 at 1:15
  • $\begingroup$ By invert I mean calculate the inverse Fourier transform. $\endgroup$ – Brannon Apr 20 '12 at 2:59
3
$\begingroup$

I haven't tried to run your code (is that C#?) or looked at it in great detail, but SlidingDFT() looks more like a Goertzel algorithm implementation to me. It should be a red flag to you that you're able to perform the calculations in place. If you were truly doing a sliding DFT, then you would have a data rate expansion by a factor of $N$, the length of the transform; for each input time-domain sample, you generate a full DFT's worth of outputs.

The Goertzel algorithm is an efficient, recursive structure for implementing a single DFT bin at an arbitrary center frequency, which can be useful for applications like DTMF decoding. However, you can't "invert" it in the sense that you're trying to; it is more of a filtering operation, not a reversible transform.

It doesn't seem to make sense to me to even try to invert a sliding DFT in the way you've shown; the outputs of a sliding DFT don't "slide" sample-by-sample the same way as the time-domain inputs do. I'd read up some more on how the DFT and its sliding version works.

$\endgroup$
  • 1
    $\begingroup$ Thanks. Yes, after further reflection, I realized that my plan to reconstruct the original signal when I only kept one of the DFT bins was very wrong. $\endgroup$ – Brannon Apr 20 '12 at 3:46

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.