소인수 FFT 알고리즘(Prime-factor FFT algorithm) 이라고도 칭하며, 발견자의 이름을 사용해 Good–Thomas algorithm (1958/1963)으로도 불리는 소인수 알고리즘(Prime-factor Algorithm, PFA) 은 중국인의 나머지 정리 (Chinese Reminder Theorem, CRT)를 활용한 고속 푸리에 변환 (FFT) 알고리즘이다. 이 알고리즘은 크기 N = N 1 N 2 인 신호를 2차원 이산 푸리에 변환 (DFT) N = N 1 X N 2 로 재표현(re-expression)한다. 다만, 여기서 N 1 과 N 2 는 서로소 (relatively prime), 즉 gcd(N 1 , N 2 )=1이어야 한다.
이 알고리즘의 장점은
회전자(twiddle factor)에 의한 곱셈 연산(multiplication)이 불필요하다. 따라서, 처리 속도가 빠르다. 변환의 in-place, in-order version을 구현하는 것이 가능하다. 즉, 변환된 값들은 메모리에서 입력값들을 바꾼다. 그리고 입력값들이 시간순서일때 출력은 주파수 순서이다. 단점은
앞서 언급한 것처럼 N 1 과 N 2 는 '서로 소'이어야한다. 중국인의 나머지 정리(CRT)를 사용하기때문에 더 복잡한 색인 재지정(re-indexing)을 수행해야한다.
Good-Thomas 알고리즘과 혼합기수 Cooley-Tukey 알고리즘을 이용한 데이터 수 N=1000의 FFT 처리 방식. 신호의 길이가 N 일 때 이산 푸리에 변환 (DFT)은
f k = ∑ n = 0 N − 1 x n e − 2 π i N n k = ∑ n = 0 N − 1 x n W N n k k = 0 , … , N − 2 , N − 1 ( 1 ) {\displaystyle f_{k}=\sum _{n=0}^{N-1}x_{n}e^{-{\frac {2\pi i}{N}}nk}=\sum _{n=0}^{N-1}x_{n}W_{N}^{nk}\qquad k=0,\dots ,N-2,N-1\qquad \qquad \qquad (1)} 신호의 길이가 N = r s 인 경우를 생각해 봅니다. r 과 s 는 서로소 인 두 정수이다(위에서 언급한 N 1 , N 2 에 해당). 입력과 출력을 재색인화(re-indexing)하고 이를 DFT 공식에 대입하여 크기가 r X s 인 2차원 DFT로 변환해 보자. 이와같이 일차원 DFT를 다차원으로 변환하는 것을 다중차원 사상(multi-dimensional mapping) 이라고 한다.
정수론의 이론을 이용하여 (식 1) DFT의 입력 쪽 색인 n 과 출력 쪽 색인 k 를 다음과 같이 나타낼 수 있다.
1) 굿스 매핑(Good's mapping)
n = ( n 1 , n 0 ) ≡ s n 0 + r n 1 ( mod N ) ( 0 ≤ n 0 ≤ r − 1 , 0 ≤ n 1 ≤ s − 1 ) ( 2 ) {\displaystyle n=(n_{1},n_{0})\equiv sn_{0}+rn_{1}({\bmod {N}})\qquad (0\leq n_{0}\leq r-1,\quad 0\leq n_{1}\leq s-1)\qquad \qquad \qquad (2)} 최대공약수 표현 정리에 의하면 각 n에 대해 위 식을 만족하는 정수 n 1 과 n 2 가 있고, 합동이론에 따라 위의 식을 만족하고 위의 범위에 속하는 정수 n 1 과 n 2 가 유일하다는 것을 쉽게 알 수 있다. 이를 루리타니안 매핑(Ruritanian mapping) 또는 굿스 매핑(Good's mapping)이라 한다.
2) CRT 매핑(CRT mapping)
k ≡ k 0 ( mod r ) ( 0 ≤ k 0 ≤ r − 1 ) {\displaystyle k\equiv k_{0}({\bmod {r}})\qquad (0\leq k_{0}\leq r-1)} k ≡ k 1 ( mod s ) ( 0 ≤ k 1 ≤ s − 1 ) ( 3 ) {\displaystyle k\equiv k_{1}({\bmod {s}})\qquad (0\leq k_{1}\leq s-1)\qquad \qquad \qquad (3)} 중국인의 나머지 정리 (CRT)에 따르면 모든 쌍 (k1 , k2 )에 대해 (식 3)의 두 방정식을 만족하는 0 과 N-1 사이에 고유한 k가 존재한다. 이를 CRT 매핑이라고 하며, 이를 수식으로 표현하면 다음과 같다.
k = ( k 1 , k 0 ) ≡ [ s s − 1 ( mod r ) k 0 + r r − 1 ( mod s ) k 1 ] mod N ( 4 ) {\displaystyle k=(k_{1},k_{0})\equiv [ss^{-1}({\bmod {r}})k_{0}+rr^{-1}({\bmod {s}})k_{1}]{\bmod {N}}\qquad \qquad \qquad (4)} n 과 k 를 위와 같이 (식 3)과 (식 4)로 나타낼 수 있는 이유는r 과 s 가 '서로 소'이기 때문이다. 위의 전개와는 반대로 입력 n 을 CRT 매핑하고 출력 k 를 Good's 매핑 할 수도 있다.
위의 두 식 (식 2)와 (식 4)를 (식 1)에 대입하면 식의 회전자(twiddle factor)는 다음 (식 5)와 같은 형태를 가지게 된다.
W r s [ s n 0 + r n 1 ] mod N [ s s − 1 ( mod r ) k 0 + r r − 1 ( mod s ) k 1 ] mod N {\displaystyle W_{rs}^{[sn_{0}+rn_{1}]{\bmod {N}}[ss^{-1}({\bmod {r}})k_{0}+rr^{-1}({\bmod {s}})k_{1}]{\bmod {N}}}} 회전자 자체가 mod N 의 성질을 가지고 있으므로 지수(exponent)에 나타나는 mod N 은 생략할 수 있다. 즉,
= W r s [ s n 0 + r n 1 ] [ s s − 1 ( mod r ) k 0 + r r − 1 ( mod s ) k 1 ] {\displaystyle =W_{rs}^{[sn_{0}+rn_{1}][ss^{-1}({\bmod {r}})k_{0}+rr^{-1}({\bmod {s}})k_{1}]}} = W r s s n 0 s s − 1 ( mod r ) k 0 + s n 0 r r − 1 ( mod s ) k 1 + r n 1 s s − 1 ( mod r ) k 0 + r n 1 r r − 1 ( mod s ) k 1 {\displaystyle =W_{rs}^{sn_{0}ss^{-1}({\bmod {r}})k_{0}+sn_{0}rr^{-1}({\bmod {s}})k_{1}+rn_{1}ss^{-1}({\bmod {r}})k_{0}+rn_{1}rr^{-1}({\bmod {s}})k_{1}}} 어떤 임의의 수 m에 대해 W r s m r s = 1 {\displaystyle W_{rs}^{mrs}=1}
= W r s s n 0 s s − 1 ( mod r ) k 0 + r n 1 r r − 1 ( mod s ) k 1 {\displaystyle =W_{rs}^{sn_{0}ss^{-1}({\bmod {r}})k_{0}+rn_{1}rr^{-1}({\bmod {s}})k_{1}}} = W r n 0 s s − 1 ( mod r ) k 0 W s n 1 r r − 1 ( mod s ) k 1 {\displaystyle =W_{r}^{n_{0}ss^{-1}({\bmod {r}})k_{0}}W_{s}^{n_{1}rr^{-1}({\bmod {s}})k_{1}}} r과 s는 '서로 소'이므로 각각의 역수 r − 1 , s − 1 {\displaystyle r^{-1},s^{-1}} s s − 1 ≡ 1 mod r , r r − 1 ≡ 1 mod s {\displaystyle ss^{-1}\equiv 1{\bmod {r}},\quad rr^{-1}\equiv 1{\bmod {s}}} W r s s − 1 mod r = W r 1 , W s r r − 1 mod s = W s 1 {\displaystyle W_{r}^{ss^{-1}{\bmod {r}}}=W_{r}^{1},\quad W_{s}^{rr^{-1}{\bmod {s}}}=W_{s}^{1}}
= W r n 0 k 0 W s n 1 k 1 {\displaystyle =W_{r}^{n_{0}k_{0}}W_{s}^{n_{1}k_{1}}} ( 5 ) {\displaystyle \qquad \qquad (5)}
f ( k 1 , k 0 ) = ∑ n 1 = 0 s − 1 ∑ n 0 = 0 r − 1 x ( n 1 , n 0 ) W r n 0 k 0 W s n 1 k 1 ( 6 ) {\displaystyle f(k_{1},k_{0})=\sum _{n_{1}=0}^{s-1}\sum _{n_{0}=0}^{r-1}x(n_{1},n_{0})W_{r}^{n_{0}k_{0}}W_{s}^{n_{1}k_{1}}\qquad \qquad (6)} 고속 푸리에 변환 위키 페이지에서 볼 수 있는 혼합 기수(mixed-radix) FFT의 다음(식 20)의 회전자 W N 1 N 2 k 1 n 2 {\displaystyle W_{N_{1}N_{2}}^{k_{1}n_{2}}}
H ( k 1 , k 2 ) = ∑ n 1 = 0 N 1 − 1 W N 1 N 2 k 1 n 2 [ ∑ n 2 = 0 N 2 − 1 h ( n 1 , n 2 ) W N 2 k 2 n 2 ] W N 1 k 1 n 1 {\displaystyle H(k_{1},k_{2})=\sum _{n_{1}=0}^{N_{1}-1}W_{N_{1}N_{2}}^{k_{1}n_{2}}{\Biggl [}\sum _{n_{2}=0}^{N_{2}-1}h(n_{1},n_{2})W_{N_{2}}^{k_{2}n_{2}}{\Biggr ]}W_{N_{1}}^{k_{1}n_{1}}} = ∑ n 1 = 0 N 1 − 1 W N 1 N 2 k 1 n 2 h ′ ( n 1 , n 2 ) W N 1 k 1 n 1 ( 20 ) {\displaystyle =\sum _{n_{1}=0}^{N_{1}-1}W_{N_{1}N_{2}}^{k_{1}n_{2}}h^{'}(n_{1},n_{2})W_{N_{1}}^{k_{1}n_{1}}\qquad \qquad (20)} r=3, s=4 인 신호의 소인수 FFT 를 적용해보자. 이 경우 (식 2)는
n = ( n 1 , n 0 ) ≡ 4 n 0 + 3 n 1 ( mod 1 2 ) ( 0 ≤ n 0 ≤ 2 , 0 ≤ n 1 ≤ 3 ) {\displaystyle n=(n_{1},n_{0})\equiv 4n_{0}+3n_{1}({\bmod {1}}2)\qquad (0\leq n_{0}\leq 2,\quad 0\leq n_{1}\leq 3)} 이며, (식 4)에서 r r − 1 ( mod s ) = 33 − 1 ( mod 4 ) = 9 , s s − 1 ( mod r ) = 44 − 1 ( mod 3 ) = 4 {\displaystyle rr^{-1}({\bmod {s}})=33^{-1}({\bmod {4}})=9,\qquad ss^{-1}({\bmod {r}})=44^{-1}({\bmod {3}})=4} k 는
k = ( k 1 , k 0 ) ≡ 4 k 0 + 9 k 1 ( mod 1 2 ) ( 0 ≤ k 0 ≤ 2 , 0 ≤ k 1 ≤ 3 ) {\displaystyle k=(k_{1},k_{0})\equiv 4k_{0}+9k_{1}({\bmod {1}}2)\qquad (0\leq k_{0}\leq 2,\quad 0\leq k_{1}\leq 3)} 이 된다. (식 6)은
f ( k 1 , k 0 ) = ∑ n 1 = 0 3 ∑ n 0 = 0 2 x ( n 1 , n 0 ) W 3 n 0 k 0 W 4 n 1 k 1 {\displaystyle f(k_{1},k_{0})=\sum _{n_{1}=0}^{3}\sum _{n_{0}=0}^{2}x(n_{1},n_{0})W_{3}^{n_{0}k_{0}}W_{4}^{n_{1}k_{1}}} 그럼 ( n 1 , n 0 ) {\displaystyle (n_{1},n_{0})} x ( n ) = x ( n 1 , n 0 ) {\displaystyle x(n)=x(n_{1},n_{0})} ( k 1 , k 0 ) {\displaystyle (k_{1},k_{0})} f ( k ) = f ( k 1 , k 0 ) {\displaystyle f(k)=f(k_{1},k_{0})}
1차원 신호 데이터 원본 n 0 1 2 3 4 5 6 7 8 9 10 11 xn x0 x1 x2 x3 x4 x5 x6 x7 x8 x9 x10 x11
Good's mapping x(n1 ,n0 ) n0 0 1 2 n1 0 x(0) x(4) x(8) 1 x(3) x(7) x(11) 2 x(6) x(10) x(2) 3 x(9) x(1) x(5)
CRT mapping f(k1 ,k0 ) k0 0 1 2 k1 0 f(0) f(4) f(8) 1 f(9) f(1) f(5) 2 f(6) f(10) f(2) 3 f(3) f(7) f(11)
![{\displaystyle k=(k_{1},k_{0})\equiv [ss^{-1}({\bmod {r}})k_{0}+rr^{-1}({\bmod {s}})k_{1}]{\bmod {N}}\qquad \qquad \qquad (4)}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a4c28ce9b031eb86e9d867fe2b11c0dde7d920fd)
![{\displaystyle W_{rs}^{[sn_{0}+rn_{1}]{\bmod {N}}[ss^{-1}({\bmod {r}})k_{0}+rr^{-1}({\bmod {s}})k_{1}]{\bmod {N}}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/92bb175b2f7f08876ae989568fe6d21632ca3a34)
![{\displaystyle =W_{rs}^{[sn_{0}+rn_{1}][ss^{-1}({\bmod {r}})k_{0}+rr^{-1}({\bmod {s}})k_{1}]}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/9d7abb1d59725e0b2050759db04facf39c315c78)
![{\displaystyle H(k_{1},k_{2})=\sum _{n_{1}=0}^{N_{1}-1}W_{N_{1}N_{2}}^{k_{1}n_{2}}{\Biggl [}\sum _{n_{2}=0}^{N_{2}-1}h(n_{1},n_{2})W_{N_{2}}^{k_{2}n_{2}}{\Biggr ]}W_{N_{1}}^{k_{1}n_{1}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/deb062a74573dbd00ab769731a517272067b07c6)
