[Commit-gnuradio] r8211 - in gnuradio/trunk/gcell/src/lib/wrapper: . spu

[eb]

Author: eb
Date: 2008-04-15 23:32:54 -0600 (Tue, 15 Apr 2008)
New Revision: 8211

Modified:
   gnuradio/trunk/gcell/src/lib/wrapper/gcp_fft_1d_r2.cc
   gnuradio/trunk/gcell/src/lib/wrapper/gcp_fft_1d_r2.h
   gnuradio/trunk/gcell/src/lib/wrapper/qa_gcp_fft_1d_r2.cc
   gnuradio/trunk/gcell/src/lib/wrapper/spu/gcs_fft_1d_r2.c
Log:
FFT inverse tranform is now correct.


Modified: gnuradio/trunk/gcell/src/lib/wrapper/gcp_fft_1d_r2.cc
===================================================================
--- gnuradio/trunk/gcell/src/lib/wrapper/gcp_fft_1d_r2.cc       2008-04-16 
03:52:57 UTC (rev 8210)
+++ gnuradio/trunk/gcell/src/lib/wrapper/gcp_fft_1d_r2.cc       2008-04-16 
05:32:54 UTC (rev 8211)
@@ -89,7 +89,7 @@
 }
 
 void
-gcp_fft_1d_r2_forward_twiddle(unsigned int log2_fft_length, 
std::complex<float> *W)
+gcp_fft_1d_r2_twiddle(unsigned int log2_fft_length, std::complex<float> *W)
 {
   unsigned int n = 1 << log2_fft_length;
 
@@ -100,17 +100,3 @@
     W[n/4 - i].imag() = -W[i].real();
   }
 }
-
-
-void
-gcp_fft_1d_r2_reverse_twiddle(unsigned int log2_fft_length, 
std::complex<float> *W)
-{
-  // FIXME this is wrong/insufficient.  inverse is still incorrect
-
-  // reverse factors are the conjugate of the forward factors
-  gcp_fft_1d_r2_forward_twiddle(log2_fft_length, W);
-
-  unsigned int n = 1 << log2_fft_length;
-  for (unsigned i=0; i < n/4; i++)
-    W[i] = conj(W[i]);
-}

Modified: gnuradio/trunk/gcell/src/lib/wrapper/gcp_fft_1d_r2.h
===================================================================
--- gnuradio/trunk/gcell/src/lib/wrapper/gcp_fft_1d_r2.h        2008-04-16 
03:52:57 UTC (rev 8210)
+++ gnuradio/trunk/gcell/src/lib/wrapper/gcp_fft_1d_r2.h        2008-04-16 
05:32:54 UTC (rev 8211)
@@ -28,7 +28,7 @@
  * \brief Submit a job that computes the forward or reverse FFT.
  *
  * \param mgr is the job manager instance
- * \param log2_fft_length is the log2 of the fft_length (4 <= x <= 13).
+ * \param log2_fft_length is the log2 of the fft_length (4 <= x <= 12).
  * \param forward is true to compute the forward xform
  * \param out is the fft_length output from FFT (must be 16-byte aligned).
  * \param in is the fft_length input to FFT (must be 16-byte aligned).
@@ -46,21 +46,12 @@
                     const std::complex<float> *W);
 
 /*!
- * \brief Compute twiddle factors for forward transform.
+ * \brief Compute twiddle factors 
  *
  * \param log2_fft_length is the log2 of the fft_length.
  * \param W is fft_length/4 twiddle factor output (must be 16-byte aligned).
  */
 void
-gcp_fft_1d_r2_forward_twiddle(unsigned int log2_fft_length, 
std::complex<float> *W);
+gcp_fft_1d_r2_twiddle(unsigned int log2_fft_length, std::complex<float> *W);
 
-/*!
- * \brief Compute twiddle factors for reverse transform.
- *
- * \param log2_fft_length is the log2 of the fft_length.
- * \param W is fft_length/4 twiddle factor output (must be 16-byte aligned).
- */
-void
-gcp_fft_1d_r2_reverse_twiddle(unsigned int log2_fft_length, 
std::complex<float> *W);
-
 #endif /* INCLUDED_GCP_FFT_1D_R2_H */

Modified: gnuradio/trunk/gcell/src/lib/wrapper/qa_gcp_fft_1d_r2.cc
===================================================================
--- gnuradio/trunk/gcell/src/lib/wrapper/qa_gcp_fft_1d_r2.cc    2008-04-16 
03:52:57 UTC (rev 8210)
+++ gnuradio/trunk/gcell/src/lib/wrapper/qa_gcp_fft_1d_r2.cc    2008-04-16 
05:32:54 UTC (rev 8211)
@@ -150,10 +150,7 @@
   std::complex<float> *cell_out = (std::complex<float> *) cell_out_void.get();
   std::complex<float> *cell_twiddle = (std::complex<float> *) 
cell_twiddle_void.get();
 
-  if (forward)
-    gcp_fft_1d_r2_forward_twiddle(log2_fft_size, cell_twiddle);
-  else
-    gcp_fft_1d_r2_reverse_twiddle(log2_fft_size, cell_twiddle);
+  gcp_fft_1d_r2_twiddle(log2_fft_size, cell_twiddle);
 
   srandom(1);          // we want reproducibility
 
@@ -206,6 +203,5 @@
   fprintf(stdout, "%s fft_size = %4d  max_rel_error = %e\n",
          forward ? "fwd" : "rev", fft_size, max_rel);
 
-  // CPPUNIT_ASSERT(max_rel <= 1e-4);
-
+  CPPUNIT_ASSERT(max_rel <= 5e-3);
 }

Modified: gnuradio/trunk/gcell/src/lib/wrapper/spu/gcs_fft_1d_r2.c
===================================================================
--- gnuradio/trunk/gcell/src/lib/wrapper/spu/gcs_fft_1d_r2.c    2008-04-16 
03:52:57 UTC (rev 8210)
+++ gnuradio/trunk/gcell/src/lib/wrapper/spu/gcs_fft_1d_r2.c    2008-04-16 
05:32:54 UTC (rev 8211)
@@ -22,7 +22,19 @@
 #include <gc_declare_proc.h>
 #include <libfft.h>
 
+/*
+ * v is really vector complex<float>
+ */
 static void
+conjugate_vector(vector float *v, int nelements)
+{
+  vector float k = {1, -1, 1, -1};
+  int i;
+  for (i = 0; i < nelements; i++)
+    v[i] *= k;
+}
+
+static void
 gcs_fft_1d_r2(const gc_job_direct_args_t *input,
              gc_job_direct_args_t *output __attribute__((unused)),
              const gc_job_ea_args_t *eaa)
@@ -31,9 +43,16 @@
   vector float *in = (vector float *) eaa->arg[1].ls_addr;
   vector float *W = (vector float *) eaa->arg[2].ls_addr;
   int log2_fft_length = input->arg[0].u32;
-  int forward = input->arg[1].u32;     // non-zero if forward xform (FIXME use)
+  int forward = input->arg[1].u32;     // non-zero if forward xform
 
-  fft_1d_r2(out, in, W, log2_fft_length);
+  if (forward){
+    fft_1d_r2(out, in, W, log2_fft_length);
+  }
+  else {
+    conjugate_vector(in, 1 << (log2_fft_length - 1));
+    fft_1d_r2(out, in, W, log2_fft_length);
+    conjugate_vector(out, 1 << (log2_fft_length - 1));
+  }
 }
 
 GC_DECLARE_PROC(gcs_fft_1d_r2, "fft_1d_r2");