ゲルツェルアルゴリズム

Goertzelアルゴリズムは、離散フーリエ変換（DFT）の各項を効率的に評価するためのデジタル信号処理（DSP）の手法です。従来のアナログ電話のキーパッドのプッシュボタンによって生成されるデュアルトーン多周波信号（DTMF）トーンの認識など、特定の実用的なアプリケーションで役立ちます。このアルゴリズムは、1958年にGerald Goertzelによって初めて記述されました。 ^{[ 1 ]}

DFT と同様に、Goertzel アルゴリズムは離散信号から選択可能な 1 つの周波数成分を分析します。^{[ 2 ] [ 3 ] [ 4 ]} 直接 DFT 計算とは異なり、Goertzel アルゴリズムは、実数値入力シーケンスに対して実数値演算を使用して、各反復で単一の実数値係数を適用します。全スペクトルをカバーする場合 (係数が後続の計算で再利用される連続データストリームを使用する場合を除く。この場合の計算複雑度はスライディング DFTと同等です)、Goertzel アルゴリズムは高速フーリエ変換(FFT) アルゴリズムよりも複雑度の次数が高いですが、少数の選択周波数成分を計算する場合は、数値的に効率的です。Goertzel アルゴリズムのシンプルな構造により、小型プロセッサや組み込みアプリケーションに適しています。

Goertzelアルゴリズムは、正弦波合成関数として「逆向きに」使用することもでき、生成されるサンプルごとに1回の乗算と1回の減算のみが必要です。^{[ 5 ]}

ゲルツェルアルゴリズムの主要な計算はデジタルフィルタの形式をとるため、このアルゴリズムはしばしばゲルツェルフィルタと呼ばれます。フィルタは入力シーケンスに対して動作します。${\displaystyle x[n]}$パラメータを持つ2段階のカスケード$\displaystyle \omega _{0}}$分析対象となる周波数を、サンプルあたりのラジアンに正規化して示す。

最初の段階では中間シーケンスを計算します。${\displaystyle s[n]}$:

第２段階では、次のフィルタを適用します。${\displaystyle s[n]}$出力シーケンスを生成する${\displaystyle y[n]}$:

最初のフィルタ段は、直接形式構造を持つ2次IIRフィルタであることが観察されます。この特定の構造は、内部状態変数がその段からの過去の出力値と等しいという特性を持っています。入力値${\displaystyle x[n]}$のために${\displaystyle n<0}$すべて0とみなされます。初期フィルタ状態を確立して、サンプルから評価を開始できるようにします。${\displaystyle x[0]}$フィルタ状態には初期値が割り当てられます${\displaystyle s[-2]=s[-1]=0}$エイリアシングの危険性を避けるため、周波数$\displaystyle \omega _{0}}$はしばしば 0 から π の範囲に制限されます (ナイキスト・シャノン標本化定理を参照)。この範囲外の値を使用することは無意味ではなく、この範囲内のエイリアシングされた周波数を使用することと同等です。指数関数は2π の周期で周期的であるためです。$\displaystyle \omega _{0}}$。

第2段階のフィルタは、過去の出力を一切使用しない計算を行うため、FIRフィルタであることがわかる。

Z変換法はフィルタカスケードの特性を調べるために適用できる。式(1)で与えられる第1フィルタ段のZ変換は次のようになる。

式（2）で示される第2フィルタ段のZ変換は

2つのフィルタ段のカスケード接続の合成伝達関数は次のようになる。

これは等価な時間領域シーケンスに変換でき、項はインデックスの最初の入力項まで展開されます。${\displaystyle n=0}$:

フィルタのZ変換の極は、${\displaystyle e^{+j\omega _{0}}}$そして${\displaystyle e^{-j\omega _{0}}}$複素Z変換平面の原点を中心とする単位半径の円上。この特性は、フィルタ処理が不安定であり、低精度の演算と長い入力シーケンスを使用して計算すると数値誤差の蓄積に対して脆弱であることを示している。 ^{[ 6 ]}数値的に安定したバージョンはChristian Reinschによって提案された。^{[ 7 ]}

DFT項を計算するという重要なケースにおいては、以下の特別な制約が適用される。

• フィルタリングはインデックスで終了します${\displaystyle n=N}$、 どこ${\displaystyle N}$は、DFTの入力シーケンスにおける項数です。
• ゲルツェル分析に選択された周波数は、特別な形式に限定される。

• インデックス番号${\displaystyle k}$DFTの「周波数ビン」がインデックス番号のセットから選択されていることを示す。

これらの代入を式（6）に代入し、項が${\displaystyle e^{+j2\pi k}=1}$すると、式（6）は次の形式になる。

式(9)の右辺は、DFT項の定義式と非常によく似ていることがわかります。${\displaystyle X[k]}$インデックス番号のDFT用語${\displaystyle k}$だが、全く同じではない。式（9）に示す総和は、${\displaystyle N+1}$入力用語のみ、ただし${\displaystyle N}$DFTを評価する際に、入力項が利用可能となる。単純だが洗練されていない方法として、入力シーケンスを拡張する方法がある。${\displaystyle x[n]}$もう一つ人工的な値${\displaystyle x[N]=0}$^{[ 8 ]} 式(9)から、最終結果に対する数学的な影響は項を削除するのと同じであることがわかる。${\displaystyle x[N]}$合計から得られるため、意図したDFT値が得られます。

しかし、余分なフィルタパスを回避する、より洗練されたアプローチがあります。式（1）から、拡張入力項の場合、${\displaystyle x[N]=0}$最終段階で使用されます。

したがって、アルゴリズムは次のように完成させることができる。

• 入力項の処理後、IIRフィルタを終了する${\displaystyle x[N-1]}$、
• 式(10)を適用して構築する${\displaystyle s[N]}$以前の出力から${\displaystyle s[N-1]}$そして${\displaystyle s[N-2]}$、
• 計算された式(2)を適用する${\displaystyle s[N]}$価値と${\displaystyle s[N-1]}$フィルターの最終的な直接計算によって生成される。

最後の2つの数学演算は、代数的に組み合わせることで簡略化されます。

フィルターの更新を終了時に停止することに注意してください${\displaystyle N-1}$また、式（11）ではなく式（2）をすぐに適用すると、最終的なフィルタ状態の更新が見逃され、位相が誤った結果が生じる。^{[ 9 ]}

Goertzelアルゴリズムに選択された特定のフィルタリング構造が、効率的なDFT計算の鍵となります。出力値は1つだけであることがわかります。${\displaystyle y[N]}$DFTの計算に使用されるため、他のすべての出力項の計算は省略されます。FIRフィルタは計算されないため、IIRステージの計算は省略されます。${\displaystyle s[0],s[1]}$などについては、第1段階の内部状態を更新した直後に破棄できます。

これは一見矛盾しているように思える。アルゴリズムを完了するには、FIRフィルタ段をIIRフィルタ段の最後の2つの出力を使用して一度評価する必要があるが、計算効率のためにIIRフィルタの反復処理ではその出力値が破棄される。ここで直接形式フィルタ構造の特性が適用される。IIRフィルタの2つの内部状態変数は、IIRフィルタ出力の最後の2つの値を提供し、これらはFIRフィルタ段を評価するために必要な項である。

• 計算複雑性理論によれば、一連の${\displaystyle M}$DFT 項を使用する${\displaystyle M}$データセットに対するGoertzelアルゴリズムの適用${\displaystyle N}$「1回の操作あたりのコスト」が${\displaystyle K}$複雑性がある${\displaystyle O(KNM)}$。

単一のDFTビンを計算するには${\displaystyle X(f)}$長さの複雑な入力シーケンスの場合${\displaystyle N}$Goertzelアルゴリズムでは、${\displaystyle 2N}$乗算と${\displaystyle 4\ N}$ループ内の加算/減算、4回の乗算、および4回の最終的な加算/減算の合計${\displaystyle 2N+4}$乗算と${\displaystyle 4N+4}$加算/減算。これは、${\displaystyle M}$周波数。

• 対照的に、FFTをデータセットに適用すると、${\displaystyle N}$値には複雑さがある${\displaystyle O(KN\log _{2}(N))}$。

これは使用するFFTアルゴリズムに依存するため直接適用するのは難しいが、典型的な例として基数2のFFTがあり、${\displaystyle 2\log _{2}(N)}$乗算と${\displaystyle 3\log _{2}(N)}$DFTビンごとの加算/減算、各${\displaystyle N}$ゴミ箱。

複雑度オーダー式では、計算された項の数${\displaystyle M}$より小さい${\displaystyle \log N}$ゴーツェルアルゴリズムの利点は明らかです。しかし、FFTコードは比較的複雑であるため、「単位作業あたりのコスト」要因が${\displaystyle K}$FFT の場合、多くの場合、より大きくなり、実用上の利点は、Goertzel アルゴリズムが有利になります。${\displaystyle M}$数倍大きい${\displaystyle \log _{2}(N)}$。

基数2のFFTとGoertzelアルゴリズムのどちらが効率的かを判断する目安として、項数を調整してください。${\displaystyle N}$データセット内で最も近い正確な 2 乗に上向きにし、これを${\displaystyle N_{2}}$、そして、Goertzel アルゴリズムは、以下の場合に高速になる可能性が高い。

${\displaystyle M\leq {\frac {5N_{2}}{6N}}\log _{2}(N_{2})}$

FFTの実装と処理プラットフォームは、相対的なパフォーマンスに大きな影響を与えます。一部のFFT実装^{[ 11 ]}は、内部で複素数計算を実行して係数をオンザフライで生成するため、「単位作業あたりのコストK」が大幅に増加します。FFTおよびDFTアルゴリズムは、数値効率を向上させるために事前に計算された係数値のテーブルを使用できますが、これには外部メモリにバッファリングされた係数値へのアクセスが増えるため、数値的な利点の一部を打ち消すキャッシュ競合の増加につながる可能性があります。

どちらのアルゴリズムも、入力データに複素数値ではなく実数値を使用すると、効率が約2倍向上します。ただし、この向上はGoertzelアルゴリズムでは当然のことながら、FFTでは実数値データの変換に特化した特定のアルゴリズム変種を使用しない限り実現できません。

• ブルースタインのFFTアルゴリズム（チャープZ）
• 周波数偏移変調（FSK）
• 位相偏移変調（PSK）

• Proakis, JG; Manolakis, DG (1996), Digital Signal Processing: Principles, Algorithms, and Applications , Upper Saddle River, NJ: Prentice Hall, pp. 480–481 , Bibcode : 1996dspp.book.....P 

• Wayback Machineに保存されているGoertzelアルゴリズム（2018年6月28日アーカイブ）
• 周波数解析のためのDSPアルゴリズム
• ゲルツェル アルゴリズム (Kevin Banks 著)
• ウーヴェ・バイスによるゲルツェルアルゴリズムの分析では、アナログの2次チェビシェフ低域通過フィルタと比較している。