[PYTHON] Pseudo code FFT

Il s'agit d'un pseudo code de FFT qui résume simplement les facteurs de rotation de DFT. Puisqu'il s'agit d'un pseudo code ou Python, il peut être exécuté pour le moment, mais la conception FFT appropriée est décrite en détail sur la page de M. Oura (attention aux caractères déformés).

Radix-2 FFT

\begin{eqnarray}
A_{2k} &=& \sum_{j=0}^{N/2-1} (a_j + a_{N/2+j}) W_{N/2}^{jk} \\
A_{2k+1} &=& \sum_{j=0}^{N/2-1} ( a_j - a_{N/2+j}) W_N^j W_{N/2}^{jk} 
\end{eqnarray}

radix_2_fft.py


from cmath import exp, pi
from itertools import chain

def flatten(l):
    return list(chain.from_iterable(l))

def fft(x):
    n = len(x)
    if(n == 1):
        return x
    else:
        # Radix-2
        t1 = (exp(-2j*pi*k/n) for k in xrange(n/2))
        x0 = fft( [(a+b)   for a, b    in zip(x, x[n/2:])] )
        x1 = fft( [(a-b)*w for a, b, w in zip(x, x[n/2:], t1)] )
        return flatten(zip(x0, x1))

Mixed-Radix FFT (Radix-4 & Radix-2)

\begin{eqnarray}
A_{4k} &=& \sum_{j=0}^{N/4-1} (a_j + a_{N/4+j} + a_{2N/4+j} + a_{3N/4+j}) W_{N/4}^{jk} \\
A_{4k+1} &=& \sum_{j=0}^{N/4-1} (a_j - a_{2N/4+j} - i(a_{N/4+j} - a_{3N/4+j}) ) W_N^j W_{N/4}^{jk} \\
A_{4k+2} &=& \sum_{j=0}^{N/4-1} (a_j - a_{N/4+j} + a_{2N/4+j} - a_{3N/4+j}) W_N^{2j} W_{N/4}^{jk} \\
A_{4k+3} &=& \sum_{j=0}^{N/4-1} (a_j - a_{2N/4+j} + i(a_{N/4+j} - a_{3N/4+j}) ) W_N^{3j} W_{N/4}^{jk}
\end{eqnarray}

mixed_radix_fft.py


from cmath import exp, pi
from itertools import chain

def flatten(l):
    return list(chain.from_iterable(l))

def fft(x):
    n = len(x)
    if(n == 1):
        return x
    elif(n == 2):
        # Radix-2
        return [x[0]+x[1], x[0]-x[1]]
    else:
        # Radix-4
        t1 = (exp(-2j*pi*k/n)   for k in xrange(n/4))
        t2 = (exp(-2j*pi*2*k/n) for k in xrange(n/4))
        t3 = (exp(-2j*pi*3*k/n) for k in xrange(n/4))
        x0 = fft( [(a+b+c+d)        for a, b, c, d    in zip(x, x[n/4:], x[n*2/4:], x[n*3/4:])] )
        x1 = fft( [(a-c-1j*(b-d))*w for a, b, c, d, w in zip(x, x[n/4:], x[n*2/4:], x[n*3/4:], t1)] )
        x2 = fft( [(a-b+c-d)*w      for a, b, c, d, w in zip(x, x[n/4:], x[n*2/4:], x[n*3/4:], t2)] )
        x3 = fft( [(a-c+1j*(b-d))*w for a, b, c, d, w in zip(x, x[n/4:], x[n*2/4:], x[n*3/4:], t3)] )
        return flatten(zip(x0, x1, x2, x3))

Split-Radix FFT

\begin{eqnarray}
A_{2k} &=& \sum_{j=0}^{N/2-1} (a_j + a_{N/2+j}) W_{N/2}^{jk} \\
A_{4k+1} &=& \sum_{j=0}^{N/4-1} (a_j - a_{2N/4+j} - i(a_{N/4+j} - a_{3N/4+j}) ) W_N^j W_{N/4}^{jk} \\
A_{4k+3} &=& \sum_{j=0}^{N/4-1} (a_j - a_{2N/4+j} + i(a_{N/4+j} - a_{3N/4+j}) ) W_N^{3j} W_{N/4}^{jk}
\end{eqnarray}

split_radix_fft.py


from cmath import exp, pi
from itertools import chain

def flatten(l):
    return list(chain.from_iterable(l))

def fft(x):
    n = len(x)
    if(n == 1):
        return x
    elif(n == 2):
        # Radix-2
        return [x[0]+x[1], x[0]-x[1]]
    else:
        # Split-Radix
        t1 = (exp(-2j*pi*k/n) for k in xrange(n/4))
        t2 = (exp(-2j*pi*3*k/n) for k in xrange(n/4))
        x0 = fft( [(a+b)            for a, b in zip(x, x[n/2:])] )
        x1 = fft( [(a-c-1j*(b-d))*w for a, b, c, d, w in zip(x, x[n/4:], x[n*2/4:], x[n*3/4:], t1)] )
        x2 = fft( [(a-c+1j*(b-d))*w for a, b, c, d, w in zip(x, x[n/4:], x[n*2/4:], x[n*3/4:], t2)] )
        return flatten(zip(x0[0::2], x1, x0[1::2], x2))

Recommended Posts

Pseudo code FFT
Code de caractère
Mauvais code