Algorytm FFT Bruuna

Algorytm Bruuna jest algorytmem szybkiej transformaty Fouriera (FFT) opartym na niezwykłym rekurencyjnym wielomianowym podejściu do faktoryzacji, zaproponowanym dla potęg dwójki przez G. Bruuna w 1978 r. I uogólnionym na dowolne, nawet złożone rozmiary przez H. Murakamiego w 1996 r. Ponieważ jego operacje obejmują tylko rzeczywiste współczynniki aż do ostatniego etapu obliczeń, początkowo zaproponowano jako sposób wydajnego obliczania dyskretnej transformaty Fouriera (DFT) rzeczywistych danych. Algorytm Bruuna nie był jednak szeroko stosowany, ponieważ podejścia oparte na zwykłym algorytmie FFT Cooleya-Tukeya zostały z powodzeniem dostosowane do rzeczywistych danych z co najmniej taką samą wydajnością. Ponadto istnieją dowody na to, że algorytm Bruuna może być z natury mniej dokładny niż algorytm Cooleya-Tukeya w obliczu skończonej precyzji numerycznej (Storn, 1993).

Niemniej jednak algorytm Bruuna ilustruje alternatywną strukturę algorytmiczną, która może wyrażać zarówno siebie, jak i algorytm Cooleya-Tukeya, a tym samym zapewnia interesujące spojrzenie na FFT, które pozwala na łączenie dwóch algorytmów i inne uogólnienia.

Wielomianowe podejście do DFT

Przypomnijmy, że DFT jest zdefiniowany wzorem:

${\ Displaystyle X_ {k} = \ suma _ {n = 0} ^ {N-1} x_{n}e^{-{\frac {2\pi i}{N}}nk}\qquad k=0,\kropki ,N-1.}$

Dla wygody oznaczmy N pierwiastków jedności przez ω _N ⁿ ( n = 0, ..., N − 1):

${\ Displaystyle \ omega _ {N} ^ {n} = e ^ {- {\ Frac {2 \ pi ja} {N}} n}}$

i zdefiniuj wielomian x ( z ), którego współczynniki to x _n :

${\ Displaystyle x (z) = \ suma _ {n = 0} ^ {N-1} x_ {n} z ^ {n}.}$

DFT można zatem rozumieć jako redukcję tego wielomianu; czyli X _k jest określone przez:

${\ Displaystyle X_ {k} = x (\ omega _ {N} ^ {k}) = x (z) \ mod (z-\omega _{N}^{k})}$

gdzie mod oznacza operację reszty z wielomianu . Kluczem do szybkich algorytmów, takich jak Bruun czy Cooley-Tukey, jest fakt, że można wykonać ten zestaw N operacji resztowych w etapach rekurencyjnych.

Rozkłady rekurencyjne i FFT

$obliczyć DFT$ , musimy ocenić resztę wielomianów N stopnia-1, jak Obliczanie tych reszt jedna po drugiej jest równoznaczne z bezpośrednim obliczaniem zwykłej formuły DFT i wymaga operacji O( N ^{2 ).} Można jednak łączyć te reszty rekurencyjnie, aby zmniejszyć koszt $}$ następującą sztuczkę: jeśli chcemy oszacować dwa wielomiany $Displaystyle U (z$ i ${\ Displaystyle V (z)}$ , możemy najpierw wziąć resztę modulo ich iloczyn ${\ Displaystyle U (z)}$${\ Displaystyle V (z)}$ , co zmniejsza stopień wielomianu i $sprawia, że$ ​​​​kolejne modulo są mniej

$Iloczyn$ wszystkich jednomianów $z^{N}-1$ = k=0..N-1 is simply (whose roots are clearly the N roots of unity). One then wishes to find a recursive factorization of $z^{N}-1$ into polynomials of few terms and smaller and smaller degree. To compute the DFT, one takes $x(z)$ modulo each level of this factorization in turn, recursively, until one arrives at the monomials and the final result. If each level of the factorization splits every polynomial into an O(1) (constant-bounded) number of smaller polynomials, each with an O(1) number of nonzero coefficients, then the modulo operations for that level take O(N) time; since there will be a logarithmic number of levels, the overall complexity is O (N log N).

Bardziej wyraźnie, załóżmy na przykład, że ${\ Displaystyle z ^ {N} -1 = F_ {1} (z) F_ {2} ( z) F_ {3} (z)}$ i że ${\ Displaystyle F_ {k} (z) = F_ {k, 1} (z)F_{k,2}(z)}$ i tak dalej. Odpowiedni algorytm FFT składałby się najpierw z obliczenia x _k ( z ) = x ( z ) mod F _k ( z ), następnie obliczenia x _{k , j} ( z ) = x _k ( z ) mod F _{k , j} ( z ), i tak dalej, rekurencyjnie tworząc coraz więcej wielomianów resztowych o coraz mniejszym stopniu, aż dojdzie się do końcowych wyników stopnia 0.

Co więcej, dopóki czynniki wielomianowe na każdym etapie są względnie pierwsze (co w przypadku wielomianów oznacza, że ​​nie mają one wspólnych pierwiastków), można skonstruować algorytm dualny, odwracając proces z chińskim twierdzeniem o resztach .

Cooley-Tukey jako faktoryzacja wielomianu

Standardowy algorytm decymacji w częstotliwości (DIF) radix- r Cooley-Tukey odpowiada ściśle rekurencyjnej faktoryzacji. Na przykład radix-2 DIF Cooley-Tukey czynniki ${\ Displaystyle z ^ {N} -1}$ na ${\ Displaystyle F_ {1} = (z ^ { N/2} -1)}$ i ${\ Displaystyle F_ {2} = (z ^ {N/2} + 1)}$ . Te operacje modulo zmniejszają stopień o $rozmiaru$ problemu przez 2. Jednak zamiast rekurencyjnego rozkładu na czynniki bezpośrednio, $}$ Zamiast tego Cooley-Tukey najpierw oblicza x ₂ ( z ω _N ), przesuwając wszystkie pierwiastki (o $)$ twiddle , aby mógł zastosować rekurencyjną faktoryzację podproblemów. Oznacza to, że Cooley – Tukey zapewnia, że ​​wszystkie podproblemy są również DFT, podczas gdy generalnie nie jest to prawdą w przypadku arbitralnej faktoryzacji rekurencyjnej (takiej jak Bruun poniżej).

Faktoryzacja Bruuna

Podstawowy algorytm Bruuna dla potęg dwóch N = 2 ⁿ rozkłada na czynniki z ^{2 n} - 1 rekurencyjnie za pomocą reguł:

${\ Displaystyle z ^ {2M} -1 = (z ^ {M} -1) (z ^ {M} + 1) \,}$

${\ Displaystyle z ^ {4M} + az ^ {2M }+1=(z^{2M}+{\sqrt {2-a}}z^{M}+1)(z^{2M}-{\sqrt {2-a}}z^{M}+ 1)}$

gdzie a jest rzeczywistą stałą z | | _ ≤ 2. Jeśli za $\ Displaystyle a = 2 \ cos (\ phi)}$${\ Displaystyle \ phi \ in (0, \ pi)}$ , , to ${\ Displaystyle {\ sqrt {2 + a}} = 2 \ sałata {\ tfrac {\ phi} {2}}}$ i ${\ Displaystyle {\ sqrt {2-a}} = 2 \ cos ({\ tfrac {\ pi} {2}} - {\ tfrac {\ phi} {2}})$ }

Na etapie s , s = 0,1,2, n -1 stan pośredni składa się z 2 ^s wielomianów ${\ Displaystyle p_ {s, 0}, \ kropki, p_{s,2^{s}-1}}$ stopnia 2 ^{n - s} - 1 lub mniej , gdzie

${\ Displaystyle {\ rozpocząć {wyrównane} p_ {s, 0} (z) i = p (z) \ mod \ left(z^{2^{ns}}-1\right)&\quad &{\text{and}}\\p_{s,m}(z)&=p(z)\mod \left(z ^{2^{ns}}-2\cos \left({\tfrac {m}{2^{s}}}\pi \right)z^{2^{n-1-s}}+1\ po prawej)&m&=1,2,\kropki ,2^{s}-1\end{wyrównane}}}$

Dzięki konstrukcji rozkładu na czynniki z ^{2 n} - 1 , wielomiany p _{s , m} ( z ) każdy kodują 2 ^{n - s} wartości

${\ Displaystyle X_ {k} = p (e ^ {2 \ pi ja {\ tfrac {k} {2 ^ {n}}}}}}$

00 transformaty Fouriera, dla m =0 indeksy objęte to k = , 2 ^k , 2∙2 ^s , 3∙2 ^s ,…, (2 ^{n - s} -1)∙2 ^s , dla m > indeksy objęte są k = m , 2 ^{s +1} - m , 2 ^{s +1} + m , 2∙2 ^{s +1} - m , 2∙2 ^{s +1} + m , …, 2 ⁿ - m .

0 Podczas przejścia do następnego etapu wielomian $} s + 1, \ ell} (z)}$ wielomianów ${\ Displaystyle p_ {s,$ i ${\ Displaystyle p_ {s + 1,2 ^ {s} - \ ell} (z)}$ przez dzielenie wielomianowe. Jeśli ktoś chce zachować wielomiany w rosnącej kolejności indeksów, ten wzorzec wymaga implementacji z dwiema tablicami. Implementacja na miejscu tworzy przewidywalną, ale wysoce nieuporządkowaną sekwencję indeksów, na przykład dla N = 16 ostateczna kolejność 8 reszt liniowych to ( , 4 , 2 , 6 , 1 , 7 , 3 , 5 ).

₀ Na końcu rekurencji dla s = n - 1 pozostaje 2 ^{n - 1} wielomianów liniowych kodujących dwa współczynniki Fouriera X i X _{2 ^{n -1}} dla pierwszego i dla każdego k -tego wielomianu współczynniki X _k i X _{2 ⁿ - k} .

Na każdym etapie rekurencyjnym wszystkie wielomiany wspólnego stopnia 4M - 1 są redukowane do dwóch części połowy stopnia 2M - 1 . Dzielnikiem tego obliczenia reszty z wielomianu jest wielomian kwadratowy z ^m , tak że wszystkie redukcje można sprowadzić do dzielenia wielomianowego wielomianów sześciennych przez kwadratowe. Na każdym etapie jest N / 2 = 2 ^{n - 1} takich małych podziałów, co prowadzi do algorytmu O ( N log N ) dla FFT.

Co więcej, ponieważ wszystkie te wielomiany mają czysto rzeczywiste współczynniki (aż do ostatniego etapu), automatycznie wykorzystują specjalny przypadek, w którym dane wejściowe x _n są czysto rzeczywiste, aby zaoszczędzić mniej więcej dwukrotność obliczeń i przechowywania. Można również bezpośrednio wykorzystać przypadek rzeczywistych danych symetrycznych do obliczenia dyskretnej transformaty kosinusowej (Chen i Sorensen, 1992).

Uogólnienie na dowolne pierwiastki

Faktoryzacja Bruuna, a tym samym algorytm FFT Bruuna, został uogólniony do obsługi dowolnych parzystych długości złożonych, tj. dzielenia stopnia wielomianu przez dowolną podstawę (czynnik), w następujący sposób. Najpierw definiujemy zbiór wielomianów φ _{N ,α} ( z ) dla dodatnich liczb całkowitych N i dla α w [0,1) przez:

${\ Displaystyle \ phi _ {N, \ alfa} (z)=\left\{{\begin{matrix}z^{2N}-2\cos(2\pi \alpha )z^{N}+1&{\mbox{if }}0<\alpha <1 \\\\z^{2N}-1&{\mbox{if}}\alpha =0\end{macierz}}\right.}$

Zauważ, że wszystkie wielomiany, które pojawiają się w powyższym rozkładzie na czynniki Bruuna, można zapisać w tej formie. $_$ zerowe ${\ Displaystyle k = 0,1, \ kropki, N-1}$ dla $w przypadku i mi$ 2 π $^ {2 \$$}$ dla ${\ Displaystyle k = 0,1, \ kropki, 2N-1$ w przypadku Stąd te wielomiany można rekurencyjnie rozłożyć na czynniki dla czynnika (podstawy) r poprzez:

${\ Displaystyle \ phi _ {rM, \ alfa} (z) = \ lewo \ {{\ rozpocząć {tablica}} ll}\prod _{\ell =0}^{r-1}\phi _{M,(\alpha +\ell)/r}&{\mbox{if }}0<\alpha \leq 0.5\\ \\\prod _{\ell =0}^{r-1}\phi _{M,(1-\alpha +\ell )/r}&{\mbox{if }}0.5<\alpha <1\ \\\\prod _{\ell =0}^{r-1}\phi _{M,\ell /(2r)}&{\mbox{if}}\alpha =0\end{tablica}}\ Prawidłowy.}$

• Georg Bruun, „ Z -Transform DFT filtry i FFT”, IEEE Trans. w sprawie akustyki, przetwarzania mowy i sygnałów (ASSP) 26 (1), 56-63 (1978).
• HJ Nussbaumer, Szybka transformata Fouriera i algorytmy splotu (Springer-Verlag: Berlin, 1990).
• Yuhang Wu, „Nowe struktury FFT oparte na algorytmie Bruuna”, IEEE Trans. ASSP 38 (1), 188-191 (1990)
• Jianping Chen i Henrik Sorensen, „Wydajny algorytm FFT dla danych rzeczywistych-symetrycznych”, Proc. ICASSP 5 , 17-20 (1992).
• Rainer Storn, „Niektóre wyniki analizy błędów punktu stałego algorytmu Bruuna-FTT [ sic ]”, IEEE Trans. Proces sygnału. 41 (7), 2371-2375 (1993).
• Hideo Murakami, „Algorytmy dziesiątkowania w czasie i częstotliwości o wartościach rzeczywistych”, IEEE Trans. Układ obwodów II: Sygnał analogowy i cyfrowy. proc. 41 (12), 808-816 (1994).
• Hideo Murakami, „Szybka dyskretna transformata Fouriera o wartościach rzeczywistych i algorytmy splotów cyklicznych o wysoce złożonej parzystej długości”, Proc. ICASSP 3 , 1311-1314 (1996).
• Shashank Mittal, Md. Zafar Ali Khan, MB Srinivas, „Studium porównawcze różnych architektur FFT dla radia definiowanego programowo”, notatki z wykładów z informatyki 4599 ( wbudowane systemy komputerowe: architektury, modelowanie i symulacja ), 375-384 (2007 ). proc. 7. międzynarodowy Warsztaty, SAMOS 2007 (Samos, Grecja, 16–19 lipca 2007).