JP2003227854A

JP2003227854A - 周波数推定装置および周波数推定方法

Info

Publication number: JP2003227854A
Application number: JP2002027013A
Authority: JP
Inventors: Wataru Horiuchi; 弥堀内
Original assignee: Mitsubishi Electric Corp
Current assignee: Mitsubishi Electric Corp
Priority date: 2002-02-04
Filing date: 2002-02-04
Publication date: 2003-08-15

Abstract

(57)【要約】【課題】高精度で入力信号中の周波数を推定できる周
波数推定装置を得る。【解決手段】本発明に係る周波数推定装置は、周波数
ω_ｒ成分を含む信号1が入力され第1周波数ω_１で極値を
有する利得特性に基づき処理する第1フィルタ2と、第1フィルタ
出力からアルコ゛リス゛ムに従い計算する第1振幅計算器4と、信
号1が入力され第2周波数ω_２で極値を有する利得特性に
基づき処理する第2フィルタ3と、第2フィルタ出力からアルコ゛リス゛ム
に従い計算する第2振幅計算器5と、第1振幅計算器の結
果6から第2振幅計算器の結果7を減算する減算器8と、減
算結果9に定数kを乗算する係数器10と、係数器の処理結
果11とn-1回周波数推定値13とを加算する加算器12とを
備え、n回周波数推定時の第1の周波数ω_１および第2の
周波数ω_２がn-1回周波数推定値13で決定され、加算器1
2の処理結果をn回周波数推定値14として、減算結果9が
所定値以下まで周波数推定を繰り返すこととした。

Description

【発明の詳細な説明】

【０００１】

【発明の属する技術分野】本発明は、入力信号中の特定
周波数を推定する周波数推定装置および周波数推定方法
に関するものである。

【０００２】

【従来の技術】ある物体がある周波数ω_ｒで振動してお
り、この物体の振動を変位センサで検出するような用途
では、周波数ω_ｒの値を正確に知ることが重要である。
この変位センサ出力信号Ｓ_ＳＥＮＳには一般に物体の振
動周波数ω_ｒの他にノイズ成分が含まれているので、不
要な周波数成分を除去するために適切なフィルタ処理が
必要となる。しかしながら、振動周波数が変動する場
合、どの周波数成分が必要でどの周波数成分が不要なの
かが判別できないため、適切なフィルタを構成すること
が困難であった。

【０００３】未知のシステムの伝達関数を推定する最も
一般的な方法の一つとして、SimonHaykin著“Introduct
ion to Adaptive Filters”，(Macmillan Publishing C
ompany New York， 1984)にあるような適応フィルタが
ある。そして、適応フィルタの各フィルタ係数が安定に
収束することが保証された非再帰型(以後、ＦＩＲ型あ
るいは単にＦＩＲと呼ぶ)フィルタを構成する方法が一
般的であった。

【０００４】しかしながら、かかるＦＩＲ型適応フィル
タでは、同等の性質を有する再帰型（以後、ＩＩＲ型あ
るいは単にＩＩＲと呼ぶ）フィルタに比べてその次元が
高くなる問題があるため、例えば、特開平１１-２０２
９４号公報に開示された振動低減方法では未知システム
の伝達関数を推定するためにＩＩＲ型適応フィルタを用
いていた。かかる周波数推定方法では、適応フィルタの
各フィルタ係数をＬＭＳアルゴリズムに基づいて更新す
ることで、ＦＩＲ型適応フィルタのような安定性は保証
できないものの、比較的低次元のフィルタを構成するこ
とが可能であった。ここで、ＬＭＳとはＬｅａｓｔＭ
ｅａｎＳｑｕａｒｅの略で、最小二乗法に基づくアル
ゴリズムを意味する。

【０００５】

【発明が解決しようとする課題】このように、従来はフ
ィルタ係数の決定には最急降下法に基づくＬＭＳアルゴ
リズムを用いていた。しかし、かかるＬＭＳアルゴリズ
ムを用いたＦＩＲ型適応フィルタは安定性に優れる反
面、複雑な周波数特性への対応が困難であり、また、多
くのタップ係数とメモリ領域が必要であった。一方、Ｉ
ＩＲ型適応フィルタでは低次元化が可能であるという利
点を有する反面、その安定性が保証されないという問題
があった。

【０００６】本発明は、上記問題を解決するためのもの
であり、ある周波数ω_ｒの信号を含む入力信号Ｓ_ＩＮか
ら周波数ω_ｒを、従来の技術に比べて高速かつ精度よく
安定に周波数推定する周波数推定装置および周波数推定
方法を得ることを目的とする。

【０００７】

【課題を解決するための手段】本発明に係る周波数推定
装置は、周波数推定の対象となる周波数ω_ｒの成分を含
む信号が入力され、第１の周波数ω_１で利得が極値を有
する利得特性に基づきフィルタ処理する第１のフィルタ
と、上記第１のフィルタの出力から所定のアルゴリズム
に従い振幅を計算する第１の振幅計算器と、上記周波数
ω_ｒの成分を含む信号が入力され、第２の周波数ω_２で
利得が極値を有する利得特性に基づきフィルタ処理する
第２のフィルタと、上記第２のフィルタの出力から上記
所定のアルゴリズムに従い振幅を計算する第２の振幅計
算器と、上記第１の振幅計算器の処理結果から上記第２
の振幅計算器の処理結果を減算する減算器と、上記減算
器の処理結果に定数を乗算する係数器と、上記係数器の
処理結果とｎ−１回目の周波数推定値とを加算する加算
器と、を備え、ｎ回目の周波数推定時における第１の周
波数ω_１および第２の周波数ω_２が上記ｎ−１回目の周
波数推定値に基づいて決定され、上記加算器の処理結果
をｎ回目の周波数推定値として、上記減算器の処理結果
が所定の値以下になるまで周波数推定を繰り返すことと
した。

【０００８】また、本発明に係る周波数推定装置は、上
記ｎ回目の周波数推定時の第１の周波数ω_１が、上記ｎ
−１回目の周波数推定値であることとした。

【０００９】また、本発明に係る周波数推定装置は、上
記ｎ回目の周波数推定時の第２の周波数ω_２が、上記第
１の周波数ω_１に所定の値を加えた数値であることとし
た。

【００１０】また、本発明に係る周波数推定装置は、上
記第１のフィルタおよび上記第２のフィルタの前段に所
定の周波数帯域のみを通過させる帯域通過特性を具備す
る第３のフィルタを設けることとした。

【００１１】また、本発明に係る周波数推定装置は、上
記第１の振幅計算器と上記減算器との間に高周波成分を
除去する第１の高周波除去フィルタが設けられ、上記第
２の振幅計算器と上記減算器との間に高周波成分を除去
する第２の高周波除去フィルタを設けることとした。

【００１２】また、本発明に係る周波数推定装置は、上
記加算器の後段に、さらに低域通過フィルタを設けるこ
ととした。

【００１３】また、本発明に係る周波数推定装置は、上
記係数器と上記加算器との間に出力制限器を設けること
とした。

【００１４】本発明に係る周波数推定方法は、ｎ−１回
目の周波数推定値を呼び出すステップと、上記ｎ−１回
目の周波数推定値を第１の周波数ω_１とし所定の帯域除
去あるいは帯域通過特性を備えた第１のフィルタを設計
するステップと、入力信号を上記第１のフィルタでフィ
ルタ処理するステップと、上記第１のフィルタの出力か
ら振幅を計算するステップと、上記第１の周波数ω_１よ
りも所定値だけ高い周波数を第２の周波数ω_２とし所定
の帯域除去あるいは帯域通過特性を備えた第２のフィル
タを設計するステップと、上記入力信号を上記第２のフ
ィルタでフィルタ処理するステップと、上記第２のフィ
ルタの出力から振幅を計算するステップと、上記第１の
フィルタの出力によって算出された振幅から上記第２の
フィルタの出力によって算出された振幅を減算するステ
ップと、上記減算結果に所定の係数を乗算した値に上記
ｎ−１回目の周波数推定値を加算することによりｎ回目
の周波数推定値を算出するステップと、を含んでなり、
上記減算結果が所定の値以下になるまで上記ステップを
繰り返すこととした。

【００１５】

【発明の実施の形態】実施の形態１．本発明の実施の形
態１の周波数推定装置を図１〜４に基づいて説明する。
ここで、図１は、本発明の実施の形態１の周波数推定装
置における装置機能を示すブロック図、図２は、第１の
フィルタＦ_１と第２のフィルタＦ_２の周波数に対する利
得特性であり、図３は入力信号Ｓ_ＩＮの周波数成分を表
す図、図４は、第１の処理結果Ｓ_ＯＵＴ１と第２の処理
結果Ｓ_ＯＵＴ２をそれぞれ表す図である。

【００１６】図中、１は入力信号Ｓ_ＩＮ、２は第１のフ
ィルタＦ_１、３は第２のフィルタＦ _２、４は第１の振幅
計算器ＡＭＰ_１、５は第２の振幅計算器ＡＭＰ_２、６は
第１の処理結果Ｓ_ＯＵＴ１、７は第２の処理結果Ｓ
_ＯＵＴ２、８は減算器ＳＵＢ、９は差ｅ、１０は係数器
Ｋ、１１は係数器Ｋの出力、１２は加算器ＡＤＤ、１３
はｎ−１回目の周波数推定結果ω（ｎ−１）、１４はｎ
回目の周波数推定結果ω（ｎ）、をそれぞれ表す。

【００１７】まず、図１に基づいて本発明の実施の形態
１の周波数推定装置の構成を説明する。ある周波数ω_ｒ
の信号を含む入力信号（Ｓ_ＩＮ）１からその周波数ω_ｒ
を推定するために、入力信号（Ｓ_ＩＮ）１を第１のフィ
ルタ（Ｆ_１）２でフィルタ処理し、さらに第１の振幅計
算器（ＡＭＰ_１）４で第１のフィルタ（Ｆ_１）２の出力
を二乗して、第１の処理結果（Ｓ_ＯＵＴ１）６を求め
る。このフィルタ処理は、具体的には、後述する所定の
帯域の周波数のみ通過および増幅させる帯域通過フィル
タ処理である。

【００１８】並行して、入力信号（Ｓ_ＩＮ）１を第２の
フィルタ（Ｆ_２）３でフィルタ処理し、第２の振幅計算
器（ＡＭＰ_２）５で第２のフィルタ（Ｆ_２）３の出力を
二乗して第２の処理結果（Ｓ_ＯＵＴ２）７を求める。減
算器（ＳＵＢ）８は振幅（Ｓ _ＯＵＴ１）６と振幅（Ｓ
_ＯＵＴ２）７の差（ｅ）９を出力し、係数器（Ｋ）１０
では差（ｅ）９が係数倍、つまり定数ｋが乗算される。
加算器（ＡＤＤ）１２では、ｎ−１回目の周波数推定結
果ω（ｎ−１）１３と係数器（Ｋ）１０の出力を加算す
る。かかる加算結果が、ｎ回目の周波数推定結果ω
（ｎ）１４となる。

【００１９】次に、図１に示す本発明の実施の形態１に
おける周波数推定装置の動作を説明する。第１のフィル
タ（Ｆ_１）２と第２のフィルタ（Ｆ_２）３の周波数に対
する利得は、それぞれ図２で表される利得特性を有す
る。一方、入力信号（Ｓ_ＩＮ）１の周波数成分は図３の
特性を有する。

【００２０】かかる特性を具備した第１のフィルタ（Ｆ
_１）２と第２のフィルタ（Ｆ_２）３で入力信号
（Ｓ_ＩＮ）１をフィルタ処理、すなわち、それぞれのフ
ィルタにおける利得特性に基づいた周波数域のみ通過お
よび増幅する処理を行った結果得られた第１の処理結果
（Ｓ_ＯＵＴ１）６と第２の処理結果（Ｓ_ＯＵＴ２）７
は、各フィルタの利得特性を反映して、図４に示したよ
うな周波数特性になる。

【００２１】図４から明らかなように、入力信号（Ｓ
_ＩＮ）１では周波数ω_ｒの成分が支配的なので、第１の
処理結果（Ｓ_ＯＵＴ１）６および第２の処理結果（Ｓ
_ＯＵＴ２）７においても周波数ω_ｒの成分が支配的とな
る。つまり、第１の処理結果（Ｓ _ＯＵＴ１）６では周波
数ω_ｒの部分で入力信号（Ｓ_ＩＮ）１中の周波数ω_ｒに
対応した振幅ａ_１が現われる。なお、周波数ω_１を中心
としたブロードな振幅は入力信号（Ｓ_ＩＮ）１に含まれ
ているホワイトノイズに対応した振幅である。一方、第
２の処理結果（Ｓ_ＯＵＴ２）７においても第１の処理結
果（Ｓ_ＯＵＴ１）６と同様、入力信号（Ｓ_ＩＮ）１中の
周波数ω_ｒに対応した振幅ａ_２が現われるが、第２のフ
ィルタ（Ｆ_２）３では周波数ω_ｒにおける利得は小さい
ので、その振幅値はａ_１に比べて格段に小さい。

【００２２】減算器（ＳＵＢ）８の出力である差（ｅ）
９は、ｅ＝ａ_１−ａ_２＞０ (式１) となるので、係数器（Ｋ）１０における係数ｋを負の数
に設定しておけば、 ω（ｎ）＝ω（ｎ−１）＋ｋ(ａ_１−ａ_２) (式２) から、ｎ回目の周波数推定結果ω（ｎ）１４はｎ−１回
目の周波数推定結果ω（ｎ−１）１３より小さな値とな
る。ここで、ｎは１から始まる自然数を意味する。

【００２３】上述の第１のフィルタ（Ｆ_１）２における
第１の周波数ω_１はω（ｎ−１）１３と一致するように
設定されるので、周波数推定が繰り返されるに従い、周
波数推定結果ω（ｎ）１４は真の周波数ω_ｒに漸近して
いく。ｎ回目の周波数推定結果ω（ｎ）１４は次回、つ
まりｎ＋１回目の周波数推定における第１の周波数ω _１
となり、この周波数推定動作は第１の周波数ω_１が十分
に真の周波数ω_ｒに近づくまで繰り返される。

【００２４】一方、第２のフィルタ（Ｆ_２）３における
第２の周波数ω_２は、上記第１の周波数ω_１に所定値を
加えた値に設定される。仮にかかる所定値をｄWとする
と、この場合のｄWは定数でも変数でもよい。ｄWとして
定数を採用すれば、本実施の形態の周波数推定装置の構
成が簡素化される利点がある。

【００２５】一方、ｄWを変数とすれば、より高精度な
周波数推定が可能となる。例えば、ｄWを周波数推定値
の変動の大きさに比例した変数とすれば、周波数推定値
が収束するに伴い、第１の周波数ω_１と第２の周波数ω
_２の間隔が狭くなり、より高精度な周波数推定が可能と
なる。さらに、ｄWの変動範囲を限定すれば、下限値ｄW
Bで必要な周波数推定精度が、また、上限値ｄWTで推定
動作の収束速度が決められる。ここでは、ｄWを周波数
推定値の変動に比例するように選択したが、周波数推定
値の変動に応じて離散的に変化する数値として選択して
も良い。例えば、周波数推定値の変動が大きい場合はｄ
Wを大きな定数ｄWTとし、周波数推定値の変動が小さい
場合はｄWを小さな定数ｄWBとすれば、ｄWを連続的に変
化する変数としたときと同様の効果が得られる。

【００２６】真の周波数ω_ｒが第２の周波数ω_２よりも
大きい場合は、式（１）で表される差（ｅ）９が負にな
るので、式（２）の右辺第２項は正となる。つまり、ｎ
回目の周波数推定結果ω（ｎ）１４はｎ−１回目の周波
数推定値ω（ｎ−１）１３よりも大きくなり、やはり真
の周波数ω_ｒに漸近していく。

【００２７】次に、上記のような周波数推定を繰り返し
た結果、周波数推定結果ω（ｎ）１４が第１の周波数ω
_１に十分に漸近すると、式（１）の差（ｅ）９の大きさ
（絶対値）は十分に小さくなっているので、周波数推定
結果ω（ｎ）１４の変動も十分に小さくなる。よって、
最終的な周波数推定結果ω（ｎ）１４、つまり、最終的
に得られた第１の周波数ω_１は、 ω_１＝ω（ｎ）＜ ω_ｒ＜ ω_２(式３) なる状態に収束するか、あるいは、式（３）で表される
状態のごく近傍で、ω（ｎ）１４が微小な振動を繰り返
す結果となる。

【００２８】したがって、第２の周波数ω_２と第１の周
波数ω_１との差を周波数推定における許容誤差よりも小
さく設定しておけば、式（３）が示すように、周波数推
定結果ω（ｎ）１４はほぼ真の周波数ω_ｒに一致する結
果が得られる。

【００２９】このようにして、簡単かつ所望の精度で入
力信号（Ｓ_ＩＮ）１中の周波数ω_ｒが推定できる。ま
た、この周波数推定結果ω（ｎ）１４をパラメータとす
るフィルタが必要であれば、ω（ｎ）１４に基づきフィ
ルタ係数を計算すればよいから、例えＩＩＲ型フィルタ
であっても、安定な適応フィルタとして構成可能とな
る。すなわち、本周波数推定装置を用いれば安定なＩＩ
Ｒ型適応フィルタが実現できる。

【００３０】さらに、係数器（Ｋ）１０における係数ｋ
の大きさによって、意図的に収束の速度を変化させるこ
とが可能となる。つまり、係数ｋが大きくすれば、周波
数推定の収束速度を速くすることができる。また、入力
信号（Ｓ_ＩＮ）１に含まれる周波数ω_ｒの成分の大きさ
が大きいほど、周波数推定の収束速度が速くなる。しか
し、収束速度が速すぎると、収束が振動的になったり式
（３）の状態の近傍での振動幅が増大するので、実際の
入力信号（Ｓ_ＩＮ）１の特性に合わせて係数ｋを調整す
ればよい。

【００３１】なお、第１の振幅計算器（ＡＭＰ_１）４お
よび第２の振幅計算器（ＡＭＰ_２）５は、第１のフィル
タ（Ｆ_１）２および第２のフィルタ（Ｆ_２）３との振幅
を計算する機能を具備するものであり、本実施の形態で
は振幅計算を二乗計算として説明した。しかしながら、
他にも同等の機能を有する計算法があり、例えば、周波
数解析に基づいて振幅を計算する、ＦＦＴ解析して最大
振幅を選択する、一定期間内でのピーク値で代表させ
る、時間波形の絶対値を計算する、ＲＭＳ(二乗平均)処
理する、等が挙げられる。

【００３２】また、上述の説明では第１のフィルタ（Ｆ
_１）２および第２のフィルタ（Ｆ_２）３を、第１の周波
数ω_１および第２の周波数ω_２において極大値を有する
と仮定したが、逆に極小値を有するものであってもよ
い。ただし、この場合には、係数器Ｋにおける係数ｋを
正の数とする必要がある。そうすれば、式（２）の計算
がそのまま使用可能となり、式（３）に示すように、周
波数推定値ω（ｎ）１４が真の周波数ω_ｒに十分漸近す
るので、有効な周波数推定が実行できる。

【００３３】実施の形態２．図5は、本発明の実施の形
態２の周波数推定装置の機能を示すブロック図である。
図5に基づいて装置構成を説明する。図中、１５は、第
３のフィルタＦ_３を示す。

【００３４】ある周波数ω_ｒの信号を含む入力信号（Ｓ
_ＩＮ）１に対してその周波数ω_ｒを推定するために、先
ず、入力信号（Ｓ_ＩＮ）１を第３のフィルタ（Ｆ_３）１
５でフィルタ処理する。第３のフィルタ（Ｆ_３）１５は
帯域通過フィルタであり、その通過帯域は真の周波数ω
_ｒの変動範囲と一致するよう予め設定しておく。なお、
かかる変動範囲は図6中の斜線部に該当する。本装置の
他の構成は実施の形態１の周波数推定装置と全く同一で
ある。

【００３５】次に、本発明の実施の形態２の周波数推定
装置の動作を説明する。入力信号（Ｓ_ＩＮ）１には、図
３で示した周波数推定の対象となる周波数ω_ｒとホワイ
トノイズの他に、さらに不要な周波数成分ω_４が含まれ
ていることもある。このような場合に、実施の形態1に
示した構成で図６のような入力信号（Ｓ_ＩＮ）１を処理
すると、周波数推定結果ω（ｎ）１４はω_ｒとω_４のど
ちらに収束するかが確定しないので、有効な周波数推定
を実行できない。そこで、かかる不具合を防止するため
に、周波数推定装置を図５のブロック図に示すような第
３のフィルタ（Ｆ_３）１５を加えた構成とした。ここ
で、第１のフィルタ（Ｆ_１）２と第２のフィルタ
（Ｆ_２）３の利得は実施の形態１と同様、図2で示され
た利得特性を具備する。一方、入力信号（Ｓ_ＩＮ）１の
周波数成分は図６に示す周波数特性を有するものとす
る。

【００３６】実施の形態２に係る周波数測定装置では、
上述したように所定の帯域通過特性を有する第３のフィ
ルタ（Ｆ_３）１５を使用している。真の周波数ω_ｒの変
動範囲が予め判明している場合は、その変動範囲と第３
のフィルタ（Ｆ_３）１５の通過帯域とを一致させてお
く。この結果、第１のフィルタ（Ｆ_１）２および第２の
フィルタ（Ｆ_２）３への入力信号（Ｓ_ＩＮ）１の周波数
特性は、図７で表される特性を呈する。

【００３７】すなわち、図６と図７の各特性を比較すれ
ば判るように、第１のフィルタ（Ｆ _１）２および第２の
フィルタ（Ｆ_２）３への入力信号（Ｓ_ＩＮ）１中の不必
要な周波数成分ω_４が除去される。つまり、不必要な周
波数成分ω_４は第３のフィルタ（Ｆ_３）１５の通過帯域
外だからである。以降の動作は実施の形態1と同様であ
り、周波数推定結果ω（ｎ）１４は真の周波数ω_ｒに十
分漸近していく。

【００３８】本実施の形態においても、第１のフィルタ
（Ｆ_１）２および第２のフィルタ（Ｆ_２）３では、第１
の周波数ω_１および第２の周波数ω_２において極大値を
有すると仮定したが、逆に極小値を有するものでもよ
い。ただし、この場合は、係数器（Ｋ）１０における係
数ｋを正の数としなければならない。そうすれば、式
（２）の計算式がそのまま使用でき、式（３）となるよ
うに周波数推定値ω（ｎ）１４が真の周波数ω_ｒに十分
漸近していく。

【００３９】このようにして、真の周波数ω_ｒの変動範
囲が予め判明している場合、その変動範囲と通過帯域と
を一致させた第３のフィルタ（Ｆ_３）１５を第１のフィ
ルタ（Ｆ_１）２および第２のフィルタ（Ｆ_２）３の前段
に設けることにより、不要な周波数成分が確実に除去さ
れ、簡単かつ所望の精度で周波数が推定できる。

【００４０】なお本実施の形態では、第３のフィルタ
（Ｆ_３）１５には帯域通過型のタイプを採用したが、不
要な周波数成分自体が予め判明している場合は、この不
要な周波数成分を選択的に除去する帯域除去型のフィル
タを第３のフィルタ（Ｆ_３）１５として採用しても、本
実施の形態と同様の効果が得られることは言うまでもな
い。

【００４１】実施の形態３．上記の実施の形態1および2
では、第１の振幅推定器（ＡＭＰ_１）４および第２の振
幅推定器（ＡＭＰ_２）５を乗算による計算で近似的に処
理したが、本実施の形態による周波数推定装置では、図
8のブロック図に示すように第１の振幅推定器（ＡＭＰ
_１）４を第１の乗算器（ＭＵＬ_１)１６と第１の高周波
除去フィルタ（Ｆ_４)１７、第２の振幅推定器（ＡＭＰ
_２)５を乗算器（ＭＵＬ_２)１８と第２の高周波除去フィ
ルタ(Ｆ_５)１９で構成している。

【００４２】第１のフィルタ（Ｆ_１）２および第２のフ
ィルタ（Ｆ_２）３の出力は時々刻々変化するものであ
り、その主要な変動周波数はω_ｒである。したがって、
振幅推定器（ＡＭＰ）を単に乗算器（ＭＵＬ）だけで構
成すると、その出力は2×ω_ｒの周波数で脈動するおそ
れがある。さらに、第１フィルタ（Ｆ_１）２と第２のフ
ィルタ（Ｆ_２）３の利得特性が相違するため、第１の処
理結果（Ｓ_ＯＵＴ１）６および第２の処理結果（Ｓ
_ＯＵＴ２）７に見られる脈動の位相がそれぞれ異なり、
図1の減算器（ＳＵＢ）８の出力である差（ｅ）９では
さらに複雑な脈動が生じうる。安定な周波数推定を妨害
するかかる脈動を抑制するために、本実施の形態の周波
数推定装置では、カットオフ周波数が周波数ω_ｒよりも
十分低い第１の高周波除去フィルタ（Ｆ_４）１７および
第２の高周波除去フィルタ（Ｆ_５）１９を新たに設けて
いる。振幅推定器（ＡＭＰ）以外の装置構成と周波数推
定動作は実施の形態1と同一である。

【００４３】かかる装置構成によって出力脈動を防止す
ることができるので、振幅計算器（ＡＭＰ_１）４および
（ＡＭＰ_２）５の出力が安定し、結果的に周波数推定動
作が安定化する。また、周波数解析に基づく振幅計算に
比べて簡単に振幅が計算できる。さらに、実施の形態２
に示す帯域通過型の第３のフィルタ（Ｆ_３）１５と組み
合わせると、より高精度に周波数を推定できる。

【００４４】実施の形態４．本発明の実施の形態４にお
ける周波数推定装置は、図1に示された実施の形態１の
周波数推定装置中の加算器（ＡＤＤ）１２の後段に低域
通過フィルタを挿入したものである。実施の形態４にお
ける周波数推定装置ブロック図を図9に示す。図中、２
０は低域通過フィルタである。

【００４５】周波数推定の収束目標は真の周波数ω_ｒで
なければならないが、ω_ｒ成分よりも大きなノイズが入
力信号（Ｓ_ＩＮ）１に混入した場合は、一時的に周波数
推定結果が望ましくない応答をすることがある。本実施
の形態では、加算器（ＡＤＤ）１２の後段に低域通過フ
ィルタ２０を追加することで、周波数推定結果が不必要
に変動するのを抑制している。かかる低域通過フィルタ
２０は、不要な高周波成分を除去する機能を果たすから
である。

【００４６】実施の形態４の周波数推定装置によると、
入力信号（Ｓ_ＩＮ）１中に単発的で大きなノイズが混入
しても、周波数推定結果における不必要な変動が抑制さ
れ、安定な周波数推定が実現できる。また、本実施の形
態の周波数推定装置の構成を実施の形態２あるいは実施
の形態３と併用することにより、より安定で高精度な周
波数推定が可能となる。

【００４７】実施の形態５．本実施の形態における周波
数推定装置は、図１０のブロック図に示すように、図１
で表された実施の形態１の周波数推定装置中の係数器
（Ｋ）１０と加算器（ＡＤＤ）１２との間に出力制限器
（Ｌ）２１を加えたものである。また、出力制限器
（Ｌ）２１の入出力特性を図１１に示す。

【００４８】周波数推定動作の概略は実施の形態1と同
様であるが、図１１に示すように、入力値が所定の値以
上になった場合に出力を一定値に制限する機能を有する
出力制限器（Ｌ）２１を用いることにより係数器（Ｋ）
１０の出力が所定の大きさ以上になることを防止してい
る。つまり、周波数推定における真の値への収束速度に
上下限を設けたことになる。この結果、周波数推定結果
ω（ｎ）１４の更新速度の最大値が制限され、不用意に
大きな更新動作を回避できる。

【００４９】実施の形態５の周波数推定装置によると、
単発的で大きなノイズの混入や、入力信号（Ｓ_ＩＮ）１
の入力レベルの増大、あるいは係数器（Ｋ）１０の係数
ｋにおける例えば不適切に大きな値の設定のような設定
ミス等があっても、周波数推定動作が過剰に反応するこ
とを防止でき、安定に周波数推定することが可能とな
る。

【００５０】さらに、本実施の形態の周波数推定装置を
実施の形態２から実施の形態４に示した装置と併用する
ことにより、より安定で高精度な周波数推定が可能とな
る。

【００５１】実施の形態６．図１２は本実施の形態１の
周波数推定装置を機械共振のある位置決め制御系に適用
した場合のブロック図を示すものである。図中、２２は
機械共振のある位置決め制御系、２３は位置指令、２４
は例えばＰＩＤのように利得と位相を調整して制御系全
体を安定化させる制御補償器、２５は必要な制御力を制
御対象に与えるアクチュエータ、２６は制御対象、２７
は制御対象の位置を検出するセンサ、２８は実際の位
置、２９は適応ノッチフィルタ、３０はセンサ信号、３
１はノッチフィルタ出力、をそれぞれ示す。

【００５２】制御帯域よりも高い周波数に機械共振があ
る位置決め制御では、共振モードを安定化するために共
振周波数あるいは共振周波数よりもやや低い周波数にノ
ッチフィルタを挿入することが多い。本実施の形態は、
このような制御系に挿入するノッチフィルタを本発明の
実施の形態１の周波数推定装置、つまり、適応フィルタ
で構成したものである。なお、ノッチフィルタとは、特
定の周波数成分のみを除去し、それ以外の周波数成分を
通過させるフィルタを指し、帯域除去フィルタとも呼ば
れる。

【００５３】ここで、制御対象２６が上記の位置決め制
御系２２の帯域よりも高い周波数の機械共振を有してい
るとする。機械共振での一巡伝達関数の利得を低減する
ために、あるいは機械共振での位相を進行させるため
に、制御ループ内のいずれかの箇所にノッチフィルタを
挿入すると、機械共振が安定化される効果がある。

【００５４】本実施の形態の周波数推定装置の動作を以
下に説明する。まず、機械共振の周波数成分を含んだセ
ンサからの信号３０を、実施の形態1における周波数推
定装置、つまり適応ノッチフィルタ２９に入力する。こ
の際、第１のフィルタ（Ｆ_１）２および第２のフィルタ
（Ｆ_２）３は、挿入したいノッチフィルタと同一形式と
する。たとえば以下の式（４）に示すように、第１のフ
ィルタ（Ｆ_１）２および第２のフィルタ（Ｆ_２）３を２
次のＩＩＲ型ノッチフィルタとする。

【数１】そして、式（４）における各フィルタ係数a、 b、 c、
dは、第１の周波数ω _１あるいは第２の周波数ω_２の関
数である。この理由としては、デジタル制御系のサンプ
リング周波数をω_Ｓとすれば、 ω_A1 = tan (πω₁ / ω_S)、 ω_A2 = tan (πω₂ / ω_S) (式５) なる中心周波数でノッチフィルタを形成するアナログフ
ィルタを構成し、これを双一次z変換したものが式
（４）の第１のフィルタ（Ｆ_１）２および第２のフィル
タ（Ｆ_２）３だからである。

【００５５】適用ノッチフィルタ２９はセンサ信号３０
に含まれる共振周波数ω_ｒに対して推定動作を繰り返
し、最終的に第１の周波数ω_１はほぼ共振周波数ω_ｒと
一致する。すると、ノッチフィルタ出力３１は自動的に
センサからの信号３０から共振周波数ω_ｒの成分を除去
した信号となる。こうして位置決め制御系における機械
共振が安定化される。適応ノッチフィルタ２９の推定周
波数が収束する速度は、係数器（Ｋ）１０の係数ｋと、
センサ信号３０に含まれる機械共振周波数ω_ｒの成分の
大きさと、第１のフィルタ（Ｆ_１）２および第２のフィ
ルタ（Ｆ_２）３の利得特性によって決定される。一方、
周波数推定精度はこれらに加えて第２の周波数ω_２と第
１の周波数ω_１との差の大きさによって決定される。

【００５６】このように、本発明の実施の形態６におけ
る周波数推定装置における第１のフィルタ（Ｆ_２）２お
よび第２のフィルタ（Ｆ_２）３を所望のフィルタ形状に
設定することで、第１のフィルタ（Ｆ_２）２の出力ある
いは第２のフィルタ（Ｆ_２）３の出力が自動的に所望の
利得特性を具備した適応フィルタとなるので、新たなフ
ィルタの設置が不必要となる結果、適応フィルタの構成
が簡素化される。また、この適応フィルタの周波数推定
が高速かつ高精度に実行できる。本実施の形態を上記の
実施の形態２から実施の形態５に示す構成と組み合わせ
れば、さらに、高速かつ高精度で安定な周波数推定が可
能となる。

【００５７】実施の形態７．本実施の形態における周波
数推定装置は、上記の実施の形態1の周波数測定装置に
おいて、第１のフィルタ（Ｆ_１）２および第２のフィル
タ（Ｆ_２）３の利得特性を図１３に示すように構成した
ものである。かかる利得特性は、実施の形態１における
第１のフィルタ（Ｆ_１）２および第２のフィルタ
（Ｆ_２）３の利得特性（図２）と比較して、中心周波数
に対して利得がより緩やかに減衰するという特徴を有し
ている。

【００５８】本実施の形態における周波数推定動作は以
下のようになる。真の周波数ω_ｒに対する周波数推定結
果ω（ｎ）１４の誤差が大きければ、式（１）で表され
る差（ｅ）９が増大し、この結果、式（２）の右辺第２
項で表される更新量も大きくなるので、ω（ｎ）１４は
速やかにω_ｒに漸近する。逆にω（ｎ）１４がω_ｒに近
づくと、差（ｅ）９が小さくなるので、ω（ｎ）１４は
ゆっくりとω_ｒに漸近してゆく。つまり、図１３に示す
第１のフィルタ（Ｆ_１）２および第２のフィルタ
（Ｆ_２）３の利得特性を反映して、推定誤差が大きいと
きは速やかな収束動作となり、逆に推定誤差が小さいと
ゆっくりとした推定動作となる。

【００５９】これにより、真の周波数ω_ｒが変動した場
合には周波数推定値が速やかに周波数ω_ｒに追従する一
方、真の周波数ω_ｒがほぼ一定の場合に周波数推定動作
が振動的になることを防止するので、高速で安定かつ高
精度な周波数推定ができる。さらに、本実施の形態を上
記の実施の形態２から実施の形態６に示す構成と組み合
わせることで、より高速で安定かつ高精度な周波数推定
ができる。

【００６０】実施の形態８．本実施の形態の周波数推定
方法のフローチャートを図１４に示す。図１４における
’はじめ’とは、適応フィルタのプログラム実行開始
を意味する。まず、ｎ−１回目の周波数推定値を呼び出
す。但し、１回目の周波数推定動作では当然前回の周波
数推定値は存在しないので、適当な初期値を用いるが、
ゼロとしても良い。続いて、第１のフィルタ（Ｆ_１）２
と第２のフィルタ（Ｆ_２）３のフィルタ係数をｎ−１回
目の周波数推定値に基づいて計算する。そして、フィル
タ係数が決定された第１のフィルタ（Ｆ_１）２および第
２のフィルタ（Ｆ_２）３に周波数推定の対象となる周波
数ω_ｒの成分を含んだ信号を入力して、それぞれフィル
タ処理する。次に、第１のフィルタ（Ｆ_１）２の出力の
振幅および第２のフィルタ（Ｆ_２）３の出力の振幅を計
算する。さらに、両者の振幅の差を計算し、最後に、こ
の振幅の差を所定の係数ｋで乗算した値をｎ−１回目の
周波数推定値に加算して、ｎ回目の周波数推定値とす
る。

【００６１】なお、上述の第１のフィルタ（Ｆ_１）２お
よび第２のフィルタ（Ｆ_２）３は、ともに帯域除去フィ
ルタとするか、あるいはともに帯域通過フィルタとす
る。

【００６２】上記の一連の処理を繰り返すことで、第１
のフィルタ（Ｆ_１）２および第２のフィルタ（Ｆ_２）３
の中心周波数ω_１、ω_２が次第に推定したい周波数ω_ｒ
に漸近し、十分な時間経過後にはこれらの中心周波数ω
_１、ω_２は入力信号の中で最も支配的な成分の周波数に
収束するか、あるいは推定対象の周波数ω_ｒのごく近傍
で振動する。そして、入力信号の支配的な周波数成分ω
_ｒが変動すると周波数推定値は直ちにこれを検知し、二
つのフィルタの中心周波数ω_１、ω_２がその支配的な周
波数成分ω_ｒに収束あるいは十分に漸近する。上記振幅
の差が所定の値以下であれば、ｎ回目の周波数推定値を
周波数ω_ｒの推定値に決定する。図１４のフローチャー
トに’おわり’がないのは、周波数推定値が入力信号に
含まれる最も支配的な成分の周波数に常に追従する状態
にあるからである。

【００６３】本実施の形態の周波数推定方法では、最急
降下法に基づくＬＭＳアルゴリズム等の複雑な計算を必
要としないので、より高速な周波数推定が可能となる。
さらに、第１のフィルタ（Ｆ_１）２および第２のフィル
タ（Ｆ_２）３をＩＩＲ型のフィルタに設定しても、本実
施の形態の周波数推定方法は何ら変更することなく使用
できる。したがって、ＬＭＳアルゴリズムによる適応フ
ィルタでは保証が不可能あるいは困難なＩＩＲ型適応フ
ィルタに対しても、本実施の形態の周波数推定方法は有
効である。また、第１のフィルタ（Ｆ_１）２と第２のフ
ィルタ（Ｆ_２）３の中心周波数ω_１、ω_２のずらし方
や、振幅の差に作用させる係数の値を適切に選択するこ
とで、安定でかつ高精度な周波数推定ができる。さら
に、最終的にはこれらのフィルタが自動的に推定したい
周波数に対して中心周波数ω_１、ω_２とする帯域除去あ
るいは帯域通過フィルタとなるので、上述の周波数推定
方法は適応フィルタの係数を決定するための適応アルゴ
リズムとしても使える。

【００６４】本実施の形態では、上記の実施の形態１と
ほぼ同じ構成の周波数推定装置に対する周波数推定方法
を示したが、上記の実施の形態２から実施の形態７に示
した構成要素による計算処理等を本実施の形態に追加す
れば、さらに高速かつ高精度で安定な周波数推定が可能
となる。

【００６５】

【発明の効果】本発明に係る周波数推定装置では、周波
数推定の対象となる周波数ω_ｒの成分を含む信号が入力
され、第１の周波数ω_１で利得が極値を有する利得特性
に基づきフィルタ処理する第１のフィルタと、上記第１
のフィルタの出力から所定のアルゴリズムに従い振幅を
計算する第１の振幅計算器と、上記周波数ω_ｒの成分を
含む信号が入力され、第２の周波数ω_２で利得が極値を
有する利得特性に基づきフィルタ処理する第２のフィル
タと、上記第２のフィルタの出力から上記所定のアルゴ
リズムに従い振幅を計算する第２の振幅計算器と、上記
第１の振幅計算器の処理結果から上記第２の振幅計算器
の処理結果を減算する減算器と、上記減算器の処理結果
に定数を乗算する係数器と、上記係数器の処理結果とｎ
−１回目の周波数推定値とを加算する加算器と、を備
え、ｎ回目の周波数推定時における第１の周波数ω_１お
よび第２の周波数ω_２が上記ｎ−１回目の周波数推定値
に基づいて決定され、上記加算器の処理結果をｎ回目の
周波数推定値として、上記減算器の処理結果が所定の値
以下になるまで周波数推定を繰り返すこととしたので、
簡単かつ所望の精度で入力信号Ｓ_ＩＮ中の周波数ω_ｒが
推定でき、また、この周波数推定結果ω（ｎ）をパラメ
ータとするフィルタが必要であれば、ω（ｎ）に基づき
フィルタ係数を計算すればよいから、例えＩＩＲ型フィ
ルタであっても、安定な適応フィルタとして構成可能と
なる。

【００６６】また、本発明に係る周波数推定装置では、
上記ｎ回目の周波数推定時の第１の周波数ω_１を、上記
ｎ−１回目の周波数推定値であることとしたので、入力
信号Ｓ_ＩＮ中の周波数ω_ｒの推定が容易になる。

【００６７】また、本発明に係る周波数推定装置では、
上記ｎ回目の周波数推定時の第２の周波数ω_２を、上記
第１の周波数ω_１に所定の値を加えた数値であることと
したので、より入力信号Ｓ_ＩＮ中の周波数ω_ｒの推定が
容易になる。

【００６８】また、本発明に係る周波数推定装置では、
上記第１のフィルタおよび上記第２のフィルタの前段に
所定の周波数帯域のみを通過させる帯域通過特性を具備
する第３のフィルタを設けることとしたので、不要な周
波数成分が確実に除去され、簡単かつ所望の精度で周波
数が推定できる。

【００６９】また、本発明に係る周波数推定装置では、
上記第１の振幅計算器と上記減算器との間に高周波成分
を除去する第１の高周波除去フィルタが設けられ、上記
第２の振幅計算器と上記減算器との間に高周波成分を除
去する第２の高周波除去フィルタを設けることとしたの
で、出力脈動を防止することができるため、振幅計算器
の出力が安定し、結果的に周波数推定動作が安定化す
る。

【００７０】また、本発明に係る周波数推定装置では、
上記加算器の後段に、さらに低域通過フィルタを設ける
こととしたので、入力信号Ｓ_ＩＮ中に単発的で大きなノ
イズが混入しても、周波数推定結果における不必要な変
動が抑制され、安定な周波数推定が実現できる。

【００７１】また、本発明に係る周波数推定装置では、
上記係数器と上記加算器との間に出力制限器を設けるこ
ととしたので、単発的で大きなノイズの混入や、入力信
号Ｓ _ＩＮの入力レベルの増大、あるいは係数器Ｋの係数
ｋにおける例えば不適切に大きな値の設定のような設定
ミス等があっても、周波数推定動作が過剰に反応するこ
とを防止でき、安定に周波数推定することができる。

【００７２】本発明に係る周波数推定方法では、ｎ−１
回目の周波数推定値を呼び出すステップと、上記ｎ−１
回目の周波数推定値を第１の周波数ω_１とし、所定の帯
域除去あるいは帯域通過特性を備えた第１のフィルタを
設計するステップと、入力信号を上記第１のフィルタで
フィルタ処理するステップと、上記第１のフィルタの出
力から振幅を計算するステップと、上記第１の周波数ω
_１よりも所定値だけ高い周波数を第２の周波数ω_２と
し、所定の帯域除去あるいは帯域通過特性を備えた第２
のフィルタを設計するステップと、上記入力信号を上記
第２のフィルタでフィルタ処理するステップと、上記第
２のフィルタの出力から振幅を計算するステップと、上
記第１のフィルタの出力によって算出された振幅から上
記第２のフィルタの出力によって算出された振幅を減算
するステップと、上記減算結果に所定の係数を乗算した
値に上記ｎ−１回目の周波数推定値を加算することによ
りｎ回目の周波数推定値を算出するステップと、を含ん
でなり、上記減算結果が所定の値以下になるまで上記ス
テップを繰り返すこととしたので、最急降下法に基づく
ＬＭＳアルゴリズム等の複雑な計算を必要としないた
め、より高速な周波数推定が可能となる。

【図面の簡単な説明】

【図１】実施の形態１の周波数推定装置における装置
機能を示すブロック図である。

【図２】実施の形態１の周波数推定装置における第１
のフィルタＦ_１と第２のフィルタＦ_２の周波数に対する
利得特性を示す図である。

【図３】実施の形態１の周波数推定装置における入力
信号Ｓ_ＩＮの周波数成分を表す図である。

【図４】実施の形態１の周波数推定装置における第１
の処理結果Ｓ_ＯＵＴ _１と第２の処理結果Ｓ_ＯＵＴ２をそ
れぞれ表す図である。

【図５】実施の形態２の周波数推定装置の機能を示す
ブロック図である。

【図６】実施の形態２の周波数推定装置における入力
信号Ｓ_ＩＮの周波数特性を示す図である。

【図７】実施の形態２の周波数推定装置における第１
のフィルタＦ_１および第２のフィルタＦ_２への入力信号
の周波数特性を示す図である。

【図８】実施の形態３の周波数推定装置の機能を示す
ブロック図である。

【図９】実施の形態４の周波数推定装置の機能を示す
ブロック図である。

【図１０】実施の形態５の周波数推定装置の機能を示
すブロック図である。

【図１１】実施の形態５の周波数推定装置における出
力制限器Ｌ２１の入出力特性を示す図である。

【図１２】実施の形態６における機械共振のある位置
決め制御系に適用した場合のブロック図を示すものであ
る。

【図１３】実施の形態７の周波数推定装置における第
１のフィルタＦ_１と第２のフィルタＦ_２の周波数に対す
る利得特性を示す図である。

【図１４】実施の形態８における周波数推定方法のフ
ローチャートである。

【符号の説明】

１入力信号Ｓ_ＩＮ、２第１のフィルタＦ_１、３
第２のフィルタＦ_２、４第１の振幅計算器ＡＭＰ
_１、５第２の振幅計算器ＡＭＰ_２、６第１の処理
結果Ｓ_ＯＵＴ１、７第２の処理結果Ｓ_ＯＵＴ２、
８減算器ＳＵＢ、９差ｅ、１０係数器Ｋ、
１１係数器Ｋの出力、１２加算器ＡＤＤ、１３
ｎ−１回目の周波数推定結果ω（ｎ−１）、１４
ｎ回目の周波数推定結果ω（ｎ）、１５第３のフィ
ルタＦ_３、１６第１の乗算器ＭＵＬ_１、１７第
１の高周波除去フィルタＦ_４、１８乗算器ＭＵ
Ｌ_２、１９第２の高周波除去フィルタＦ_５、２０
低域通過フィルタ、２１出力制限器Ｌ、２２機
械共振のある位置決め制御系、２３位置指令、２４
制御補償器、２５アクチュエータ、２６制御
対象、２７センサ、２８実際の位置、２９
適応ノッチフィルタ、３０センサ信号、３１ノ
ッチフィルタ出力。

Claims

【特許請求の範囲】

【請求項１】周波数推定の対象となる周波数ω_ｒの成
分を含む信号が入力され、第１の周波数ω_１で利得が極
値を有する利得特性に基づきフィルタ処理する第１のフ
ィルタと、前記第１のフィルタの出力から所定のアルゴ
リズムに従い振幅を計算する第１の振幅計算器と、前記
周波数ω_ｒの成分を含む信号が入力され、第２の周波数
ω_２で利得が極値を有する利得特性に基づきフィルタ処
理する第２のフィルタと、前記第２のフィルタの出力か
ら前記所定のアルゴリズムに従い振幅を計算する第２の
振幅計算器と、前記第１の振幅計算器の処理結果から前
記第２の振幅計算器の処理結果を減算する減算器と、前
記減算器の処理結果に定数を乗算する係数器と、前記係
数器の処理結果とｎ−１回目の周波数推定値とを加算す
る加算器と、を備え、ｎ回目の周波数推定時における第
１の周波数ω_１および第２の周波数ω_２が前記ｎ−１回
目の周波数推定値に基づいて決定され、前記加算器の処
理結果をｎ回目の周波数推定値として、前記減算器の処
理結果が所定の値以下になるまで周波数推定を繰り返す
ことを特徴とする周波数推定装置。
【請求項２】前記ｎ回目の周波数推定時の第１の周波
数ω_１が、前記ｎ−１回目の周波数推定値であることを
特徴とする請求項１記載の周波数推定装置。
【請求項３】前記ｎ回目の周波数推定時の第２の周波
数ω_２が、前記第１の周波数ω_１に所定の値を加えた数
値であることを特徴とする請求項２記載の周波数推定装
置。
【請求項４】前記第１のフィルタおよび前記第２のフ
ィルタの前段に所定の周波数帯域のみを通過させる帯域
通過特性を具備する第３のフィルタが設けられているこ
とを特徴とする請求項１ないし３のいずれか１項記載の
周波数推定装置。
【請求項５】前記第１の振幅計算器と前記減算器との
間に高周波成分を除去する第１の高周波除去フィルタが
設られ、前記第２の振幅計算器と前記減算器との間に高
周波成分を除去する第２の高周波除去フィルタが設られ
ていることを特徴とする請求項１ないし４のいずれか１
項記載の周波数推定装置。
【請求項６】前記加算器の後段に、さらに低域通過フ
ィルタが設られていることを特徴とする請求項１ないし
５のいずれか１項記載の周波数推定装置。
【請求項７】前記係数器と前記加算器との間に出力制
限器が設られていることを特徴とする請求項１ないし６
のいずれか１項記載の周波数推定装置。
【請求項８】ｎ−１回目の周波数推定値を呼び出すス
テップと、前記ｎ−１回目の周波数推定値を第１の周波
数ω_１とし、所定の帯域除去あるいは帯域通過特性を備
えた第１のフィルタを設計するステップと、入力信号を
前記第１のフィルタでフィルタ処理するステップと、前
記第１のフィルタの出力から振幅を計算するステップ
と、前記第１の周波数ω_１よりも所定値だけ高い周波数
を第２の周波数ω_２とし、所定の帯域除去あるいは帯域
通過特性を備えた第２のフィルタを設計するステップ
と、前記入力信号を前記第２のフィルタでフィルタ処理
するステップと、前記第２のフィルタの出力から振幅を
計算するステップと、前記第１のフィルタの出力によっ
て算出された振幅から前記第２のフィルタの出力によっ
て算出された振幅を減算するステップと、前記減算結果
に所定の係数を乗算した値に前記ｎ−１回目の周波数推
定値を加算することによりｎ回目の周波数推定値を算出
するステップと、を含んでなり、前記減算結果が所定の
値以下になるまで上記ステップを繰り返すことを特徴と
する周波数推定方法。