/[marvin]/marvin/src/libsip/algo/fft.cc
ViewVC logotype

Diff of /marvin/src/libsip/algo/fft.cc

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.5 by exos, Fri Sep 12 01:50:01 2003 UTC revision 1.6 by exos, Fri Sep 12 03:06:18 2003 UTC
# Line 40  Line 40 
40  *******************************************************************/  *******************************************************************/
41  #include "fft.hh"  #include "fft.hh"
42    
43    static unsigned get_puiss(int n)
44    {
45      int res = 1;
46      while (res < n)
47        res *= 2;
48      return res;
49    }
50    
51  fft::fft ( DiscSignal &signal, bool inv){  fft::fft ( DiscSignal &signal, bool inv){
52    
53    double angle = -2.0 * PI;    double angle = -2.0 * PI;
# Line 47  fft::fft ( DiscSignal &signal, bool inv) Line 55  fft::fft ( DiscSignal &signal, bool inv)
55    unsigned ns = signal.length();    unsigned ns = signal.length();
56    double *result_r_tmp = (double *)malloc(sizeof(double) * ns);    double *result_r_tmp = (double *)malloc(sizeof(double) * ns);
57    double *result_i_tmp = (double *)malloc(sizeof(double) * ns);    double *result_i_tmp = (double *)malloc(sizeof(double) * ns);
58      DiscSignal &tmpsig = signal;
59      
60    if ( !puiss_de_2(ns) )    if ( !puiss_de_2(ns) )
61      {                                    {
62        std::cerr << "Not a two power length: "<<ns << std::endl;        unsigned l = get_puiss(ns);
63        exit(1);        double *tmp = (double *)malloc(sizeof(double) * ns);
64          int i;
65          for (i = 0; i < ns; i++)
66            tmp[i] = signal[i];
67          while (i < l)
68            tmp[i++] = 0;
69          tmpsig = *new DiscSignal(tmp, l);
70      }      }
   
71    if (inv)    if (inv)
72      angle *= -1;      angle *= -1;
73        
# Line 66  fft::fft ( DiscSignal &signal, bool inv) Line 79  fft::fft ( DiscSignal &signal, bool inv)
79    for (i = 0; i < ns; i++ )    for (i = 0; i < ns; i++ )
80      {      {
81        j = ReverseBits (i, nbi);        j = ReverseBits (i, nbi);
82        result_r_tmp[j] = signal[i];        result_r_tmp[j] = tmpsig[i];
83        result_i_tmp[j] = 0.0;        result_i_tmp[j] = 0.0;
84        }        }
85        

Legend:
Removed from v.1.5  
changed lines
  Added in v.1.6

savannah-hackers-public@gnu.org
ViewVC Help
Powered by ViewVC 1.1.26