changeset 48:f97abcda094f

Implement and test inverse fft

--- a/complex.yeti	Sat Dec 29 15:36:05 2012 +0000
+++ b/complex.yeti	Sat Dec 29 23:06:22 2012 +0000
@@ -9,7 +9,19 @@
     int getImag()
         imag,
     String toString()
-        "\(real) + \(imag)i",
+        if real == int real and imag == int imag then
+            if imag < 0 then
+                " \(int real) - \(int (-imag))i"
+            else 
+                " \(int real) + \(int imag)i"
+            fi
+        else
+            if imag < 0 then
+                " \(real) - \((-imag))i"
+            else 
+                " \(real) + \(imag)i"
+            fi
+        fi,
     int hashCode()
         Double#valueOf(real)#hashCode() + Double#valueOf(imag)#hashCode(),
     boolean equals(Object other)

--- a/fft.yeti	Sat Dec 29 15:36:05 2012 +0000
+++ b/fft.yeti	Sat Dec 29 23:06:22 2012 +0000
@@ -7,7 +7,7 @@
 vec = load fvector;
 complex = load complex;
 
-unpack p =
+packedToComplex p =
    (n = (vec.length p) / 2;
     array
        (map do i:
@@ -16,15 +16,38 @@
             complex.complex re im;
         done [0..n]));
 
-fft n = 
+complexToPacked arr =
+   (n = length arr;
+    v = vec.vector
+       (map do i:
+            ix = int (i/2);
+            if i == ix*2 then
+                complex.real arr[ix]
+            else 
+                complex.imaginary arr[ix] 
+            fi;
+            done [0..(n-1)*2-1]);
+    v[1] := complex.real arr[n-1];
+    v);
+
+realForward n = 
    (d = new DoubleFFT_1D(n);
     do bl:
         v = b.vector bl;
         d#realForward(v);
-        unpack v;
+        packedToComplex v;
+    done);
+
+realInverse n = 
+   (d = new DoubleFFT_1D(n);
+    do arr:
+        v = complexToPacked arr;
+        d#realInverse(v, true);
+        b.block v;
     done);
 
 {
-fft,
+realForward,
+realInverse,
 }
 

--- a/test/test_fft.yeti	Sat Dec 29 15:36:05 2012 +0000
+++ b/test/test_fft.yeti	Sat Dec 29 23:06:22 2012 +0000
@@ -1,43 +1,44 @@
 
 module test.test_fft;
 
-{ fft } = load fft;
-{ fromList } = load block;
-{ list } = load fvector;
+{ realForward, realInverse } = load fft;
+{ list, fromList } = load block;
 { complex } = load complex;
 
 { declare, compare } = load test.test;
 
+testFFT orig reals imags =
+   (out = realForward (length orig) (fromList orig);
+    back = realInverse (length orig) out;
+    compare out (array (map2 complex reals imags)) and compare orig (list back));
+
 declare [
 
 "dc": \(
-    out = fft 4 (fromList [1,1,1,1]);
-    compare out (array (map2 complex [4,0,0] [0,0,0]))
+    testFFT [1,1,1,1] [4,0,0] [0,0,0];
 ),
 
 "sine": \(
-    out = fft 4 (fromList [0,1,0,-1]);
-    compare out (array (map2 complex [0,0,0] [0,-2,0]))
+    testFFT [0,1,0,-1] [0,0,0] [0,-2,0];
 ),
 
 "cosine": \(
-    out = fft 4 (fromList [1,0,-1,0]);
-    compare out (array (map2 complex [0,2,0] [0,0,0]))
+    testFFT [1,0,-1,0] [0,2,0] [0,0,0];
 ),
 
 "sineCosine": \(
-    out = fft 4 (fromList [0.5,1,-0.5,-1]);
-    compare out (array (map2 complex [0,1,0] [0,-2,0]))
+    testFFT [0.5,1,-0.5,-1] [0,1,0] [0,-2,0];
 ),
 
 "nyquist": \(
-    out = fft 4 (fromList [1,-1,1,-1]);
-    compare out (array (map2 complex [0,0,4] [0,0,0]))
+    testFFT [1,-1,1,-1] [0,0,4] [0,0,0];
 ),
 
 "dirac": \(
-    out = fft 4 (fromList [1,0,0,0]);
-    compare out (array (map2 complex [1,1,1] [0,0,0]))
+    testFFT [1,0,0,0] [1,1,1] [0,0,0] and
+        testFFT [0,1,0,0] [1,0,-1] [0,-1,0] and
+        testFFT [0,0,1,0] [1,-1,1] [0,0,0] and
+        testFFT [0,0,0,1] [1,0,-1] [0,1,0];
 ),
 
 ];