Merge pull request #74 from gsgou/inverse_real

implements InverseReal

src/FftSharp.Tests/Inverse.cs

@@ -51,5 +51,24 @@ public void Test_IFFT_MatchesOriginal()
                 Assert.AreEqual(original[i].Imaginary, ifft[i].Imaginary, 1e-6);
             }
         }
+
+        [Test]
+        public void Test_InverseReal_MatchesOriginal()
+        {
+            Random rand = new Random(0);
+            double[] original = new double[1024];
+            for (int i = 0; i < original.Length; i++)
+                original[i] = rand.NextDouble() - .5;
+
+            System.Numerics.Complex[] rfft = FftSharp.FFT.ForwardReal(original);
+            double[] irfft = FftSharp.FFT.InverseReal(rfft);
+
+            Assert.AreEqual(original.Length, irfft.Length);
+
+            for (int i = 0; i < irfft.Length; i++)
+            {
+                Assert.AreEqual(original[i], irfft[i], 1e-10);
+            }
+        }
     }
 }

src/FftSharp/FFT.cs

@@ -1,4 +1,4 @@
-using System;
+using System;
 
 namespace FftSharp;
 
@@ -73,6 +73,37 @@ public static System.Numerics.Complex[] ForwardReal(double[] samples)
         return real;
     }
 
+    /// <summary>
+    /// Calculate IFFT and return just the real component of the spectrum
+    /// </summary>
+    public static double[] InverseReal(System.Numerics.Complex[] spectrum)
+    {
+        int spectrumSize = spectrum.Length;
+        int windowSize = (spectrumSize - 1) * 2;
+
+        if (!FftOperations.IsPowerOfTwo(windowSize))
+            throw new ArgumentException($"{nameof(spectrum)} length must be a power of two plus 1");
+
+        System.Numerics.Complex[] buffer = new System.Numerics.Complex[windowSize];
+
+        for (int i = 0; i < spectrumSize; i++)
+            buffer[i] = spectrum[i];
+
+        for (int i = spectrumSize; i < windowSize; i++)
+        {
+            int iMirrored = spectrumSize - (i - (spectrumSize - 2));
+            buffer[i] = System.Numerics.Complex.Conjugate(spectrum[iMirrored]);
+        }
+
+        Inverse(buffer);
+
+        double[] result = new double[windowSize];
+        for (int i = 0; i < windowSize; i++)
+            result[i] = buffer[i].Real;
+
+        return result;
+    }
+
     /// <summary>
     /// Apply the inverse Fast Fourier Transform (iFFT) to the Complex array in place.
     /// Length of the array must be a power of 2.