lib/fft.lib

   1 /******************************************************************
   2  *  FFT
   3  *  Implementation contributed by Remy Muller
   4  *****************************************************************/
   5
   6
   7 // twiddle_mult(n) : n parallel cables
   8
   9 W(k, n) = 1, (0, ( k, ((2 * PI) / n) : *) : -) : polar_cplx;
  10
  11 twiddle_mult(k, n) = _, W(k, n) : pcplx_mul;
  12
  13 // butterfly(n) : addition then substraction of interleaved signals :
  14 xbutterfly(n) = (bus(n/2), par(k, n/2, twiddle_mult(k, n))) <: interleave(n/2,2), interleave(n/2,2) : par(i, n/2, pcplx_add), par(i, n/2, pcplx_sub);
  15
  16 //btf_downside(n) = bus(n) : interleave(n/2,2);
  17
  18 // fft(n) : fft matrix function of size n = 2^k
  19 //fft(2) = butterfly(2);
  20 //fft(n) = butterfly(n) : (fft(n/2) , fft(n/2));
  21
  22 xbutterflies(2) = xbutterfly(2);
  23 xbutterflies(n) = (xbutterflies(n/2) , xbutterflies(n/2)) : xbutterfly(n);
  24
  25
  26 evens = case {
  27         (2) => _ , ! ;
  28         (n) => _ , ! , evens(n - 2);
  29 };
  30 odds = case {
  31         (2) => ! , _ ;
  32         (n) => ! , _ , odds(n - 2);
  33 };
  34
  35 eo(n) = evens(n), odds(n);
  36
  37 shuffling = case {
  38         (2) => eo(2);
  39         (n) => (evens(n) <: shuffling(n/2)), (odds(n) <: shuffling(n/2));
  40 };
  41
  42 shuffle(n) = bus(n) <: shuffling(n);
  43
  44 real2pcplx(n) = par(i, n, (sca2pcplx));
  45
  46 //fft(n) =  shuffle(n) : xbutterflies(n);
  47 fft(n) = _ <: picks(n) : real2pcplx(n) : shuffle(n) : xbutterflies(n);
  48 fftc(n) = _ <: picks(n) : shuffle(n) : xbutterflies(n) : pcplx_moduls(n); // already complex input
  49
  50
  51 picks(n) = par(i, n, [i]);
  52
  53 concats = case {
  54         (1) => vectorize(1);
  55         (n) => concats(n-1) # vectorize(1);
  56 };
  57
  58 nconcat(n) = concats(n); //fake name for svg block encapsulation
  59
  60 pack = case {
  61         (1) => _;
  62         (n) => pack(n-1), _ : #;
  63 };
  64
  65 delays(m) = _ <: par(i, m, @(i));
  66
  67 overlap(n,m) = vectorize(n/m) : delays(m) : pack(m);
  68
  69 stops(n) = par(i, n, !);