JP2016090670A

JP2016090670A - 周波数分析装置およびプログラム

Info

Publication number: JP2016090670A
Application number: JP2014221867A
Authority: JP
Inventors: 陽前澤; Akira Maezawa
Original assignee: Yamaha Corp
Current assignee: Yamaha Corp
Priority date: 2014-10-30
Filing date: 2014-10-30
Publication date: 2016-05-23

Abstract

【課題】音響信号などの時系列ディジタルデータの周波数分析を効率よく短時間で実行できる周波数分析装置およびプログラムを提供することを目的とする。
【解決手段】入力した時系列ディジタルデータに対し所定のホップ長Ｈ毎に所定のフレーム長Ｎの長さのフレームの周波数分析を行う周波数分析装置であって、既に分析済みのｍ−１番目のフレームの周波数分析結果に所定の位相差を加えることによりｍ番目のフレームの分析結果の第１項を取得し、ｍ−１番目のフレームのサンプルデータとｍ番目のフレームの新たなサンプルデータとの差分に対して所定の離散フーリエ変換の処理を施してｍ番目のフレームの分析結果の第２項を取得し、前記第１項に第２項を加えることにより今回の分析対象であるｍ番目のフレームの周波数分析結果を取得する。
【選択図】図２

Description

この発明は、入力した時系列ディジタルデータの周波数成分を分析する周波数分析装置およびプログラムに関する。

従来より、音響信号などの時系列ディジタルデータの周波数分析等を行う手法として離散フーリエ変換（Discrete Fourier Transform：ＤＦＴ）が知られている。通常は無限の区間に渡って積分を実行するわけにはいかないので、分析対象のディジタル信号から、所定のホップ長毎に所定のフレーム長の区間のデータを切り出し窓関数を掛けて離散フーリエ変換する処理を繰り返して分析する、いわゆる短時間フーリエ変換（Short-Time Fourier Transform：ＳＴＦＴ）の手法が使われる。

図７は、ＳＴＦＴの概念図を示す。７０１は時系列データであるディジタル信号ｘ（ｔ）の時間変化を示すグラフである。７０２は最初に切り出したフレーム長がＮの区間のフレーム（データ列）を示す。このフレーム７０２に対してＤＦＴの演算を実行し、最初の分析結果X₀(f)を得る。X₀(f)の下付の添え字は分析したフレームを特定するフレーム番号（フレーム０を先頭とする）、パラメータｆは周波数ビンを特定する番号（あるいは周波数）を示す。分析結果X₀(f)は、各周波数ビンにおける分析対象の信号の成分を複素数表現したデータ列である。７０３は次のフレーム１を示す。フレーム１は、フレーム０からホップ長Ｈだけ進んだ位置からフレーム長Ｎの区間を切り出したものである。フレーム１に対してＤＦＴの演算を実行し分析結果X₁(f)を得る。同様にして、ホップ長Ｈ毎にフレーム長Ｎの区間をフレームとして順次切り出しＤＦＴの演算を実行していく。なお、ＤＦＴの計算においては、演算を高速に行うことができる高速フーリエ変換（Fast Fourier Transform：ＦＦＴ）の手法を用いることが多い。

さらに、ＤＦＴの循環性を用いることによるスライディング離散フーリエ変換（Sliding-DFT：ＳＤＦＴ）という手法も提案されている（非特許文献１，２参照）。図８は、ＳＤＦＴの概念図を示す。８０１は、７０１と同様の信号ｘ（ｔ）のグラフを示す。８０２は最初に切り出したフレーム長がＮの区間のフレームを示す。このフレーム８０２に対してＤＦＴの演算を実行し、最初の分析結果X₀(f)を得る。次に、フレーム８０２の処理対象データを１サンプル分スライドする。すなわち、先頭サンプル８０３を削除するとともにフレーム８０２の最終サンプルの次のサンプル８０５を追加し、この網掛け８０４と８０５に示す一連のデータの周波数分析結果を求める。この計算は、フレーム８０２の分析結果X₀(f)を補正することで得ることができる。この補正の計算は各周波数ビン毎の計算であるが、ＤＦＴの計算を行うのに比べれば簡単で計算量が少ない。以下同様にして、８０６，８０７に示すように分析対象データをスライドさせていき、補正計算で各周波数ビンの分析結果を求めていく。ホップ長Ｈまでスライドさせることで、次のフレーム８０８の分析結果が得られる。

IEEE SIGNAL PROCESSING MAGAZINE 2003年3月第74-80頁、dsp tips & tricks "The Sliding DFT"、Eric Jacobsen and Richard Lyons IEEE SIGNAL PROCESSING MAGAZINE 2004年1月第110-111頁、dsp tips & tricks "An Update to the Sliding DFT"、Eric Jacobsen and Richard Lyons

音響信号の周波数分析等に適用する短時間フーリエ変換の計算は、できるだけ効率よく短時間で行うことが求められる。

一般にＦＦＴの計算量は O[N log₂(N)] である。これを踏まえ、上述のＳＴＦＴおよびＳＤＦＴにおける各フレームでの計算量を考える。まず、図７のＳＴＦＴでは、ＦＦＴを１回行うので計算量は O[N log₂(N)] となる。図８のＳＤＦＴでは、スライド処理をホップ長Ｈだけ行うので計算量は O(NH) となる。

本発明は、音響信号などの時系列ディジタルデータの周波数分析を効率よく短時間で実行できる周波数分析装置およびプログラムを提供することを目的とする。

上記目的を達成するため、請求項１に係る発明は、入力した時系列ディジタルデータに対し所定のホップ長Ｈ毎に所定のフレーム長Ｎの長さのフレームの周波数分析を行う周波数分析装置であって、既に分析済みのｍ−１番目のフレームの周波数分析結果に所定の位相差を加えることにより今回の分析対象であるｍ番目のフレームの分析結果の第１項を取得する手段と、ｍ−１番目のフレームのサンプルデータと今回の分析対象であるｍ番目のフレームの新たなサンプルデータとの差分をとる手段と、前記差分に対して所定の離散フーリエ変換の処理を施すことにより今回の分析対象であるｍ番目のフレームの分析結果の第２項を取得する手段と、前記第１項に第２項を加えることにより今回の分析対象であるｍ番目のフレームの周波数分析結果を取得する手段とを備えることを特徴とする。

請求項２に係る発明は、入力した時系列ディジタルデータに対し所定のホップ長Ｈ毎に所定のフレーム長Ｎの長さのフレームの周波数分析を行う周波数分析装置であって、既に分析済みのｍ−１番目のフレームの周波数分析結果X_m-1(k)に所定の位相差を加える演算

を施すことにより今回の分析対象であるｍ番目のフレームの分析結果の第１項を取得する手段と、ｍ−１番目のフレームのサンプルデータx(n+mN)と今回の分析対象であるｍ番目のフレームの新たなサンプルデータx(n+mN-N)との差分を

にて求める手段と、前記差分に対して離散フーリエ変換の処理を施すことにより

でｍ番目のフレームの分析結果の第２項を取得する手段と、前記第１項に第２項を加えることにより今回の分析対象であるｍ番目のフレームの周波数分析結果X_m(k)を取得する手段とを備えることを特徴とする。

請求項３に係る発明は、請求項１または２に記載の周波数分析装置において、前記ホップ長Ｈが２のべき乗であり、前記差分に対する離散フーリエ変換の処理が高速フーリエ変換の処理であることを特徴とする。

請求項４に係る発明は、請求項１から３の何れか１つに記載の周波数分析装置において、ホップ長Ｈとフレーム長Ｎとを、ＨＱ＝Ｎ（ＨはＮの約数）を満たすように設定し、前記差分に対する離散フーリエ変換の処理をＱ回の離散フーリエ変換のループ処理で実行することを特徴とする。

請求項５に係る発明は、請求項１から３の何れか１つに記載の周波数分析装置において、ホップ長Ｈとフレーム長Ｎとを、ＨＱ＝Ｎ（ＨはＮの約数）を満たすように設定し、前記差分に対する離散フーリエ変換の処理は、Ｈ点の差分に対する離散フーリエ変換の処理を実行し、得られたＨ点の離散フーリエ変換の結果を用いて補間演算によりＮ点の離散フーリエ変換の結果を取得するものであることを特徴とする。

請求項６に係る発明は、請求項１から５の何れか１つに記載の周波数分析装置において、前記求めたｍ番目のフレームの周波数分析結果に、所定の窓関数の周波数応答を畳み込む演算を行うことにより最適化された周波数分析結果を求めることを特徴とする。

請求項７に係る発明は、入力した時系列ディジタルデータに対し所定のホップ長Ｈ毎に所定のフレーム長Ｎの長さのフレームの周波数分析を行う周波数分析装置を実現する周波数分析プログラムであって、処理装置で実行させることにより、既に分析済みのｍ−１番目のフレームの周波数分析結果に所定の位相差を加えることにより今回の分析対象であるｍ番目のフレームの分析結果の第１項を取得する手段と、ｍ−１番目のフレームのサンプルデータと今回の分析対象であるｍ番目のフレームの新たなサンプルデータとの差分をとる手段と、前記差分に対して所定の離散フーリエ変換の処理を施すことにより今回の分析対象であるｍ番目のフレームの分析結果の第２項を取得する手段と、前記第１項に第２項を加えることにより今回の分析対象であるｍ番目のフレームの周波数分析結果を取得する手段とを備える周波数分析装置として機能させる周波数分析プログラムである。

本発明によれば、時系列ディジタルデータに対する周波数分析を従来よりも高速に効率よく行うことができる。

この発明の実施の形態である周波数分析装置のハードウェア構成を示すブロック図本実施形態における周波数分析の概念図本実施形態における周波数分析の演算の原理を説明するための数式（その１）を示す図本実施形態における周波数分析の演算の原理を説明するための数式（その２）を示す図本実施形態における周波数分析の演算の原理を説明するための数式（その３）を示す図本実施形態の周波数分析装置における周波数分析処理の流れを示すフローチャート（その１）本実施形態の周波数分析装置における周波数分析処理の流れを示すフローチャート（その２）周波数分析処理の変形例を示すフローチャート短時間フーリエ変換（ＳＴＦＴ）の概念図スライディング離散フーリエ変換（ＳＤＦＴ）の概念図

以下、図面を用いて本発明の実施の形態を説明する。

図１は、この発明の実施の形態である周波数分析装置のハードウェア構成を示すブロック図である。ＣＰＵ１０１は、この装置全体の動作を制御する処理装置である。メモリ１０２は、ＣＰＵ１０１が実行する各種のプログラムや各種のデータなどを格納した記憶装置であり、ＲＡＭ、ＲＯＭ、フラッシュメモリ、ハードディスクなどを適宜組み合わせて構成することができる。表示器１０３は、この装置の操作パネル上に設けられた各種の情報を表示するためのディスプレイである。操作子１０４は、この装置の操作パネル上に設けられたユーザが操作するための各種の操作子である。信号処理部１０５は、例えばＤＳＰであり、ＣＰＵ１０１の指示に基づいて各種の信号処理プログラムを実行する。特に、信号処理部１０５は、後述する図５Ａ，５Ｂの周波数分析処理プログラムを実行することにより、波形Ｉ／Ｏ１０６経由で入力した音響信号の周波数分析処理を行い、その分析結果をＣＰＵ１０１に出力する。バス１１０は、これら各部を接続するバスラインであり、コントロールバス、データバス、およびアドレスバスを総称したものである。

図１の装置では、信号処理装置１０５で周波数分析処理を実行したが、ＣＰＵ１０１で実行するように構成してもよい。また、汎用のパーソナルコンピュータ（ＰＣ）で、以下で説明する本実施形態の装置のプログラムを実行させることにより周波数分析装置を実現することもできる。さらに、図１の装置は単体の周波数分析装置として構成したが、例えば図１の装置を楽音信号を発生する電子楽器や音響信号を処理するミキサとして構成し、それらの装置の機能の一部として周波数分析の機能を実現するようにしてもよい。

図２は、図１の装置における周波数分析の概念図を示す。２０１は分析対象の音響信号ｘ（ｔ）の時間的な変化を示すグラフである。２０２は処理対象として最初に切り出したフレーム０を示す。最初のフレーム０については、従来技術の図７の７０２や図８の８０２と同様にしてＤＦＴの演算を行い、分析結果Ｘ₀（ｋ）を得る。

２０６は次の処理対象のフレーム１を示す。フレーム１に対するＤＦＴの演算結果を得るため、本実施形態では、前のフレーム０と今回のフレーム１との差分に着目する。フレーム０と１とでは網掛け２０４の波形サンプルが重複している。また、フレーム０から網掛け２０３の波形サンプルを削除し、網掛け２０５の波形サンプルを追加すればフレーム１となる。そこで、本実施形態では、前のフレームの分析結果に所定の位相差を加えた項と、前のフレームと今回分析対象のフレームとの差分に対するＤＦＴ演算結果を元に求めた補正項とから、今回分析対象のフレームの分析結果を取得する。

以下、そのような本発明に係る分析手法の原理を説明する。図３、図４Ａ、および図４Ｂは、本実施形態における周波数分析の演算の原理を説明するための数式である。

式（１）は、フレーム０に対して施すＤＦＴの演算式を示す。式（１）において、x(n)は入力信号であり時刻nにおける波形サンプル値（複素数表現）を示している。eは自然対数の底、jは虚数単位、πは円周率である。Nはフレーム長、kは周波数ビンを特定する番号（ビン番号と呼ぶものとする）である。計算結果のX₀(k)は、各周波数ビン毎の成分を示す複素数であり、これから各周波数ビンにおける入力信号の振幅と位相が計算できる。

次の式（２）は、フレーム１に対して施すＤＦＴの演算式を示す。式（１）と同様の演算式であるが、フレーム１に対する演算でありN=0からN-1の総和（Σ）をとる演算であるので、式（１）のx(n)の部分はx(n+H)に置き換えている。式（３）は、式（２）の総和の計算を、図２の網掛け２０４の範囲に対する計算に相当する第１項と、網掛け２０５の範囲に対する計算に相当する第２項に分けたものである。式（４）の第１項は、式（３）の第１項のnの範囲をn=0〜N-H-1からn=H〜N-1に書き換え、x(n+H)の+Hの部分を調整した。式（４）の第２項は、式（３）の第２項の指数部の定数をΣの前に出した。式（５）は、式（４）の第１項の指数部の定数をΣの前に出した。式（６）は、式（５）の第１項についてn=0〜N-1の範囲の総和の形とし、その総和の余分を第２項で引いている。式（７）は、式（６）の第２項にexp(-j2π/N)kN=exp(-j2πk)=１を書き加えた。式（７）を（８）（９）と変形していき、最終的に式（１０）を得ることができる。

式（１０）の第１項は、いま求めようとしているフレーム１の直前のフレーム０のＤＦＴ結果X₀(k)に位相差を加えたものである。式（１０）の第２項に含まれているx(n+N)-x(n)は、フレーム１の波形サンプル値からフレーム０の波形サンプル値を引いた差であるから、式（１０）の第２項はフレーム０の分析結果を元にフレーム１の分析結果を求める際の補正項といえる。この補正項は、図２の網掛け２０３の範囲を削除し網掛け２０５の範囲を追加することに相当しているといえる。

図４Ａの式（１１）は、式（１０）の第２項の補正項のΣの部分を取り出したものである。式（１１）のY(k)を用いて、式（１０）は式（１２）のように変形できる。式（１０）の形から式（１１）のY(k)が計算量が多い部分といえるが、式（１１）は差分に対するＤＦＴ演算であり、まとめて計算することができるので、従来の手法より高速に計算できる。

式（１３）は、HQ=N（HはNの約数）と仮定して、式（１１）を変形したものである。HはNの約数であるので、式（１１）のkを^kQ+iに置き換えている（^kは、図４Ａ中の上部に^（サーカムフレックスあるいはハット記号）が付いているkを示すものとする）。この式（１３）は、ＤＦＴをＱ回実行すれば計算できる。従って、式（１３）を利用して計算した場合の計算量は、O[Q*H*log(H)]=O[N*log(H)]である。

図４Ｂは、m-1フレーム目とmフレーム目の表現を示す。式（１５）と（１６）は、m-1フレーム目とmフレーム目のＤＦＴの演算式である。式（１７）は、式（１６）を変形したものであり、図３の式（１０）に相当する。式（１８）は、式（１９）のように定義したY_m(k)を用いて式（１７）を変形したものである。式（２０）は、HQ=N（HはNの約数）と仮定した場合の式（１３）に相当する式である。m-1フレーム目の分析結果X_m-1(k)を元にmフレーム目を計算する際には式（１８）を適用すれば良い。また、HQ=N（HはNの約数）の場合は、式（２０）を適用して高速に計算できる。

図５Ａおよび図５Ｂは、本実施形態の周波数分析装置における周波数分析処理の流れを示すフローチャートである。この処理プログラムは、図１の信号処理部１０５で実行される。なお、このフローチャートの処理は、HQ=N（HはNの約数）と仮定した場合の処理である。

ステップ５０１〜５０４は、最初のフレーム０の演算を行う部分である。ステップ５０１で、最初のフレーム０のN個（Nはフレーム長）の波形サンプルをバッファに取り込む。mはフレーム番号を格納する変数であり、ここでは初期値m=0に設定される。ステップ５０２で、バッファ内のN個の波形サンプルに対して、ＤＦＴを施し、X₀(k)を得る。ステップ５０３で演算結果X₀(k)を出力する。この演算結果はＣＰＵ１０１に出力され、周波数分析結果として利用される。ステップ５０４で、次のフレームの処理のために、フレーム０の波形サンプルx(n)とＤＦＴ結果X₀(k)を一時記憶する。

ステップ５０５から５１７は、mを１から順に歩進しながら、フレームm-1の分析結果を元にフレームmの分析結果を求めるループ処理である。まずステップ５０６では、一時記憶しておいた前のフレームm-1のＤＦＴ結果の各周波数ビンの値に位相差を加える。この計算は、図４Ｂの式（１８）の第１項の計算である。ステップ５０７で、今回の処理対象であるフレームmの新しいH個の波形サンプルを取得する。これは図２で言えば網掛け２０５の範囲の波形サンプルに相当する。次にステップ５０８で、１つ前のフレームm-1の古いH個の波形サンプル（図２の２０３に相当）とフレームmの新しいH個の波形サンプル（図２の２０５に相当）との差をとる。これは図４Ｂの式（１７）のΣの中の差分をとる部分に相当する。

次のステップ５０９〜５１４は、iを0〜Q-1の範囲で歩進させながら、式（１８）の第２項の補正項を求めて第１項に順次加えていくループ処理である。まずステップ５１０では、ステップ５０８で求めた差分に位相調整を施す前処理を行う。これは、１つのiの値での式（２０）のΣの括弧の中の計算である。ステップ５１１で、ステップ５１０の結果である前処理されたH個のサンプルにＤＦＴを施す。これは式（２０）の計算である。ステップ５１２で、ステップ５１１の結果を調整する。これは式（１８）の第２項の計算である。ステップ５１３で、算出しておいた前のフレームm-1のＤＦＴ結果に位相差を加えたもの（すなわちステップ５０６の結果）に、ステップ５１２の結果を加える。以上をiを歩進させながら処理し、Q回処理したらループから抜ける。この時点でのステップ５１３の結果が、フレームmの分析結果X_m(k)である。

ステップ５１５で、分析結果X_m(k)を出力する。ステップ５１６で、次のフレームの処理のために、今回のフレームの波形サンプルと分析結果を一時記憶する。ステップ５１７で、最終フレームまで処理したらループから抜けて処理を終了する。

上記実施形態では、ステップ５１１や式（２０）でＤＦＴの計算をすると説明したが、ホップ長Hが２のべき乗であれば高速フーリエ変換（ＦＦＴ）の手法が適用できる。この場合、Q回のＦＦＴで計算でき、一般にＤＦＴよりＦＦＴのほうが高速に計算できる。そこで、上記実施形態の第１の変形例では、ホップ長Hを２のべき乗とし、ステップ５１１の計算をＦＦＴとする。これにより、さらに高速な分析を実現することができる。

上記実施形態の第２の変形例を説明する。通常、ＤＦＴやＦＦＴの計算では、窓関数を適用して分析結果を最適化することができる。ただし、通常の窓関数はＤＦＴやＦＦＴ演算を施す対象のフレームのデータ列に対して窓関数を掛けて最適化するものであるが、上記実施形態の手法では今回のフレームの分析結果を前のフレームの分析結果を利用して求めているので、前のフレームに窓関数を掛けてしまうと今回のフレームとの差分を求めるときに窓関数を掛けた前フレームの値を元に戻すことなどを考慮しなければならなくなり計算が煩雑になる。そのため本実施形態では、従前の窓関数の手法をそのまま適用できない。

そこで、時間領域での窓関数の掛け算が周波数領域では畳み込み演算に相当することに鑑み、予め適用する窓関数の周波数応答を求めておき、ステップ５１５の直前で、求めた結果X_m(k)に対して該窓関数の周波数応答を畳み込む演算を行う。ステップ５１５では畳み込み演算の結果を分析結果として出力し、ステップ５１６では畳み込み演算を行う前の計算結果を次のフレームのために一時記憶する。そのようにすることで、窓関数が適用され最適化された周波数分析結果を出力することができる。

なお、適切な窓関数を選ぶことにより、その窓関数の周波数応答の振幅特性はメインローブに集中しサイドローブは無視できるようにできる。従って、各周波数ビンのうちこの辺りの周波数ビンについて際だたせたいという範囲があれば、その範囲を際だたせるような周波数応答を持つ窓関数を選んで適用することで、それが実現できる。

また、全周波数ビンに対して上記畳み込みをそれぞれ計算しても良いが、ユーザが注目している周波数帯域の周波数ビンについてだけ上記畳み込みを求めるようにしてもよい。例えば、ユーザにとって不要な低すぎる周波数帯や高すぎる周波数帯については上記畳み込み演算を省略することもできる。

さらに、周波数ビン毎に、違う特性を持つ窓関数の周波数応答を畳み込んでもよい。

上記実施形態の第３の変形例を説明する。第３の変形例では、上述のステップ５０９〜５１４を図６のステップ６０１〜６０４の処理に置き換える。その他の処理は、図５Ａ，５Ｂと同じである。上記実施形態のステップ５０９〜５１４では、Ｈ点の差分に対し、１回毎にＷ_H ^(i/Q)nとのかけ算を実行して、全部でＱ回のＤＦＴ計算を行っている。これに対し、本変形例では、Ｈ点の差分に対し１回のＤＦＴ計算を実行して周波数軸上のＨ点のＤＦＴ結果を求め、そのＨ点のＤＦＴ結果を使って補間することによりＮ点ののＤＦＴ結果を求める。

具体的には、ステップ５０８でＨ点の差分を求めた後、ステップ６０１に進み、求めたＨ点の差分に対してＤＦＴ計算を実行する。次にステップ６０２で、そのＤＦＴ結果をＹ_m(^kQ)とするとともに、それを補間することによりＮ点からなるＹ_m(^kQ+i)を求める。ステップ６０３で、求めたＹ_m(^kQ+i)=Ｙ_m(k)を調整する。ステップ６０４で、算出しておいた前のフレームm-1のＤＦＴ結果に位相差を加えたものに、上記ステップ６０３の結果である調整されたＹ_m(k)を加える。以上でフレームmの分析結果X_m(k)が取得でき、ステップ５１５でそれを出力する。

上記第３の変形例によれば、厳密な計算でなく、簡単な補間計算によって大雑把なＤＦＴ結果を求めることができる。なお、補間は例えば直線補間や４点補間などの従前より知られている補間方法を用いれば良い。

１０１…中央処理装置（ＣＰＵ）、１０２…メモリ、１０３…表示器、１０４…操作子、１０５…信号処理部（ＤＳＰ）、１０６…波形Ｉ／Ｏ。

Claims

入力した時系列ディジタルデータに対し所定のホップ長Ｈ毎に所定のフレーム長Ｎの長さのフレームの周波数分析を行う周波数分析装置であって、
既に分析済みのｍ−１番目のフレームの周波数分析結果に所定の位相差を加えることにより今回の分析対象であるｍ番目のフレームの分析結果の第１項を取得する手段と、
ｍ−１番目のフレームのサンプルデータと今回の分析対象であるｍ番目のフレームの新たなサンプルデータとの差分をとる手段と、
前記差分に対して所定の離散フーリエ変換の処理を施すことにより今回の分析対象であるｍ番目のフレームの分析結果の第２項を取得する手段と、
前記第１項に第２項を加えることにより今回の分析対象であるｍ番目のフレームの周波数分析結果を取得する手段と
を備えることを特徴とする周波数分析装置。
入力した時系列ディジタルデータに対し所定のホップ長Ｈ毎に所定のフレーム長Ｎの長さのフレームの周波数分析を行う周波数分析装置であって、
既に分析済みのｍ−１番目のフレームの周波数分析結果X_m-1(k)に所定の位相差を加える演算

を施すことにより今回の分析対象であるｍ番目のフレームの分析結果の第１項を取得する手段と、
ｍ−１番目のフレームのサンプルデータx(n+mN)と今回の分析対象であるｍ番目のフレームの新たなサンプルデータx(n+mN-N)との差分を

にて求める手段と、
前記差分に対して離散フーリエ変換の処理を施すことにより

でｍ番目のフレームの分析結果の第２項を取得する手段と、
前記第１項に第２項を加えることにより今回の分析対象であるｍ番目のフレームの周波数分析結果X_m(k)を取得する手段と
を備えることを特徴とする周波数分析装置。
請求項１または２に記載の周波数分析装置において、
前記ホップ長Ｈが２のべき乗であり、
前記差分に対する離散フーリエ変換の処理が高速フーリエ変換の処理である
ことを特徴とする周波数分析装置。
請求項１から３の何れか１つに記載の周波数分析装置において、
ホップ長Ｈとフレーム長Ｎとを、ＨＱ＝Ｎ（ＨはＮの約数）を満たすように設定し、前記差分に対する離散フーリエ変換の処理をＱ回の離散フーリエ変換のループ処理で実行する
ことを特徴とする周波数分析装置。
請求項１から３の何れか１つに記載の周波数分析装置において、
ホップ長Ｈとフレーム長Ｎとを、ＨＱ＝Ｎ（ＨはＮの約数）を満たすように設定し、前記差分に対する離散フーリエ変換の処理は、Ｈ点の差分に対する離散フーリエ変換の処理を実行し、得られたＨ点の離散フーリエ変換の結果を用いて補間演算によりＮ点の離散フーリエ変換の結果を取得するものである
ことを特徴とする周波数分析装置。
請求項１から５の何れか１つに記載の周波数分析装置において、
前記求めたｍ番目のフレームの周波数分析結果に、所定の窓関数の周波数応答を畳み込む演算を行うことにより最適化された周波数分析結果を求める
ことを特徴とする周波数分析装置。
入力した時系列ディジタルデータに対し所定のホップ長Ｈ毎に所定のフレーム長Ｎの長さのフレームの周波数分析を行う周波数分析装置を実現する周波数分析プログラムであって、
処理装置で実行させることにより、
既に分析済みのｍ−１番目のフレームの周波数分析結果に所定の位相差を加えることにより今回の分析対象であるｍ番目のフレームの分析結果の第１項を取得する手段と、
ｍ−１番目のフレームのサンプルデータと今回の分析対象であるｍ番目のフレームの新たなサンプルデータとの差分をとる手段と、
前記差分に対して所定の離散フーリエ変換の処理を施すことにより今回の分析対象であるｍ番目のフレームの分析結果の第２項を取得する手段と、
前記第１項に第２項を加えることにより今回の分析対象であるｍ番目のフレームの周波数分析結果を取得する手段と
を備える周波数分析装置として機能させる周波数分析プログラム。