examples/test_fourier/test_fourier.c

Go to the documentation of this file.
00001 /*
00002    libit - Library for basic source and channel coding functions
00003    Copyright (C) 2005-2005 Vivien Chappelier, Herve Jegou
00004 
00005    This library is free software; you can redistribute it and/or
00006    modify it under the terms of the GNU Library General Public
00007    License as published by the Free Software Foundation; either
00008    version 2 of the License, or (at your option) any later version.
00009 
00010    This library is distributed in the hope that it will be useful,
00011    but WITHOUT ANY WARRANTY; without even the implied warranty of
00012    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00013    Library General Public License for more details.
00014 
00015    You should have received a copy of the GNU Library General Public
00016    License along with this library; if not, write to the Free
00017    Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
00018 */
00019 
00020 /*
00021   Test program for the Fourier transform
00022   Copyright (C) 2005 Vivien Chappelier
00023 */
00024 
00025 #include <it/io.h>
00026 #include <it/mat.h>
00027 #include <it/transform.h>
00028 #include <it/fourier.h>
00029 
00030 int cvec_eq_tol (cvec v1, cvec v2, double tol)
00031 {
00032   idx_t i;
00033   assert (v1);
00034   assert (v2);
00035 
00036   if (cvec_length (v1) != cvec_length (v2))
00037     return 0;
00038 
00039   for (i = 0; i < cvec_length (v1); i++)
00040     if (cnorm (csub (v1[i], v2[i])) > tol)
00041       return 0;
00042   return 1;
00043 }
00044 
00045 //#define FFT_TIMING
00046 //#define DFT_TIMING
00047 
00048 int main ()
00049 {
00050   cvec v, r1, r2;
00051   cvec vt;
00052   vec  a;
00053   cplx z;
00054 
00055 #ifdef FFT_TIMING
00056   cvec t = cvec_new_range (8 * 1024 * 1024);
00057   vt = it_fft (t);
00058   cvec_delete (vt);
00059   return 0;
00060 #endif
00061 
00062 #ifdef DFT_TIMING
00063   int  i;
00064   cvec t = cvec_new_range (15);
00065   for (i = 0; i < 1000000; i++) {
00066     vt = it_dft (t);
00067     cvec_delete (vt);
00068   }
00069   return 0;
00070 #endif
00071 
00072   /* -- real vector -- */
00073   v = cvec_new_string ("1 2 3");
00074   a = vec_new_string ("1 2 3");
00075   z.r = 0.5;
00076   z.i = 1.0;
00077   it_printf ("z = %z rho = %f theta = %f\n", z, cnorm (z), cang (z));
00078   vt = it_fzt (v, z);
00079   it_printf ("vt[fzt] = $.04z\n", vt);
00080   cvec_delete (vt);
00081 
00082   vt = it_dft (v);
00083   it_printf ("vt[dft] = $.04z\n", vt);
00084   cvec_delete (vt);
00085 
00086   vt = it_dft_real (a);
00087   it_printf ("vt[dft_real] = $.04z\n", vt);
00088 
00089 
00090   v = it_idft (vt);
00091   it_printf ("v[idft] = $.04z\n", v);
00092 
00093   a = it_idft_real (vt);
00094   it_printf ("a[idft_real] = $.04f\n", a);
00095   cvec_delete (vt);
00096   cvec_delete (v);
00097   vec_delete (a);
00098 
00099   /* -- complex vector -- */
00100   v = cvec_new_string ("1+0.5i 2-0.5i 3+0.25i 4-0.25i 5");
00101   it_printf ("v = $z\n", v);
00102 
00103   /* Transform the vector */
00104   vt = it_dft (v);
00105   cvec_delete (v);
00106   it_printf ("vt = $.04z\n", vt);
00107 
00108   /* Reconstruct the vector */
00109   v = it_idft (vt);
00110   cvec_delete (vt);
00111   it_printf ("v = $.04z\n", v);
00112 
00113   a = vec_new_1N (8);
00114   v = vec_to_cvec (a);
00115   it_printf ("v = $.04z\n", v);
00116 
00117   /* Transform the vector */
00118   vt = it_dft (v);
00119   cvec_delete (v);
00120   it_printf ("vt = $.04z\n", vt);
00121 
00122   /* Reconstruct the vector */
00123   v = it_idft (vt);
00124   cvec_delete (vt);
00125   it_printf ("v = $.04z\n", v);
00126 
00127   /* Transform the vector */
00128   /*  vt = __it_dft_real_fast(a); */
00129   /*  it_printf("vt = $z\n", vt); */
00130 
00131   vec_delete (a);
00132 
00133 
00134   /* -- convolution and correlation by FFT -- */
00135   r1 = cvec_new_string ("1 1 1");
00136   r2 = cvec_new_string ("2 2 2");
00137   vt = cvec_fft_conv (r1, r2);
00138 
00139   it_printf ("fft_conv($z,$z) = $z\n", r1, r2, vt);
00140   cvec_delete (vt);
00141   vt = cvec_fft_corr (r1, r2);
00142   it_printf ("fft_corr($z,$z) = $z\n", r1, r2, vt);
00143   cvec_delete (vt);
00144   vt = cvec_fft_autocorr (r1);
00145   it_printf ("fft_autocorr($z) = $z\n", r1, vt);
00146   cvec_delete (vt);
00147 
00148   return 0;
00149 }

Hosted by
Copyright (C) 2005-2006 Hervé Jégou
Vivien Chappelier
Francois Cayre
libit logo courtesy of Jonathan Delhumeau