格策尔算法

格策尔算法把离散傅立叶转换看成是一组滤波器，将输入的讯号与滤波器中的脉冲响应做卷积运算，求得滤波器的输出，即得到频率域其中一点的频率${\displaystyle X(k)}$ 。此算法利用旋转因子${\displaystyle {\omega }_{N}^{k}}$ 的周期性，将离散傅立叶转换转换为线性的滤波运算。

因为旋转因子

${\displaystyle {\omega }_{N}^{-kN}=e^{-j(2{\pi /N})(-kN)}=1}$

(1)

可得转换后第k点的频率为

${\displaystyle X(k)={\omega }_{N}^{-kN}\sum _{m=0}^{N-1}x(m){\omega }_{N}^{km}=\sum _{m=0}^{N-1}x(m){\omega }_{N}^{-k(N-m)}\qquad \quad ,k=0,1,2,...,N-1}$

(2)

定义${\displaystyle y_{k}(n)}$ 为

${\displaystyle y_{k}(n)=\sum _{m=0}^{N-1}x(m){\omega }_{N}^{-k(n-m)}}$

(3.A)

可将${\displaystyle y_{k}(n)}$ 理解为由两个讯号的卷积运算得出的结果

${\displaystyle y_{k}(n)=x(n)\otimes h_{k}(n)}$

(3.B)

其中${\displaystyle x(n)}$ 式输入的N点讯号，另外一个${\displaystyle h_{k}(n)}$ 则被看作是IIR滤波器的脉冲频率响应

${\displaystyle h_{k}(n)={\omega }_{N}^{-kn}\ u(n)}$

(4)

对比(2)和(3)式，可推知(3.A)进行卷积运算，当n=N时，滤波器的输出${\displaystyle y_{k}(N)}$ 即为${\displaystyle X(k)}$ :

${\displaystyle X(k)=y_{k}(n)\lfloor _{n=N}}$

(5)

对(4)进行Z转换，可得一阶IIR转移函数

${\displaystyle H_{k}(z)={\frac {1}{1-{{\omega }^{-k}\ z^{-1}}}}}$

(6)

图一为此系统的流程图，其对应的差分方程式为:

${\displaystyle y_{k}(n)={\omega }_{N}^{-k}\ y_{k}(n-1)+x(n),\qquad \ y(-1)=0}$

(7)

依照此差分方程进行迭代运算，迭代到${\displaystyle n=N}$ 时即可依据(5)式得到${\displaystyle X(k)}$ 。而依照转移函数(6)式进行运算时，可以先将旋转因子${\displaystyle {\omega }_{N}^{k}}$ 储存起来，每次迭代包含一次复数乘法，则按照(1)式计算N点离散傅立叶转换时则需要${\displaystyle 4N^{2}}$ 次实数乘法运算和${\displaystyle N(4N-2)}$ 次加/减法^[6]，加/减法与乘法运算皆为${\displaystyle 4N^{2}}$ 次，当N不大时运算效率不佳，若改为接下来改进的的格策尔算法(二阶)，所需的实数乘法次数约为原本的一半。

将式(6)上下同乘以${\displaystyle 1-{\omega }_{N}^{k}\ z^{-1}}$ ，可得第k点的频率响应转移函数为

${\displaystyle {\begin{alignedat}{2}H_{k}(z)&={\frac {1-{{\omega }_{N}^{k}\ z^{-1}}}{(1-{{\omega }_{N}^{-k}\ z^{-1}})(1-{{\omega }_{N}^{k}\ z^{-1})}}}\\&={\frac {1-{{\omega }_{N}^{k}\ z^{-1}}}{1-2\ \cos((2{\pi }/N)k)z^{-1}+z^{-2})}}\\\end{alignedat}}}$

(8)

此转移函数所对应的系统流程图如图二所示，复数分析(8)式，可得知此二阶滤波器有一对共轭的极点与一个零点。图二中在计算${\displaystyle x(n)}$ 的转换结果${\displaystyle X(k)}$ 时，会有两个步骤:

1. 共轭极点迭代计算 依序将输入讯号${\displaystyle x(0),x(1),x(2),...,x(n-1)}$ 放入滤波器做迭代运算，共作N次迭代，计算量是${\displaystyle 2N}$ 次实数乘法与${\displaystyle 4N}$ 次实数加/减法
2. 零点迭代计算 输入讯号${\displaystyle x(n)}$ 是N点的讯号从${\displaystyle n=0,1,2,3,...,N-1}$ 。加入${\displaystyle x(N)=0}$ 的边界条件，可以按照图二的流程图计算出${\displaystyle y_{k}(N)}$ ，此即为所求的${\displaystyle x(n)}$ 离散傅立叶转换${\displaystyle X(k)}$ ，此步骤的计算量为4次实数乘法与4次实数加/减法。

综合以上步骤，总共的计算量为${\displaystyle 2N+4}$ 次实数乘法运算以及${\displaystyle 4N+4}$ 次实数加法运算，而使用此计算算法只需储存${\displaystyle {\omega }_{N}^{k}}$ 与${\displaystyle \cos((2{\pi }/N)k)}$ 两个参数^[7]。

• 双音多频
• Chirp-Z转换
• 频率偏移调变(FSK)
• 相位偏移调变(PSK)

1. ^ Goertzel, G. (January 1958), "An Algorithm for the Evaluation of Finite Trigonometric Series", American Mathematical Monthly, 65 (1): 34–35, doi:10.2307/2310304
2. ^ Mock, P. (March 21, 1985), "Add DTMF Generation and Decoding to DSP-μP Designs （页面存档备份，存于互联网档案馆）" (PDF), EDN, ISSN 0012-7515 （页面存档备份，存于互联网档案馆）; also found in DSP Applications with the TMS320 Family, Vol. 1, Texas Instruments, 1989.
3. ^ Chen, Chiouguey J. (June 1996), Modified Goertzel Algorithm in DTMF Detection Using the TMS320C80 DSP (PDF), Application Report, Texas Instruments, SPRA066
4. ^ Schmer, Gunter (May 2000), DTMF Tone Generation and Detection: An Implementation Using the TMS320C54x （页面存档备份，存于互联网档案馆） (PDF), Application Report, Texas Instruments, SPRA096a
5. ^ http://haskell.cs.yale.edu/wp-content/uploads/2011/01/AudioProc-TR.pdf.
6. ^ http://www.docin.com/p-577391532.html
7. ^ 格策爾介紹網站(英文). [2017-06-22]. （原始内容存档于2017-06-22）. 

• Proakis, J. G.; Manolakis, D. G. (1996), Digital Signal Processing: Principles, Algorithms, and Applications, Upper Saddle River, NJ: Prentice Hall, pp. 480–481

• 格策尔算法网站介绍(英文)
• 频域分析数字信号数理算法 （页面存档备份，存于互联网档案馆）(英文)
• 格策尔算法网站介绍 （页面存档备份，存于互联网档案馆）(英文) 作者: Kevin Bank,2002