Merge pull request #47 from swharden/test-filters

Update window functions to match Python/Numpy/Scipy

dev/python/window-functions/make-window-tests.py

@@ -0,0 +1,25 @@
+import numpy as np
+import scipy.signal
+import scipy.signal.windows
+
+
+def printCsArray(name: str, values: np.ndarray, precision: int = 8):
+    arrayContents = ", ".join([str(np.round(x, precision)) for x in values])
+    assignmentString = f"public double[] Known_{name}_{len(values)} = {{ {arrayContents} }};"
+    print(assignmentString)
+
+
+def printCsArrays(count: int):
+    printCsArray("bartlett", scipy.signal.windows.bartlett(count))
+    printCsArray("blackman", scipy.signal.windows.blackman(count))
+    printCsArray("cosine", scipy.signal.windows.cosine(count))
+    printCsArray("flattop", scipy.signal.windows.flattop(count))
+    printCsArray("hamming",  scipy.signal.windows.hamming(count))
+    printCsArray("hanning",  scipy.signal.windows.hann(count))
+    printCsArray("kaiser14", scipy.signal.windows.kaiser(count, beta=14))
+    printCsArray("tukey",  scipy.signal.windows.tukey(count))
+
+
+if __name__ == "__main__":
+    printCsArrays(13)
+    printCsArrays(14)

src/FftSharp.Tests/Window.cs

@@ -72,6 +72,27 @@ public void Test_NormalizedWindows_SumIsOne()
             }
         }
 
+        [Test]
+        public void Test_OddLength_CenterIndexIsBiggest()
+        {
+            foreach (IWindow window in FftSharp.Window.GetWindows())
+            {
+                double[] values = window.Create(13);
+                Assert.GreaterOrEqual(values[6], values[5], window.Name);
+                Assert.GreaterOrEqual(values[6], values[7], window.Name);
+            }
+        }
+
+        [Test]
+        public void Test_EvenLength_CenterTwoAreSame()
+        {
+            foreach (IWindow window in FftSharp.Window.GetWindows())
+            {
+                double[] values = window.Create(12);
+                Assert.AreEqual(values[5], values[6], 1e-5, window.Name);
+            }
+        }
+
         [Test]
         public void Test_Plot_AllWindowKernels()
         {
@@ -100,35 +121,5 @@ public void Test_Window_Reflection()
             IWindow[] window = FftSharp.Window.GetWindows();
             Assert.IsNotEmpty(window);
         }
-
-        [Test]
-        public void Test_Kaiser_MatchesPython()
-        {
-            /* expected values calculated with python:
-                >>> import numpy as np
-                >>> np.kaiser(50, 14)
-            */
-            double[] expected = {
-                7.72686684e-06, 8.15094846e-05, 3.26000767e-04, 9.42588751e-04, 2.26624847e-03,
-                4.80567914e-03, 9.27621459e-03, 1.66164301e-02, 2.79789657e-02, 4.46873500e-02,
-                6.81537432e-02, 9.97574012e-02, 1.40689805e-01, 1.91778970e-01, 2.53311387e-01,
-                3.24874218e-01, 4.05241756e-01, 4.92328143e-01, 5.83222700e-01, 6.74315433e-01,
-                7.61509399e-01, 8.40504954e-01, 9.07130390e-01, 9.57685605e-01, 9.89261639e-01,
-                1.00000000e+00, 9.89261639e-01, 9.57685605e-01, 9.07130390e-01, 8.40504954e-01,
-                7.61509399e-01, 6.74315433e-01, 5.83222700e-01, 4.92328143e-01, 4.05241756e-01,
-                3.24874218e-01, 2.53311387e-01, 1.91778970e-01, 1.40689805e-01, 9.97574012e-02,
-                6.81537432e-02, 4.46873500e-02, 2.79789657e-02, 1.66164301e-02, 9.27621459e-03,
-                4.80567914e-03, 2.26624847e-03, 9.42588751e-04, 3.26000767e-04, 8.15094846e-05,
-            };
-
-            var window = new FftSharp.Windows.Kaiser(14);
-            double[] actual = window.Create(51);
-
-            for (int i = 0; i < expected.Length; i++)
-            {
-                double allowableError = .00001 * expected[i];
-                Assert.AreEqual(expected[i], actual[i], allowableError);
-            }
-        }
     }
 }

src/FftSharp.Tests/WindowValueTests.cs

@@ -0,0 +1,109 @@
+using NUnit.Framework;
+using System;
+using System.Collections.Generic;
+using System.Linq;
+using System.Text;
+using System.Threading.Tasks;
+
+namespace FftSharp.Tests
+{
+    internal class WindowValueTests
+    {
+        /* values in these test functions were calculated using SciPy */
+
+        private void AssertEqual(double[] expected, double[] actual, double precision = 1e-5)
+        {
+            Assert.IsNotNull(expected);
+            Assert.IsNotNull(actual);
+            Assert.AreEqual(expected.Length, actual.Length);
+
+            Console.WriteLine($"Expected\tActual");
+            for (int i = 0; i < expected.Length; i++)
+                Console.WriteLine($"{expected[i]}\t{actual[i]}");
+
+            for (int i = 0; i < expected.Length; i++)
+                Assert.AreEqual(expected[i], actual[i], precision);
+        }
+
+        [Test]
+        public void Test_Bartlett()
+        {
+            double[] Known_bartlett_13 = { 0.0, 0.16666667, 0.33333333, 0.5, 0.66666667, 0.83333333, 1.0, 0.83333333, 0.66666667, 0.5, 0.33333333, 0.16666667, 0.0 };
+            double[] Known_bartlett_14 = { 0.0, 0.15384615, 0.30769231, 0.46153846, 0.61538462, 0.76923077, 0.92307692, 0.92307692, 0.76923077, 0.61538462, 0.46153846, 0.30769231, 0.15384615, 0.0 };
+            AssertEqual(Known_bartlett_13, new Windows.Bartlett().Create(13));
+            AssertEqual(Known_bartlett_14, new Windows.Bartlett().Create(14));
+        }
+
+        [Test]
+        public void Test_Blackman()
+        {
+            double[] Known_blackman_13 = { -0.0, 0.0269873, 0.13, 0.34, 0.63, 0.8930127, 1.0, 0.8930127, 0.63, 0.34, 0.13, 0.0269873, -0.0 };
+            double[] Known_blackman_14 = { -0.0, 0.02271717, 0.10759924, 0.28205631, 0.53742158, 0.80389831, 0.97630739, 0.97630739, 0.80389831, 0.53742158, 0.28205631, 0.10759924, 0.02271717, -0.0 };
+            AssertEqual(Known_blackman_13, new Windows.Blackman(.42, .5, .08).Create(13));
+            AssertEqual(Known_blackman_14, new Windows.Blackman(.42, .5, .08).Create(14));
+        }
+
+        [Test]
+        public void Test_Cosine()
+        {
+            double[] Known_cosine_13 = { 0.12053668, 0.35460489, 0.56806475, 0.74851075, 0.88545603, 0.97094182, 1.0, 0.97094182, 0.88545603, 0.74851075, 0.56806475, 0.35460489, 0.12053668 };
+            double[] Known_cosine_14 = { 0.11196448, 0.33027906, 0.53203208, 0.70710678, 0.8467242, 0.94388333, 0.99371221, 0.99371221, 0.94388333, 0.8467242, 0.70710678, 0.53203208, 0.33027906, 0.11196448 };
+            AssertEqual(Known_cosine_13, new Windows.Cosine().Create(13));
+            AssertEqual(Known_cosine_14, new Windows.Cosine().Create(14));
+        }
+
+        [Test]
+        public void Test_FlatTop()
+        {
+            double[] Known_flattop_13 = { -0.00042105, -0.01007669, -0.05126316, -0.05473684, 0.19821053, 0.71155038, 1.0, 0.71155038, 0.19821053, -0.05473684, -0.05126316, -0.01007669, -0.00042105 };
+            double[] Known_flattop_14 = { -0.00042105, -0.00836447, -0.04346352, -0.06805774, 0.08261602, 0.5066288, 0.9321146, 0.9321146, 0.5066288, 0.08261602, -0.06805774, -0.04346352, -0.00836447, -0.00042105 };
+            AssertEqual(Known_flattop_13, new Windows.FlatTop().Create(13));
+            AssertEqual(Known_flattop_14, new Windows.FlatTop().Create(14));
+        }
+
+        [Test]
+        public void Test_Hamming()
+        {
+            double[] Known_hamming_13 = { 0.08, 0.14162831, 0.31, 0.54, 0.77, 0.93837169, 1.0, 0.93837169, 0.77, 0.54, 0.31, 0.14162831, 0.08 };
+            double[] Known_hamming_14 = { 0.08, 0.13269023, 0.27869022, 0.48455313, 0.70311825, 0.88431494, 0.98663324, 0.98663324, 0.88431494, 0.70311825, 0.48455313, 0.27869022, 0.13269023, 0.08 };
+            AssertEqual(Known_hamming_13, new Windows.Hamming().Create(13));
+            AssertEqual(Known_hamming_14, new Windows.Hamming().Create(14));
+        }
+
+        [Test]
+        public void Test_Hanning()
+        {
+            double[] Known_hanning_13 = { 0.0, 0.0669873, 0.25, 0.5, 0.75, 0.9330127, 1.0, 0.9330127, 0.75, 0.5, 0.25, 0.0669873, 0.0 };
+            double[] Known_hanning_14 = { 0.0, 0.05727199, 0.21596763, 0.43973166, 0.67730244, 0.87425537, 0.98547091, 0.98547091, 0.87425537, 0.67730244, 0.43973166, 0.21596763, 0.05727199, 0.0 };
+            AssertEqual(Known_hanning_13, new Windows.Hanning().Create(13));
+            AssertEqual(Known_hanning_14, new Windows.Hanning().Create(14));
+        }
+
+        [Test]
+        public void Test_Rectangular()
+        {
+            double[] Known_rectangular_13 = Enumerable.Range(0, 13).Select(x => 1.0).ToArray();
+            double[] Known_rectangular_14 = Enumerable.Range(0, 14).Select(x => 1.0).ToArray();
+            AssertEqual(Known_rectangular_13, new Windows.Rectangular().Create(13));
+            AssertEqual(Known_rectangular_14, new Windows.Rectangular().Create(14));
+        }
+
+        [Test]
+        public void Test_Tukey()
+        {
+            double[] Known_tukey_13 = { 0.0, 0.25, 0.75, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 0.75, 0.25, 0.0 };
+            double[] Known_tukey_14 = { 0.0, 0.21596763, 0.67730244, 0.98547091, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 0.98547091, 0.67730244, 0.21596763, 0.0 };
+            AssertEqual(Known_tukey_13, new Windows.Tukey().Create(13));
+            AssertEqual(Known_tukey_14, new Windows.Tukey().Create(14));
+        }
+
+        [Test]
+        public void Test_Kaiser()
+        {
+            double[] Known_kaiser14_13 = { 7.73e-06, 0.00258844, 0.03288553, 0.16493219, 0.4627165, 0.82808941, 1.0, 0.82808941, 0.4627165, 0.16493219, 0.03288553, 0.00258844, 7.73e-06 };
+            double[] Known_kaiser14_14 = { 7.73e-06, 0.00199846, 0.02397752, 0.12057536, 0.35483775, 0.69493873, 0.96081903, 0.96081903, 0.69493873, 0.35483775, 0.12057536, 0.02397752, 0.00199846, 7.73e-06 };
+            AssertEqual(Known_kaiser14_13, new Windows.Kaiser(beta: 14).Create(13));
+            AssertEqual(Known_kaiser14_14, new Windows.Kaiser(beta: 14).Create(14));
+        }
+    }
+}

src/FftSharp/Windows/Bartlett.cs

@@ -13,8 +13,12 @@ public override double[] Create(int size, bool normalize = false)
         {
             double[] window = new double[size];
 
+            bool isOddSize = size % 2 == 1;
+
+            double halfSize = isOddSize ? size / 2 : (size - 1) / 2.0;
+
             for (int i = 0; i < size; i++)
-                window[i] = 1 - Math.Abs((double)(i - (size / 2)) / (size / 2));
+                window[i] = 1 - Math.Abs((double)(i - halfSize) / halfSize);
 
             if (normalize)
                 NormalizeInPlace(window);

src/FftSharp/Windows/Blackman.cs

@@ -17,6 +17,7 @@ public Blackman()
         {
         }
 
+        //TODO: 5-term constructor to allow testing Python's flattop
         public Blackman(double a, double b, double c)
         {
             (A, B, C) = (a, b, c);
@@ -27,7 +28,10 @@ public override double[] Create(int size, bool normalize = false)
             double[] window = new double[size];
 
             for (int i = 0; i < size; i++)
-                window[i] = A - B * Math.Cos(2 * Math.PI * i / size) + C * Math.Cos(4 * Math.PI * i / size);
+            {
+                double frac = (double)i / (size - 1);
+                window[i] = A - B * Math.Cos(2 * Math.PI * frac) + C * Math.Cos(4 * Math.PI * frac);
+            }
 
             if (normalize)
                 NormalizeInPlace(window);

src/FftSharp/Windows/Cosine.cs

@@ -1,4 +1,5 @@
 using System;
+using System.Linq;
 
 namespace FftSharp.Windows
 {
@@ -11,10 +12,7 @@ public class Cosine : Window, IWindow
 
         public override double[] Create(int size, bool normalize = false)
         {
-            double[] window = new double[size];
-
-            for (int i = 0; i < size; i++)
-                window[i] = Math.Sin(i * Math.PI / (size - 1));
+            double[] window = Enumerable.Range(0, size).Select(x => Math.Sin(Math.PI / (size) * (x + .5))).ToArray();
 
             if (normalize)
                 NormalizeInPlace(window);

src/FftSharp/Windows/FlatTop.cs

@@ -4,7 +4,7 @@
 
 namespace FftSharp.Windows
 {
-    public class FlatTop : Blackman
+    public class FlatTop : Window, IWindow
     {
         public override string Name => "Flat Top";
         public override string Description =>
@@ -14,9 +14,42 @@ public class FlatTop : Blackman
             "The flat top window crosses the zero line causing a broader peak in the frequency domain, " +
             "which is closer to the true amplitude of the signal than with other windows";
 
-        public FlatTop() : base(0.2810639, 0.5208972, 0.1980399)
+        public readonly double A0 = 0.21557895;
+        public readonly double A1 = 0.41663158;
+        public readonly double A2 = 0.277263158;
+        public readonly double A3 = 0.083578947;
+        public readonly double A4 = 0.006947368;
+
+        public FlatTop()
         {
+        }
+
+        public FlatTop(double a0, double a1, double a2, double a3, double a4)
+        {
+            A0 = a0;
+            A1 = a1;
+            A2 = a2;
+            A3 = a3;
+            A4 = a4;
+        }
+
+        public override double[] Create(int size, bool normalize = false)
+        {
+            double[] window = new double[size];
+
+            for (int i = 0; i < size; i++)
+            {
+                window[i] = A0
+                    - A1 * Math.Cos(2 * Math.PI * i / (size - 1))
+                    + A2 * Math.Cos(4 * Math.PI * i / (size - 1))
+                    - A3 * Math.Cos(6 * Math.PI * i / (size - 1))
+                    + A4 * Math.Cos(8 * Math.PI * i / (size - 1));
+            }
+
+            if (normalize)
+                NormalizeInPlace(window);
 
+            return window;
         }
     }
 }

src/FftSharp/Windows/Hamming.cs

@@ -15,7 +15,7 @@ public override double[] Create(int size, bool normalize = false)
             double[] window = new double[size];
 
             for (int i = 0; i < size; i++)
-                window[i] = 0.54 - 0.46 * Math.Cos(2 * Math.PI * i / size);
+                window[i] = 0.54 - 0.46 * Math.Cos(2 * Math.PI * i / (size - 1));
 
             if (normalize)
                 NormalizeInPlace(window);

src/FftSharp/Windows/Hanning.cs

@@ -15,7 +15,7 @@ public override double[] Create(int size, bool normalize = false)
             double[] window = new double[size];
 
             for (int i = 0; i < size; i++)
-                window[i] = 0.5 - 0.5 * Math.Cos(2 * Math.PI * i / size);
+                window[i] = 0.5 - 0.5 * Math.Cos(2 * Math.PI * i / (size - 1));
 
             if (normalize)
                 NormalizeInPlace(window);

src/FftSharp/Windows/Tukey.cs

@@ -29,13 +29,29 @@ public override double[] Create(int size, bool normalize = false)
         {
             double[] window = new double[size];
 
-            double m = 2 * Math.PI / (Alpha * size);
+            double m = 2 * Math.PI / (Alpha * (size - 1));
+
             int edgeSizePoints = (int)(size * Alpha / 2);
 
+            if (size % 2 == 0)
+                edgeSizePoints += 1;
+
             for (int i = 0; i < size; i++)
             {
-                bool isEdge = (i < edgeSizePoints) || (i > size - edgeSizePoints);
-                window[i] = isEdge ? (1 - Math.Cos(i * m)) / 2 : 1;
+                if (i < edgeSizePoints)
+                {
+                    // left edge
+                    window[i] = (1 - Math.Cos(i * m)) / 2;
+                }
+                else if (i >= size - edgeSizePoints)
+                {
+                    // right edge
+                    window[i] = (1 - Math.Cos(i * m)) / 2;
+                }
+                else
+                {
+                    window[i] = 1;
+                }
             }
 
             if (normalize)