格策尔算法

格策尔算法或格兹尔算法( 英语:Goertzel algorithm )是数字信号处理的一种运算技巧，此运算技巧提供一个有效率的方式来估计部分区域的离散傅立叶转换，广泛的运用在数字电话中的的双音多频信号(每个拨号的数字键由两个频率的音所组成，一个低频，一个高频)，此算法在1958年被杰拉德 · 格策尔（英语：Gerald_Goertzel）所提出^[1]。

格策尔算法与离散傅立叶转换的相似处在于他们都可以分析某个特定频段的离散讯号^[2][3][4]；相反的，它们的不同处在于，格策尔算法每次迭代的运算都是使用实数的乘法。虽然说在全频域的计算上，格策尔算法会比其他的傅立叶转换快速算法的复杂度来的高，但是它能区段式的分析每个小区段的频率组成，因此可以编写成较简单的运算架构，实际应用在处理器内的数值计算会更有效率。

格策尔算法逆向操作生成出弦波，而这个过程只需花费一个乘法和一个加法运算^[5]。

算法

格策尔算法把离散傅立叶转换看成是一组滤波器，将输入的讯号与滤波器中的脉冲响应做卷积运算，求得滤波器的输出，即得到频率域其中一点的频率${\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)