WO2016006079A1

WO2016006079A1 - ピーク周波数検出装置、方法およびプログラム

Info

Publication number: WO2016006079A1
Application number: PCT/JP2014/068432
Authority: WO
Inventors: 一史阿久津
Original assignee: 一史阿久津
Priority date: 2014-07-10
Filing date: 2014-07-10
Publication date: 2016-01-14
Also published as: CN106415310A; KR101785714B1; JP5766896B1; US20170205448A1; EP3168638A1; EP3168638A4; US9857399B2; CN106415310B; KR20170009878A; JPWO2016006079A1

Abstract

【課題】周波数分解能（ｆ_０）と時間窓長（Ｔ_０）との間にある相反関係を回避し、所望の周波数分解能と所望の時間窓長でもって、信号波のピーク周波数の検出を可能にする。【解決手段】所定の周波数帯域（ｆ_ｃｌ～ｆ_ｃｈ）においてパワースペクトルが最大となるピーク周波数を検出するピーク周波数検出装置であって、デジタルデータ列の各要素をｎ乗（ｎは２以上の整数）するｎ乗部と、サンプリング周波数ｆ_ｓと周波数分解能ｆ_ｔｇと時間窓長Ｔ_ｔｇに応じて決まるＮ（Ｎは２のべき乗の整数）個のｎ乗されたサンプリング周波数ｆ_ｓのデジタルデータ列に対して高速フーリエ変換を行って得られるパワースペクトルの最大値に対応する周波数を仮想ピーク周波数として導出するＦＦＴ部と、前記仮想ピーク周波数を１／ｎ倍した値をデジタルデータ列のピーク周波数として出力する１／ｎ倍部と、を備え、ｎ≧１／（ｆ_ｔｇ×Ｔ_ｔｇ）、ｆ_ｓ／（ｎ×ｆ_ｔｇ）≦Ｎ≦ｆ_ｓ×Ｔ_ｔｇ、ｆ_ｓ＞２×ｎ×ｆ_ｃｈを満たすピーク周波数検出装置。

Description

ピーク周波数検出装置、方法およびプログラム

　本発明は、パワースペクトルが最大となるピーク周波数の検出装置、方法およびプログラムに関する。

　高速フーリエ変換（ＦＦＴ）の時間窓長をＴ_０（ｓ）、サンプリング周波数をｆ_ｓ（Ｈｚ）、サンプリング数をＮ（Ｎは２のべき乗の整数）とすると、

の関係がある。時間窓長の逆数ｆ_０（Ｈｚ）は周波数分解能と呼ばれ、

となる。

　信号波をＦＦＴで周波数解析をする場合、時間窓長はＦＦＴを行う信号波の時間長さ、周波数分解能は信号波周波数の最小分解能、つまり周波数の検出精度となる。以下、この時間窓長のことを、「ＦＦＴにかける信号波の切り出し時間」「ＦＦＴにかける信号波の時間」などともいう。［数２］式の意味するところは、ＦＦＴの時間窓長と周波数分解能は相反関係にあるということである。この相反関係は、周波数分解能が小さい周波数のところで影響が顕著になる。例えば、０．０１（Ｈｚ）の周波数分解能が必要であれば、その逆数の１００（ｓ）間の信号波をＦＦＴにかけなければならない。また、ＦＦＴにかける信号波の時間が０．０１（ｓ）であれば、その逆数の１００（Ｈｚ）の周波数分解能しか得られない。

　この相反関係は、ＦＦＴを利用するにあたって制約となることがあり、相反関係の制約を回避すべく様々な工夫がなされてきた。例えば、受信した信号波をＦＦＴにかけ、周波数領域に現れるスペクトルを構成するピーク周波数を求めるに際して、時間窓長Ｔ_０を、予め決めた値に固定し、受信した信号波からＴ_０間隔で相互に重複しないようにずらしながら切り出した複数の期間のそれぞれについてＦＦＴを適用し（即ち、複数回に亘ってＦＦＴを適用し）、求めた複数のピーク周波数を平均化する、という方法が考えられる。この方法では、１回のＦＦＴで求められる周波数分解能ｆ_０よりも高い分解能（細かい周波数間隔）でピーク周波数を計算できることもある。しかし、ＦＦＴの操作回数を増やせば、ＦＦＴにかける信号波の切り出し時間がその回数分長くなると言う問題がある。

　なお、デジタルデータ列にＦＦＴを適用すると、ＦＦＴに供したサンプル数Ｎに対し、０からＮ－１までの一定の周波数間隔（周波数分解能（ｆ_０）に等しい）での振幅の値が与えられるパワースペクトルが得られる。これらのパワースペクトルのうち、最大値のパワースペクトルに対応する周波数のことを、ピーク周波数と呼ぶ。最大値のパワースペクトルは、公知の方法により、例えば、１からＮ／２までのパワースペクトルを１個ずつ比較していけば特定できる。このとき、ｐの点でパワースペクトルが最大であれば、ピーク周波数ｆ_ｐｋはｆ_ｐｋ＝ｐ×ｆ_０となる。

　特許文献１には、周波数分解能が１２（Ｈｚ）以下、ＦＦＴにかける信号波の切り出し時間が１０（ｍｓ）（水中の位置分解能７．５ｍに相当）以下を同時に満たしてピーク周波数を求める方法が説明されている（同文献の［００２３］～［００２６］、［００８９］～［００９０］参照。なお［００９７］以降の数値例では、入力信号波からの切り出し時間が５（ｍｓ）に変わっているが、理由は不明である。）。周波数分解能１２（Ｈｚ）を優先すると、［数２］式から、切り出し時間は、１／１２（Ｈｚ）＝８３．３（ｍｓ）となる。一方、ＦＦＴを適用する場合にはそのサンプル数Ｎは２のべき乗であることが要請されるので、同文献では、信号切り出し時間を、Ｎ＝１０２４に相当する１０２．４（ｍｓ）としている（同文献の［００５５］参照）。ＦＦＴにかける入力信号波の切り出し時間を１０２．４（ｍｓ）とするためは、実際のデータが得られる時間５（ｍｓ）では不足するため、不足分の１０２．４－５＝９７．４（ｍｓ）間の部分にゼロ値データを付加し、このゼロ値データを加えて得られたデータ列（時間長１０２．４（ｍｓ）の部分）をＦＦＴにかけて、ピーク周波数を求めている（同文献の［０１０３］参照）。

　同文献記載の方法では、ＦＦＴにかけるデータ列の中に含まれる有効データ列は、全体の約１／１０（同文献［０１０３］によれば、約１／２０）に過ぎないので、入力信号波がアナログ信号波のときに一般に適用されるＡ／Ｄ変換をして得られるデジタル信号波にかけるデジタルバンドパスフィルタ（同文献図７では、デジタルＢＰＦ６２）を、次数を高くかつ狭帯域にすると、デジタルバンドパスフィルタ通過後の信号が鈍ってさらに有効データ数が減るため、デジタルバンドパスフィルタのフィルタ特性を緩やかにしなければならない。つまり、同文献には、時間長の短い入力信号波の場合の外乱ノイズによるネガティブな影響を回避することが困難となるという潜在的な課題は何ら認識されておらず、従って、そのような潜在的な課題解決の手段は何ら示唆も教示もなされていない。

　また、同文献の場合、ＦＦＴにかけるデジタル信号波（デジタルデータ列）の中には、ＦＦＴを適用する際に要請される時間窓長の約１／１０（もしくは、約１／２０）の長さの部分にしか実際のデータは存在しないために、多くのパワースペクトルが発生せざるを得ない。ＦＦＴを適用する実際の信号波に外乱ノイズが含まれている場合には、ピーク周波数の周辺にさらに余分なスペクトルが現れることになるので、ピーク周波数を特定することが困難になる可能性がある。

　このように、入力信号のピーク周波数を求めるにあたり、相反関係の制約を回避するために種々の試みがなされているが、一般的な解決方法には未だ至っていない。

特開２０１２－２４７３０２号公報

　本発明の目的は、ＦＦＴによる周波数解析の際に制約となる周波数分解能（ｆ_０）と時間窓長（Ｔ_０）との間にあるｆ_０＝１／Ｔ_０という相反関係を回避し、所望の周波数分解能と所望の時間窓長でもって、信号波のピーク周波数の検出を可能にする装置、方法およびプログラムを提供することにある。

　（１）上記目的を達成するためのピーク周波数検出装置は、所定の周波数帯域（ｆ_ｃｌ～ｆ_ｃｈ）においてパワースペクトルが最大となるピーク周波数を検出するピーク周波数検出装置であって、デジタルデータ列の各要素をｎ乗（ｎは２以上の整数）するｎ乗部と、サンプリング周波数ｆ_ｓと周波数分解能ｆ_ｔｇと時間窓長Ｔ_ｔｇに応じて決まるＮ（Ｎは２のべき乗の整数）個のｎ乗されたサンプリング周波数ｆ_ｓのデジタルデータ列に対して高速フーリエ変換を行って得られるパワースペクトルの最大値に対応する周波数を仮想ピーク周波数として導出するＦＦＴ部と、前記仮想ピーク周波数を１／ｎ倍した値をデジタルデータ列のピーク周波数として出力する１／ｎ倍部と、を備え、
ｎ≧１／（ｆ_ｔｇ×Ｔ_ｔｇ）
ｆ_ｓ／（ｎ×ｆ_ｔｇ）≦Ｎ≦ｆ_ｓ×Ｔ_ｔｇ
ｆ_ｓ＞２×ｎ×ｆ_ｃｈ
　を満たす。なお本明細書において、ｆ_ｃｌ～ｆ_ｃｈはｆ_ｃｌ以上ｆ_ｃｈ以下を意味するものとする。

　周波数解析の対象となる信号波ｙが、単一の周波数f（Ｈｚ）を有する[数３]式で表される正弦関数であると仮定すると、ｙのｎ乗（ｎは２以上の正の整数）は［数４］式および［数５］式で表される。［数４］式においてｎは奇数であって、［数５］においてｎは偶数である。

　［数４］式および［数５］式は、ｙをｎ乗（ｎ＝２,３,４,・・・）すると、ｎが奇数か偶数かに関わらず、周波数成分ｎ×ｆ，（ｎ－２）×ｆ，（ｎ－４）×ｆ，…が表れるということを示している。

　したがって、ｙ＝ｓｉｎ（２πｆｔ）をｎ乗して、そのｎ乗して得られた信号波から、周波数成分（ｎ×ｆ）の信号波を抽出し、その信号波のピーク周波数を求めることができれば、そこで求められた周波数を１／ｎ倍すれば、それが、元の信号波の周波数ｆに対応することが理解できる。図１は、この関係を模式的に示したものである。

　ここで、周波数解析に使用するＦＦＴのサンプリング周波数をｆ_ｓ、サンプル数をＮとする。ｙ＝ｓｉｎ（２πｆｔ）をこのＦＦＴにかけてピーク周波数を求めると、これがｆの計算値であり、その時の周波数分解能ｆ_０は、

となる。

　ここで、周波数成分が（ｎ×ｆ）の正弦関数ｙ_ｎ＝ｓｉｎ（２πｎｆｔ）について、ｙのピーク周波数を求めたときと同じ値のｆ_ｓとＮを有するＦＦＴにかけてピーク周波数を求めると、求めたｙ_ｎの周波数ｆ_ｎ'は、（ｎ×ｆ）の計算値となる。ここで注目すべきは、このｙ_ｎの周波数ｆ_ｎ'は、ｙについての周波数分解能と同じ周波数分解能ｆ_０のもとで求められた、ということである。つまり、ｆ_ｎ'の周波数分解能もｆ_０である。したがって、ｆ_ｎ'を１／ｎ倍したｆ_ｎ'／ｎは、ｙの周波数ｆを、周波数分解能ｆ_０／ｎで求めたものとなる。

　したがって、ｙの周波数ｆを直接ＦＦＴにかけて求めたときの周波数分解能をｆ_０、時間窓長をＴ_０とし、ｙの周波数ｆをｙ_ｎから求めたときの周波数分解能をｆ_ｎ、時間窓長をＴ_ｎとすると、

であり、

の関係が成立する。

　［数９］式について、ｎ＝２，３，４のときと、べき乗をしないときのｆ_０×Ｔ_０＝１のグラフを図２に示す。

　周波数分解能ｆ_ｔｇは、元の信号波のピーク周波数を求めるにあたって、ユーザが希望するピーク周波数の周波数分解能である。時間窓長Ｔ_ｔｇは、ユーザが希望するＦＦＴの時間窓長である。ｆ_ｔｇとＴ_ｔｇとは独立に定めることができる。

　つまり、本発明は、元の信号波のピーク周波数を、以下の条件を満たして計算する方法である。
・ピーク周波数の周波数分解能ｆ_ｎ≦ｆ_ｔｇ
・ＦＦＴの時間窓長Ｔ_ｎ≦Ｔ_ｔｇ

　このためには［数９］式から、

を満たす必要がある。

　つまり、乗数ｎは、

を満たさなければならない。

　なお、ｎは２のべき乗、２，４，８，１６，３２，・・・に制限されず、２以上の整数から選べばよい。したがって、用途に応じて、必要で十分な最適な値を選択することができる。すなわち、本発明によると、ＦＦＴによる周波数解析の際に制約となる周波数分解能（ｆ_０）と時間窓長（Ｔ_０）との間にあるｆ_０＝１／Ｔ_０という相反関係を回避し、所望の周波数分解能と所望の時間窓長でもって、信号波のピーク周波数を検出することが可能になる。具体的には、ｆ_ｔｇ×Ｔ_ｔｇ＜１の範囲においても所望の周波数分解能と所望の時間窓長で信号波のピーク周波数を検出することができる。

　信号波をｎ乗した後のデジタルデータ列をＦＦＴにかける場合のピーク周波数の周波数分解能ｆ_ｎ、時間窓長Ｔ_ｎには、[数９]式で示したとおりｆ_ｎ×Ｔ_ｎ＝１／ｎの関係があり、これを図示すると図４のようになる。この図４に描かれたｆ_ｎ×Ｔ_ｎ＝１／ｎのグラフにおける太線部分が、ｆ_ｔｇ、Ｔ_ｔｇの要求を満たすｆ_ｎ、Ｔ_ｎとなる。

　ｆ_ｎ＝ｆ_ｔｇのときのＴ_ｎをＴ_ｍｉｎとすると、Ｔ_ｎの取りうる範囲は、

となる。一方、

だから、

となる。

　ＦＦＴの場合、Ｎは２のべき乗でなければならないため、Ｎの取りうる範囲は、
ｆ_ｓ／（ｎ×ｆ_ｔｇ）≦Ｎ≦ｆ_ｓ×Ｔ_ｔｇ（Ｎは２のべき乗）
となる。

　もし、Ｎが存在しない場合は、Ｎが存在するまでｎを大きくするか、ｆ_ｓを大きくするか、もしくはその両方を行う。但し、ＦＦＴのサンプリング定理から、

を満たしていなければならない。

　該当するＮが複数存在する場合、どのＮを採用しても所望の周波数分解能ｆ_ｔｇと所望の時間窓長Ｔ_ｔｇを満たしている。但し、Ｎが小さい方がＦＦＴの計算量は少なくなるため、通常は、

に最も近いＮを採用するのが良い。

　このように、ｆ_ｔｇ×Ｔ_ｔｇ＜１の範囲においても、所望の周波数分解能ｆ_ｔｇ、時間窓長Ｔ_ｔｇでピーク周波数を検出できる乗数ｎ、ＦＦＴのサンプリング周波数ｆ_ｓ、ＦＦＴのサンプル数Ｎを決めることができる。なお、乗数ｎ、ｆ_ｓ、Ｎの決め方は、上記の手順によらなくても良い。適当に数値を当てはめて条件を満たすまで試行錯誤を繰り返したり、簡単なプログラムを作成して導出したり、どのような方法を用いても良い。

　（２）上記目的を達成するためのピーク周波数検出装置は、前記所定の周波数帯域に含まれる周波数のデジタルデータ列を抽出する第一デジタルバンドパスフィルタを備え、前記ｎ乗部には、前記第一デジタルバンドパスフィルタの出力が入力されてもよい。

　この構成を採用すると、様々な周波数成分を含む信号から所定の周波数帯域の信号を抽出するために前段のアナログフィルタを簡略化または省略できるため、回路規模を小さくすることができる。

　（３）上記目的を達成するためのピーク周波数検出装置は、サンプリング周波数ｆ_ｉｓのデジタルデータ列を１／ｒ（ｒは２以上の整数）に間引いてサンプリング周波数をｆ_ｓにする間引き部を備え、前記第一デジタルバンドパスフィルタには前記間引き部の出力が入力されてもよい。

　この構成を採用すると、ピーク周波数検出装置に入力されるデジタルデータ列のサンプリング周波数がｆ_ｓより大きい場合であっても、所望の周波数分解能ｆ_ｔｇ、時間窓長Ｔ_ｔｇでピーク周波数を検出できる。

　（４）上記目的を達成するためのピーク周波数検出装置は、受信デジタルデータ列をｇ倍（ｇは２以上の整数）に補間してサンプリング周波数をｆ_ｓにする補間部を備え、前記第一デジタルバンドパスフィルタには前記補間部の出力が入力されてもよい。

　この構成を採用すると、ピーク周波数検出装置に入力されるデジタルデータ列のサンプリング周波数がｆ_ｓより小さい場合であっても、所望の周波数分解能ｆ_ｔｇ、時間窓長Ｔ_ｔｇでピーク周波数を検出できる。

　（５）上記目的を達成するためのピーク周波数検出装置は、ｎ乗されたＮ個の前記デジタルデータ列から第二の周波数帯域に含まれるデジタルデータ列を抽出する第二デジタルバンドパスフィルタを備え、前記ＦＦＴ部には、前記第二デジタルバンドパスフィルタによって抽出されたデジタルデータ列が入力され、前記第二の周波数帯域は、ほぼｎ×ｆ_ｃｌ～ｎ×ｆ_ｃｈであってもよい。

　仮に単一の周波数成分しかもたないデジタルデータ列であっても、ｎ乗すると複数の周波数成分を持つことになるため、パワースペクトルに複数のピークが現れるが、検出しようとするピークは、帯域ｆ_ｃｌ～ｆ_ｃｈに対応するｎ×ｆ_ｃｌ～ｎ×ｆ_ｃｈに現れる。したがって、ｎ×ｆ_ｃｌ～ｎ×ｆ_ｃｈの成分を抽出することにより、所定の周波数帯域（ｆ_ｃｌ～ｆ_ｃｈ）においてパワースペクトルが最大となるピーク周波数を検出することが可能になる。なお、複数のピーク周波数から検出しようとするピーク周波数を特定して選択できる場合には、第二デジタルバンドパスフィルタを用いなくともよい。

　（６）上記目的を達成するためのピーク周波数検出装置は、前記第二デジタルバンドパスフィルタで抽出されたデジタルデータ列を１／ｒ（ｒは２以上の整数）に間引いてサンプリング周波数をｆ_ｓとする間引き部を備え、前記ＦＦＴ部には、前記間引き部の出力が入力されてもよい。

　（７）上記目的を達成するためのピーク周波数検出装置において、前記特定の周波数帯域に含まれる周波数のデジタルデータ列を抽出する第一デジタルバンドパスフィルタと、前記ｎ乗部の出力から第二の周波数帯域に含まれるデジタルデータ列を抽出する第二デジタルバンドパスフィルタとを備え、前記ｎ乗部には、前記第一デジタルバンドパスフィルタの出力が入力され、前記ＦＦＴ部には、前記第二デジタルバンドパスフィルタの出力が入力され、前記第二の周波数帯域は、ほぼｎ×ｆ_ｃｌ～ｎ×ｆ_ｃｈであってもよい。

　（８）上記目的を達成するためのピーク周波数検出装置は、前記ｎ乗部に代えて、入力されるデジタルデータ列をｍ_ｊ乗するべき乗部（ｊ）（ｊ＝１，２，…，ｋ）と前記べき乗部（ｊ）の出力から特定の周波数帯域ｆ_ｃｌ（ｊ）～ｆ_ｃｈ（ｊ）の信号を抽出するデジタルバンドパスフィルタ（ｊ）とを備えるべき乗ブロック（ｊ）をｋ（ｋは２以上の整数）段備える多重べき乗部を備え、
ｎ＝ｍ_１×ｍ_２×…×ｍ_ｋ
ｆ_ｃｌ（ｊ）≒（ｍ_１×ｍ_２×…×ｍ_ｊ）×ｆ_ｃｌ
ｆ_ｃｈ（ｊ）≒ｍ_１×ｍ_２×…×ｍ_ｊ）×ｆ_ｃｈ
であってもよい。

　デジタルデータ列をべき乗する場合、乗数が小さいほど、べき乗することによって発生する余分な周波数成分を減らすことができる。余分な周波数成分が少ないほど、デジタルバンドパスフィルタで余分な周波数成分を除去しやすくなる。すなわち、１段のべき乗部と第二デジタルバンドパスフィルタで対応できる乗数ｎの上限よりも、べき乗部とデジタルバンドパスフィルタをユニット化して多段構成にして対応できる乗数ｎの上限の方が大きい。したがって、この構成を採用することにより、乗数として設定できるｎの範囲を広くすることができる。

　（９）上記目的を達成するためのピーク周波数検出装置は、ユーザの指示を受け付ける操作部と、前記指示に応じたｎ、ｆ_ｓ、Ｎの少なくともいずれかを設定するパラメータ設定部と、を備えてもよい。

　この構成を採用することにより、次式を満たすｎ、ｆ_ｓ、Ｎをユーザが試行錯誤を繰り返して設定する必要がなくなる。
ｎ≧１／（ｆ_ｔｇ×Ｔ_ｔｇ）
ｆ_ｓ／（ｎ×ｆ_ｔｇ）≦Ｎ≦ｆ_ｓ×Ｔ_ｔｇ
ｆ_ｓ＞２×ｎ×ｆ_ｃｈ

　なお請求項に記載された各手段の機能は、構成自体で機能が特定されるハードウェア資源、プログラムにより機能が特定されるハードウェア資源、又はそれらの組み合わせにより実現される。また、これら各手段の機能は、各々が物理的に互いに独立したハードウェア資源で実現されるものに限定されない。さらに、本発明は、方法としても、コンピュータプログラムとしても、コンピュータプログラムの記録媒体としても成立する。むろん、そのコンピュータプログラムの記録媒体は、磁気記録媒体であってもよいし光磁気記録媒体であってもよいし、今後開発されるいかなる記録媒体であってもよい。

本発明の実施形態にかかるスペクトルを示す図。本発明の実施形態にかかるグラフ。本発明の実施形態にかかるブロック図。本発明の実施形態にかかるグラフ。本発明の実施形態にかかるスペクトルを示す図。本発明の実施形態にかかるブロック図。本発明の実施形態にかかる波形図。本発明の実施形態にかかるブロック図。本発明の実施形態にかかるブロック図。本発明の実施形態にかかるブロック図。本発明の実施形態にかかる波形図。本発明の実施形態にかかる波形図。本発明の実施形態にかかるブロック図。本発明の実施形態にかかるブロック図。本発明の実施形態にかかるブロック図。本発明の実施形態にかかる画面構成図。

　以下、本発明の実施の形態を添付図面を参照しながら説明する。尚、各図において対応する構成要素には同一の符号が付され、重複する説明は省略される。
　１．第一実施形態
　第一実施形態では、サンプリング周波数ｆ_ｓでサンプリングされた受信デジタルデータ列の、下限値ｆ_ｃｌと上限値ｆ_ｃｈとで定められた周波数帯域ｆ_ｃｌ～ｆ_ｃｈにおけるピーク周波数ｆを、所望の周波数分解能ｆ_ｔｇと所望の時間窓長Ｔ_ｔｇを満たして検出するピーク周波数検出装置とこれを用いたドップラー計測器について説明する。

　本発明の第一実施形態としてのピーク周波数検出装置１は、図３に示すように、第一デジタルバンドパスフィルタ（ＢＰＦ）部１１と、ｎ乗部１２と、第二デジタルバンドパスフィルタ（ＢＰＦ）部１３と、ＦＦＴ部１４と、１／ｎ倍部１５とを具備している。

　はじめに乗数ｎと、周波数解析に使用するＦＦＴのサンプリング周波数ｆ_ｓと、周波数解析に使用するＦＦＴのサンプリング数Ｎの決め方の手順を説明する。
　（ステップ１．乗数ｎを定める）
　まず、ｎ≧１／（ｆ_ｔｇ×Ｔ_ｔｇ）の条件を満たすｎ（ｎは、２以上の整数）を定める。例えば、ｎを、ｎ≧１／（ｆ_ｔｇ×Ｔ_ｔｇ）を満たす最小の整数とし、後述するステップ３でＮが存在しなければ、ｎを１増加させて再計算をしても良いし、最初から大きめにしておいても良い。
　（ステップ２．ＦＦＴのサンプリング周波数ｆ_ｓを選択する）

　ＦＦＴのサンプリング周波数ｆ_ｓは、サンプリング定理を満たすため、
ｆ_ｓ＞２×ｎ×ｆ_ｃｈ
となるように選択する。なお、本実施形態では、ピーク周波数検出装置１に入力されるデジタルデータ列のサンプリング周波数がＦＦＴのサンプリング周波数ｆ_ｓとなる。ｆ_ｓが回路上の制約等で上式を満たせない場合は、他の実施形態を使用する。
　（ステップ３．ＦＦＴのサンプル数Ｎを定める）
　次に、ｆ_ｓ／（ｎ×ｆ_ｔｇ）≦Ｎ≦ｆ_ｓ×Ｔ_ｔｇ（Ｎは２のべき乗）を満たすＮを選択する。もし、Ｎが存在しない場合は、Ｎが存在するまでｎを大きくするか、ｆ_ｓを大きくするか、もしくはその両方を行う。ただし、ＦＦＴのサンプリング定理から、
ｆ_ｓ＞２×ｎ×ｆ_ｃｈ
を満たしていなければならない。

　以上説明した方法で選択される乗数ｎ、ｆ_ｓ、Ｎは、ｎ乗部１２およびＦＦＴ部１４に設定される。また、ピーク周波数を検出しようとする帯域の下限値ｆ_ｃｌと上限値ｆ_ｃｈは第一デジタルＢＰＦ１１のカットオフ周波数として設定される。また第二デジタルＢＰＦ１３のカットオフ周波数には、ピーク周波数を検出しようとする帯域の下限値ｆ_ｃｌと上限値ｆ_ｃｈのｎ倍の値が設定される。

　（デジタルデータ列からピーク周波数を計算する）
　乗数ｎ、ｆ_ｓ、Ｎが設定されたピーク周波数検出装置１に対象とするデジタルデータ列を入力すると、以下のように、所望の周波数分解能ｆ_ｔｇと所望の時間窓長Ｔ_ｔｇを満たしてピーク周波数が検出される。

　対象とするデジタルデータ列がピーク周波数検出装置１に入力されると、第一デジタルＢＰＦ１１は、デジタルデータ列の幅広い周波数成分から、設定されている帯域外の余分な直流成分、低周波成分、高周波成分を排除し、単一の周波数ｆに近い周波数成分を抽出する。

　ピーク周波数を検出しようとする帯域としてｆ_ｃｌ～ｆ_ｃｈが広すぎて、その間に大きなパワースペクトルが複数存在する場合、ｆ_ｃｌ～ｆ_ｃｈを狭く設定し直して、図５のように、何回かに分けてピーク周波数を求めるのが良い。図５の例では、ｆ_ｃｌ～ｆ_ｃｈが最初の周波数帯域、ｆ_１、ｆ_２が、大きなパワースペクトルを持つ２つの周波数、ｆ_ｃｌ'～ｆ_ｃｈ'が狭めた周波数帯域を示す。なお、ｆ_ｃｌ'～ｆ_ｃｈ'の間隔は一定である必要はなく、図５に示すように、周波数が高くなるほど、ｆ_ｃｌ'～ｆ_ｃｈ'の間隔を広げてもよい。

　ここで第一デジタルＢＰＦ１１を通過した後のデジタルデータ列を、Ａ（１）：ａ_０，ａ_１，ａ_２，…とする。Ａ（１）には、周波数ｆのほかに、パワースペクトルの小さい若干の余分な周波数成分が含まれる。

　次に、このＡ（１）をｎ乗部１２に入力する。ｎ乗部１２では、Ａ（１）の各要素を、上で定めたｎを乗数としてべき乗する。ｎ乗部１２通過後のデジタルデータ列をＢ（ｎ）：ｂ_０，ｂ_１，ｂ_２，…とすると、ｂ_ｉ＝（ａ_ｉ）^ｎ（ｉ＝０，１，２，３，４，…）となる。Ｂ（ｎ）には、周波数成分、ｎ×ｆ,（ｎ－２）×ｆ,（ｎ－４）×ｆ，…,とパワースペクトルの小さい、余分な低周波成分と高周波成分が含まれる。

　次に、Ｂ（ｎ）を第二デジタルＢＰＦ１３に入力し、ｎ×ｆ_ｃｌ～ｎ×ｆ_ｃｈの帯域の周波数成分を抽出する。第二デジタルＢＰＦ１３によって、Ｂ（ｎ）の幅広い周波数成分から、余分な低周波、高周波成分を排除し、単一の周波数に近い周波数成分を抽出することができる。第二デジタルＢＰＦ１３の周波数帯域は、第一デジタルＢＰＦ１１の周波数帯域のｎ倍のｎ×ｆ_ｃｌ～ｎ×ｆ_ｃｈとするのが良いが、用途によっては、多少変えてもかまわない。第二デジタルＢＰＦ１３を通った後のデジタルデータ列を、Ｃ（ｎ）：ｃ_０，ｃ_１，ｃ_２，…とする。Ｃ（ｎ）には、周波数ｎ×ｆのほかに、若干の余分な周波数成分が含まれる。

　ここで、ピーク周波数検出装置１に入力されるデジタルデータ列のサンプリング周波数ｆ_ｓを、ｆ_ｓ＞２×ｎ×ｆ_ｃｈとした理由は、Ｃ（ｎ）の周波数帯域の上限がｎ×ｆ_ｃｈだからである。つまり、受信デジタルデータ列のサンプリング周波数ｆ_ｓは、ＦＦＴのサンプリング定理を満たすことを条件とした。

　次に、Ｃ（ｎ）を、ＦＦＴ部１４に入力してピーク周波数を算出する。ＦＦＴ部１４では、Ｃ（ｎ）に対し、設定されたサンプリング周波数ｆ_ｓとサンプル数ＮでＦＦＴを実行し、ピーク周波数を算出する。ＦＦＴ部１４から出力されるピーク周波数を（ｎ×ｆ'）とする。

　次に、ＦＦＴ部１４から出力されるピーク周波数（ｎ×ｆ'）を１／ｎ倍部１５に通し、ｆ'を求める。求めたｆ'が、ピーク周波数検出装置１に入力されたデジタルデータ列のピーク周波数ｆの計算値となる。このｆ'は、ｆ_ｔｇ、Ｔ_ｔｇを満たして算出されたものである。

　以上説明したとおり、ピーク周波数検出装置１は、受信デジタルデータ列から、所望の周波数分解能ｆ_ｔｇと所望の時間窓長Ｔ_ｔｇでピーク周波数を検出することができる。そしてピーク周波数検出装置１は、近似計算を行わないため、算出されるピーク周波数の精度劣化がない。したがって、高精度にピーク周波数を検出することができる。

　図６は、ピーク周波数検出装置１を組み込んだドップラー計測器２を示すブロック図である。ドップラー計測器２には、媒体中に信号波を発信する機能と媒体中の対象物からの反射波を受信する機能とを兼ね備えている送受信器２１が具備されている。送受信器２１は、送受信切替回路２２を介し、送信回路２３の出力と受信増幅器２４の入力に接続している。送信回路２３は、送信周波数ｆ_ｔｘの信号を生成するものである。受信増幅器２４の出力には、受信信号を中間周波数信号に周波数変換する変調器２５が設けられており、変調器２５には、局部発振回路２６から局部発振周波数ｆ_ｌｏｃの信号が供給されている。変調器２５の出力すなわち中間周波数信号は、アナログフィルタ２７を介してアナログ／デジタル変換器（Ａ／Ｄ）２８に入力されてＦＦＴに要請されるサンプリング周波数でデジタル信号に変換される。次にピーク周波数検出装置１に入力されてＦＦＴにより周波数解析され、ピーク周波数検出装置１から、ピーク周波数が出力される。

　本実施例では、以下の条件下で海水中の音波のピーク周波数を検出する方法を説明する。　海水中の音波伝搬速度Ｃ：１５００ｍ／ｓ
　送信周波数ｆ_ｔｘ：１２０ｋＨｚ
　局部発信周波数ｆ_ｌｏｃ：１３７ｋＨｚ
　検出最大速度（水平方向）Ｖ：１５ｍ／ｓ
　検出速度精度（水平方向）Ｖ_０：０．１５ｍ／ｓ
　検出対象物の位置の精度（位置分解能）Ｄ_０：７．５ｍ

　また、音波は水平方向から斜め方向（θ＝６０°）に送受信するものとする。Ｃ≫Ｖの場合、ドップラー周波数ｆ_ｄｏｐは、

となる。検出最大速度Ｖが１５ｍ／ｓだから、

となり、ドップラー信号は１２０±１．２ｋＨｚの範囲であり、観測される周波数帯域の幅Δｆ_ｐは、
　Δｆ_ｐ＝２×１２００Ｈｚ＝２４００Ｈｚ
となる。

　ここで中間周波数ｆ_ｍｉｄは、ｆ_ｌｏｃ－ｆ_ｔｘ＝１３７－１２０＝１７ｋＨｚとする。アナログフィルタ２７は、
　ｆ_ｍｉｄ±（△ｆ_ｐ／２）＝１７０００±１２００Ｈｚ
の信号を通し、次段のＡ／Ｄ変換器２８でエイリアシングが発生しないような特性とする。

　検出速度精度Ｖ_０が０．１５ｍ／ｓだから、

となり、周波数分解能ｆ_０は１２Ｈｚとなる。

　検出対象物の位置の精度（位置分解能）は７．５ｍだから、時間窓長は、７．５ｍを往復する時間である、７．５×２／１５００＝１０ｍｓとなる。

　この条件は、アナログ／デジタル変換器（Ａ／Ｄ）２８の出力をそのままＦＦＴにかけても、満たすことができない。なぜなら、周波数分解能が１２Ｈｚのときには、時間窓長は、１／１２＝８３．３ｍｓ（＞１０ｍｓ）となり、位置分解能は０．０８３３×１５００／２＝６２．５ｍ（＞７．５ｍ）となってしまう。逆に、時間窓長を１０ｍｓとすると、周波数分解能は、１／０．０１＝１００Ｈｚ（＞１２Ｈｚ）となってしまう。

　そこで、ｎ乗部１２と１／ｎ倍部１５とを備えるピーク周波数検出装置１によってピーク周波数を求める。所望の周波数分解能ｆ_ｔｇを１２Ｈｚ、所望の時間窓長Ｔ_ｔｇを１０ｍｓとすると、乗数ｎは、
ｎ≧１／（ｆ_ｔｇ×Ｔ_ｔｇ）＝１／（１２×０．０１）＝８．３
となる。ここでは、ｎ＝１２とする。

　次に、Ａ／Ｄ変換器２８のサンプリング周波数ｆ_ｓを決定する。第一デジタルＢＰＦ１１の周波数帯域ｆ_ｃｌ～ｆ_ｃｈは、１７０００±１２００Ｈｚに合わせると、ｆ_ｃｌ＝１５８００Ｈｚ、ｆ_ｃｈ＝１８２００Ｈｚとなる。したがって、Ａ／Ｄ変換器２８のサンプリング周波数ｆ_ｓは、
ｆ_ｓ＞２×ｎ×ｆ_ｃｈ＝２×１２×１８２００＝４３６８００Ｈｚ
となるため、ｆ_ｓ＝５１０ｋＨｚとする。

　ＦＦＴのサンプル数Ｎは、ｆ_ｓ／（ｎ×ｆ_ｔｇ）≦Ｎ≦ｆ_ｓ×Ｔ_ｔｇより、５１００００Ｈｚ／（１２×１２Ｈｚ）＝３５４１．７≦Ｎ≦５１００００Ｈｚ×０．０１ｓ＝５１００を満たす、２のべき乗の整数となる。すなわち、
Ｎ＝４０９６
となる。

　このようにしてサンプリング周波数ｆｓとサンプル数Ｎを定め、アナログ／デジタル変換器（Ａ／Ｄ）２８から出力されるデジタルデータ列をピーク周波数検出装置１に入力する。

　なお、第一デジタルＢＰＦ１１には、例えば、次数が８次でカットオフ周波数をｆ_ｃｌ＝１５.８ｋＨｚ、ｆ_ｃｈ＝１８.２ｋＨｚに設定したバタワース型のＩＩＲ（Infinite impulse response）フィルタを好適に使用できる。また、第二デジタルＢＰＦ１３には、例えば、次数が８次でカットオフ周波数を１８９.６ｋＨｚ（１２×ｆ_ｃｌ）および２１８.４ｋＨｚ（１２×ｆ_ｃｈ）に設定したバタワース型のＩＩＲフィルタを好適に使用できる。

　ここで、対象物の相対移動速度（水平方向）ＶをＶ＝１０ｍ／ｓとし、音波を水平方向から斜め方向（α＝６０°）に送受信するものとする。
　受信信号波のドップラー周波数ｆ_ｄｏｐは、

となる。

　Ｓｉｎ｛２π（１７０００＋８００）ｔ｝＝Ｓｉｎ（２π１７８００ｔ）をサンプリング周波数５１０ｋＨｚのＡ／Ｄ変換器２８への入力信号と仮定し、擬似的なデジタルデータ列を作成した上で、上記の設定でピーク周波数を実際に求めると、（１２×ｆ'）≒２１３６６２．１Ｈｚが得られた。求めるｆ'は、ｆ'＝（１２×ｆ'）／１２＝２１３６６２．１／１２≒１７８０５．２Ｈｚとなる。ｆ'は、ｆ_ｔｇを満たして計算されたものである。実際に、このデジタルデータ列に対する計測誤差εは、ε＝ｆ'－ｆ_ｄｏｐ＝１７８０５．２－１７８００＝５．２Ｈｚであり、±（ｆ_ｔｇ）／２＝±６Ｈｚ以下に収まっている。なお、周波数分解能ｆ_０のときのピーク周波数は、ｆ_０の間隔で最も近い点となるため、周波数分解能ｆ_０のときのピーク周波数の誤差は、±（ｆ_０／２）以下となる。

　本実施例において、第二デジタルＢＰＦ１３から出力されるデジタルデータ列の実例を図７に示す。図７において、Ｔ_１は、ＦＦＴにかける区間である。Ｔ_２は、デジタルデータ列の振幅が落ち着くまでの区間である。Ｔ_１の要素数は４０９６個であるので、Ｔ_１の区間長は、４０９６／５１００００≒８．０ｍｓとなり、所望の時間窓長Ｔ_ｔｇ＝１０ｍｓを満たしている。

　なお、区間Ｔ_２で振幅が落ち着いた後のデジタルデータ列をＦＦＴにかける方がより好ましい。Ｔ_２の要素数は約６００個であるので、Ｔ_２の区間長は６００／５１００００≒１．２ｍｓとなる。その場合でも、乗数ｎを大き目の１２に設定したためピーク周波数を計算するために必要なデジタルデータ列の時間長は、１．２＋８．０＝９．２ｍｓであり、Ｔ_ｔｇ＝１０ｍｓ以下となっている。なお、本実施例では、Ｔ_２以降の４０９６個のデジタルデータ列をＦＦＴにかけても、計算結果は、ｆ'≒１７８０５．２Ｈｚとなり、同じ値が得られた。

　このように、乗数ｎを上げることにより、余裕をもってＴ_ｔｇを満たすことができるようになり、第一デジタルＢＰＦ１１、第二デジタルＢＰＦ１３の次数を上げて急峻にし、外来ノイズに強くすることができる。なお、Ｔ_ｔｇが５ｍｓであっても、乗数ｎを上げることにより、所望の周波数分解能ｆ_ｔｇでもって周波数解析が可能であることは容易に理解できる。

　ここで、帯域通過フィルタである第一デジタルＢＰＦ１１を設けた理由について説明する。［数４］式および［数５］式を適用するためには、ｎ乗をする前に、周波数帯域ｆ_ｃｌ～ｆ_ｃｈ以外の余分な直流分、低周波分、高周波分を除いておく必要がある。これをアナログフィルタ２７で実現できれば、第一デジタルＢＰＦ１１は不要にできる。しかしこのようなアナログフィルタは、高次で高精度で回路規模が大きくなりコストがかかる。したがって、アナログフィルタ２７は、Ａ／Ｄ変換器２８でエイリアシングが発生しない程度の設計とし、周波数帯域ｆ_ｃｌ～ｆ_ｃｈの成分の抽出は、高次で高精度の設計が容易にできコストもかからないデジタルバンドパスフィルタを採用するのが得策である。

　２．第二実施形態
　図８は、本発明の第二実施形態としてのピーク周波数検出装置３の構成を示すブロック図である。ピーク周波数検出装置３は、ピーク周波数検出装置１の第一デジタルＢＰＦ１１の前に、間引き部１６が追加された構成となっている。これは、Ａ／Ｄ変換後のデジタルデータ列を間引いて減らすことによりサンプリング周波数を下げるためである。ピーク周波数検出装置３に入力されるデジタルデータ列のサンプリング周波数ｆ_ｓが高く、乗数ｎと第一デジタルＢＰＦ１１のカットオフ周波数ｆ_ｃｈがｆ_ｓ＞４×ｎ×ｆ_ｃｈを満たす場合、第二実施形態を適用することが好ましい。

　間引き部１６では、デジタルデータ列を１／ｒ（ｒは２以上の整数）に間引き、間引いた後のデジタルデータ列が、
ｆ_ｓ＞２×ｎ×ｆ_ｃｈ
となるようにする。ｆ_ｓは間引き後のサンプリング周波数である。

　１／ｒに間引くには、例えば、次のようにする。間引く前のデジタルデータ列を、
Ｐ（１）：ｐ_０，ｐ_１，ｐ_２，…とし、間引いた後のデジタルデータ列を、Ｑ（１）：ｑ_０，ｑ_１，ｑ_２，…としたとき、ｑ_ｉ＝ｐ_{（ｒ×ｉ）}（ｉ＝０，１，２，３，４，・・・）を行う。ｒ＝２のときは、ｑ_０＝ｐ_０、ｑ_１＝ｐ_２、ｑ_２＝ｐ_４、ｑ_３＝ｐ_６となる。
　間引く方法はこれ以外であっても良い。例えば、ｒ＝２のときに、ｑ_０＝（ｐ_０＋ｐ_１）／２、ｑ_１＝（ｐ_２＋ｐ_３）／２、ｑ_２＝（ｐ_４＋ｐ_５）／２、…などとしても良い。

　なお、ｒ＝２としたときに、
ｆ_ｓ／（ｎ×ｆ_ｔｇ）≦Ｎ≦ｆ_ｓ×Ｔ_ｔｇ（Ｎは２のべき乗）
を満たすＮが存在せず、ｎを大きくすることもできない場合、本実施例は適用できない。この場合は第一実施形態を適用する。

　以下、図６に示したドップラー計測器２に本実施例のピーク周波数検出装置３を適用した場合について説明する。なお、下記の数値データ以外の必要な数値データは、第一実施形態の数値例と同じとする。また、後の説明を簡略化するため、第一実施形態同様、ｆ_ｓ＝５１０ｋＨｚとなるようＡ／Ｄ変換器２８のサンプリング周波数を１０．２ＭＨｚとすると、間引き部１６に入力される受信デジタルデータ列のサンプリング周波数ｆ_ｉｓは、１０．２ＭＨｚとなる。

　ここでｆ_ｉｓ＝１０．２ＭＨｚ＞４×ｎ×ｆ_ｃｈ＝４×１２×１８．２ｋＨｚ＝８７３．６ｋＨｚである。一方、２×ｎ×ｆ_ｃｈ＝２×１２×１８．２ｋＨｚ＝４３６．８ｋＨｚである。そこで、受信デジタルデータを、１／ｒ＝１／２０に間引くものとする。間引いた後のデジタルデータ列のサンプリング周波数をｆ_ｓとすると、
ｆ_ｓ＝ｆ_ｉｓ／ｒ＝１０．２ＭＨｚ／２０＝５１０ｋＨｚ＞４３６．８ｋＨｚ
となる。

　間引く前のデジタルデータ列を、Ｐ（１）：ｐ_０，ｐ_１，ｐ_２，…とし、間引いた後のデジタルデータ列を、Ｑ（１）：ｑ_０，ｑ_１，ｑ_２，…とすると、間引き部１６では、ｑ_０＝ｐ_０、ｑ_１＝ｐ_２０、ｑ_２＝ｐ_４０、ｑ_３＝ｐ_６０、…を行うこととなる。

　このＱ（１）は、入力信号を、サンプリング周波数ｆ_ｓ＝５１０ｋＨｚでＡ／Ｄ変換した信号列と同等であることが理解できる。従って、ＦＦＴのサンプリング周波数ｆ_ｓは、ｆ_ｓ＝５１０ｋＨｚとなる。

　第一実施形態と同様、ｙ＝ｓｉｎ（２π１７８００ｔ）をサンプリング周波数１０．２ＭＨｚのＡ／Ｄ変換器２８への入力信号と仮定し、擬似的なデータ列を作成した上で、上記の実施例の設定で実際に周波数解析を行なった結果は、第一実施形態の数値例と同じとなり、ｆ_ｔｇ、Ｔ_ｔｇを満たす。

　３．第三実施形態
　図９は、本発明の第三実施形態を示すピーク周波数検出装置４の構成を示すブロック図である。ピーク周波数検出装置４は、第一実施形態の第二デジタルＢＰＦ１３の後に間引き部１７が追加された構成となっている。第一デジタルＢＰＦ部１１、ｎ乗部１２、第二デジタルＢＰＦ１３は、入力されるデジタルデータのサンプリング周波数が、第一実施形態ではｆ_ｓだったものが、ｆ_ｉｓとなっただけで、動作は同じである。受信デジタルデータ列のサンプリング周波数ｆ_ｓが高く、乗数ｎと第一デジタルＢＰＦ１１のカットオフ周波数ｆ_ｃｈに対して、ｆ_ｓ＞４×ｎ×ｆ_ｃｈの場合、計算量を削減するため、一般的には第二実施形態を適用するのが良い。しかし、本実施形態であっても、所望の時間窓長ｆ_ｔｇ、所望の時間窓長Ｔ_ｔｇを満たした処理が可能である。

　間引き部１７では、受信デジタルデータ列を１／ｒ（ｒは２以上の整数）に間引いて、間引いた後のデジタルデータ列のサンプリング周波数をｆ_ｓとしたとき、
ｆ_ｓ＞２×ｎ×ｆ_ｃｈ
となるようにする。間引いた後のデジタルデータ列は、ＦＦＴ部１４に入力され、その後の処置は、第一実施形態と同一である。

　なお、ｒ≧２でｆ_ｓ／（ｎ×ｆ_ｔｇ）≦Ｎ≦ｆ_ｓ×Ｔ_ｔｇ（Ｎは２のべき乗）を満たせない場合には本実施形態は適用できない。この場合は、第一実施形態を適用する。

　以下、図６に示したドップラー計測器２に本実施例のピーク周波数検出装置４を適用した場合について説明する。なお、下記の数値データ以外の必要な数値データは、第一実施形態の数値例と同じとする。Ａ／Ｄ変換器２８のサンプリング周波数が１０．２ＭＨｚである場合、間引き部１７に入力される受信デジタルデータ列のサンプリング周波数ｆ_ｉｓは、１０．２ＭＨｚとなる。

　ここでｆ_ｉｓ＝１０．２ＭＨｚ＞４×ｎ×ｆ_ｃｈ＝４×１２×１８．２ｋＨｚ＝８７３．６ｋＨｚである。一方、２×ｎ×ｆ_ｃｈ＝２×１２×１８．２ｋＨｚ＝４３６．８ｋＨｚである。そこで、第二デジタルＢＰＦ１３通過後のデジタルデータを、１／ｒ＝１／２０に間引くものとする。間引いた後のデジタルデータ列のサンプリング周波数をｆ_ｓとすると、
ｆ_ｓ＝ｆ_ｉｓ／ｒ＝１０．２ＭＨｚ／２０＝５１０ｋＨｚ＞４３６．８ｋＨｚ
となる。ＦＦＴのサンプリング周波数ｆ_ｓは、ｆ_ｓ＝５１０ｋＨｚとなる。

　第一実施形態の数値例と同様、ｙ＝ｓｉｎ（２π１７８００ｔ）をサンプリング周波数１０．２ＭＨｚのＡ／Ｄ変換器２８への入力信号と仮定し、擬似的なデジタルデータ列作成した上で、上記の実施例の設定で実際に周波数解析を行なった結果は、第一実施形態の数値例と同じとなり、ｆ_ｔｇ、Ｔ_ｔｇを満たす。

　４．第四実施形態
（以下、出願時に段落番号を繰り上げます）　図１０は、本発明の第四実施形態としてのピーク周波数検出装置５の構成を示すブロック図である。ピーク周波数検出装置５は、第一実施形態のデジタルＢＰＦ１１の前に、補間部１８が追加された構成となっている。補間部１８を追加することにより、Ａ／Ｄ変換後のデジタルデータ列を補間して増やすことによりサンプリング周波数を上げることができる。ピーク周波数検出装置５に入力されるデジタルデータ列のサンプリング周波数がｆ_ｓが低く、乗数ｎと第一デジタルＢＰＦ１１の高いほうのカットオフ周波数ｆ_ｃｈに対して、２×ｆ_ｃｈ＜ｆ_ｓ＜２×ｎ×ｆ_ｃｈの場合、本実施形態を適用する。

　補間部１８では、受信デジタルデータ列をｇ倍（ｇは２以上の整数）に補間して、補間した後のデジタルデータ列のサンプリング周波数ｆ_ｓが、ｆ_ｓ＞２×ｎ×ｆ_ｃｈとなるようにする。例えば、補間前のデジタルデータ列を、Ｕ（１）：ｕ_０，ｕ_１，ｕ_２，…とし、補間した後のデジタルデータ列を、Ｖ（１）：ｖ_０，ｖ_１，ｖ_２，…としたとき、補間部１８は、
ｖ_ｉ＝ｕ_０　（ｉ＝０，１，２，３，・・・，（ｇ－１））
ｖ_ｉ＝ｕ_１　（ｉ＝ｇ，ｇ＋１，ｇ＋２，ｇ＋３，・・・，（２ｇ－１））
ｖ_ｉ＝ｕ_２　（ｉ＝２ｇ，２ｇ＋１，２ｇ＋２，２ｇ＋３，・・・，（３ｇ－１））
・・・
となるように補間する。ｇ＝２のときは、ｖ_０＝ｕ_０、ｖ_１＝ｕ_０、ｖ_２＝ｕ_１、ｖ_３＝ｕ_１、ｖ_４＝ｕ_２、ｖ_５＝ｕ_２、…となる。補間する方法はこれ以外であっても良い。例えば、ｇ＝２のときに、
ｖ_０＝ｕ_０
ｖ_１＝（ｕ_０＋ｕ_１）／２
ｖ_２＝ｕ_１
ｖ_３＝（ｕ_１＋ｕ_２）／２
ｖ_４＝ｕ_２
ｖ_５＝（ｕ_２＋ｕ_３）／２
・・・
となるように補間しても良い。

　以下、図６に示したドップラー計測器２に本実施例のピーク周波数検出装置５を適用した場合について説明する。なお、下記の数値データ以外の必要な数値データは、第一実施形態の数値例と同じとする。説明を簡略化するため第一実施形態同様、ｆ_ｓ＝５１０ｋＨｚとなるようにＡ／Ｄ変換器２８のサンプリング周波数を４２．５ｋＨｚとする。このとき、補間部１８に入力される受信デジタルデータ列のサンプリング周波数ｆ_ｉｓは、４２．５ｋＨｚとなる。この場合、
２×ｆ_ｃｈ＝２×１８．２ｋＨｚ＝３６．４ｋＨｚ
２×ｎ×ｆ_ｃｈ＝２×１２×１８．２ｋＨｚ＝４３６．８ｋＨｚ
２×ｆ_ｃｈ＜ｆ_ｉｓ＜２×ｎ×ｆ_ｃｈ
であるから、本実施形態が適用される。

　受信デジタルデータ列のサンプリング周波数を４３６．８ｋＨｚ以上に拡張するため、ｇ＝１２倍に補間するものとする。補間した後のデジタルデータ列のサンプリング周波数をｆ_ｓとすると、
ｆ_ｓ＝ｆ_ｉｓ×ｇ＝４２．５ｋＨｚ×１２＝５１０ｋＨｚ＞４３６．８ｋＨｚ
となる。

　補間する前のデジタルデータ列を、Ｕ（１）：ｕ_０，ｕ_１，ｕ_２，…とし、補間した後のデジタルデータ列を、Ｖ（１）：ｖ_０，ｖ_１，ｖ_２，…としたとき、補間部１８は、
ｖ_０～ｖ_１１＝ｕ_０
ｖ_１２～ｖ_２３＝ｕ_１
ｖ_２４～ｖ_３５＝ｕ_２
ｖ_３６～ｖ_４７＝ｕ_３
・・・
となるように補間を行う。

　補間後のデジタルデータ列Ｖ（１）の例を図で示すと、図１１のように階段状になる。図１１は、ｙ＝ｓｉｎ（２π１７０００ｔ）をサンプリング周波数４２．５ｋＨｚのＡ／Ｄ変換器２８への入力信号と仮定し、擬似的なデジタルデータ列作成した上で、補間を行った例である。このＶ（１）のサンプリング周波数ｆ_ｓは、５１０ｋＨｚなっている。従って、その後の処理は、図１０の構成の場合と同様に、Ｖ（１）を、第一デジタルＢＰＦ１１へ入力されるサンプリング周波数ｆ_ｓの受信デジタルデータ列とみなして、同じ手法で周波数解析を行えばよい。

　第一実施形態の数値例と同様、ｙ＝ｓｉｎ（２π１７８００ｔ）をサンプリング周波数４２．５ｋＨｚのＡ／Ｄ変換器２８への入力信号と仮定し、擬似的なデジタルデータ列作成した上で、上記の実施例の設定で実際に周波数解析を行なった結果は、第一実施形態の数値例と同じとなり、ｆ_ｔｇ、Ｔ_ｔｇを満たす。

　ここで、ｎ乗部１２の前のデジタルフィルタを、帯域通過フィルタであるデジタルバンドパスフィルタにした理由を、図１２を参照しながら説明をする。
　図１２Ａは、Ａ／Ｄ変換器２８への入力信号で、
ｙ＝ｓｉｎ（２π１７０００ｔ）＋２
とする。つまり、アナログフィルタで取りきれなかった直流分２が残っている、振幅１、周波数１７ｋＨｚのｓｉｎ波形とする。この入力信号を、サンプリング周波数が１７ｋＨｚの２．５倍である４２．５ｋＨｚでＡ／Ｄ変換したものが、ピーク周波数検出装置１に入力されるデジタルデータ列に相当する。

　このデジタルデータ列を、補間部１８で１２倍に拡張する。図１２Ｂは、補間部１８から出力されるデジタルデータ列の例である。このデジタルデータ列のサンプリング周波数は、１２×４２．５ｋＨｚ＝５１０ｋＨｚ相当になっている。

　図１２Ｃは、次数が８次でカットオフ周波数をｆ_ｃｌ＝１５.８ｋＨｚ、ｆ_ｃｈ＝１８.２ｋＨｚに設定したバタワース型のＩＩＲフィルタを用いた場合に、第一デジタルＢＰＦ部１１から出力されるデジタルデータ列である。第一デジタルＢＰＦ部１１の出力は、補間部１８から出力される余分な直流分、余分な低周波分、余分な高周波分が削除され、収束した波形はｓｉｎ波形に近い波形となっている。

　このように第一デジタルＢＰＦ部１１で周波数帯域を制限しておけば、次のｎ乗部１２でｎ乗した時に、［数４］式、［数５］式以外の周波数成分の発生が抑制される。すなわち、ｎ乗部１２の前段にデジタルバンドパスフィルタを備えることにより、デジタルデータ列に余分な周波数成分が含まれていても対応できることになる。したがって、前段のアナログフィルタ２７の性能を次段のＡ／Ｄ変換器２８でエイリアシングが発生しない程度に落とすことができ、Ａ／Ｄ変換器２８のサンプリング周波数も落とすことができるので、回路の規模とコストを削減できる。

　なお、カットオフ周波数ｆ_ｃｌ、ｆ_ｃｈのデジタルバンドパスフィルタを、カットオフ周波数ｆ_ｃｌのデジタルハイパスフィルタと、カットオフ周波数ｆ_ｃｈのデジタルローパスフィルタを組み合わせた構成としてもよい。

　５．第五実施形態
　図１３は、本発明の第五実施形態としてのピーク周波数検出装置６の構成を示すブロック図である。ピーク周波数検出装置６は、第一実施形態のｎ乗部１２と第二デジタルＢＰＦ１３を多重べき乗部１９に置き換えたものである。

　多重べき乗部１９は、図１４に示すように、デジタルデータ列をｍ_ｊ乗（ｍ_ｊは２以上の整数）するべき乗部（ｊ）（ｊ＝１，２，３，・・・，ｋ）と、べき乗部（ｊ）の出力から特定の周波数帯域ｆ_ｃｌ（ｊ）～ｆ_ｃｈ（ｊ）の信号を抽出するデジタルバンドパスフィルタ（ｊ）（ＢＰＦ（ｊ））とから構成されているべき乗ブロック（ｊ）がｋ段（ｋは２以上の整数）番号順に従属接続されたものである。

　べき乗ブロック（ｊ）のべき乗部（ｊ）の乗数ｍ_ｊは、置き換え前のｎ乗部の乗数ｎについて、
ｎ＝ｍ_１×ｍ_２×・・・×ｍ_ｋ
が成立するように選択される。
また、
ｆ_ｃｌ（ｊ）≒（ｍ_１×ｍ_２×…×ｍ_ｊ）×ｆ_ｃｌ
ｆ_ｃｈ（ｊ）≒（ｍ_１×ｍ_２×…×ｍ_ｊ）×ｆ_ｃｈ
と設定する。

　乗数ｍ_ｊが偶数の場合、［数５］式より、ｙをｍ_ｊ乗すると直流分が発生する。したがって、デジタルＢＰＦ（ｊ）は、直流分を削除できる必要がある。

　［数４］式、［数５］式で明らかなように、ｓｉｎ（２πｆｔ）をｍ乗する場合、乗数ｍを小さくすれば、発生する（ｍ－２）×ｆ以下の周波数成分の数を減らすことができる。そのため、その後のデジタルＢＰＦで、必要な帯域以外の周波数成分を、より抑制しやすくなる。したがって、本実施形態を適用すると、より高い乗数ｎまで対応できるようになる。目安としては、ｎ乗部の乗数ｎが１６を越える場合は、本実施形態の採用を検討するのがよい。

　本実施形態においても、第一実施形態同様、所望の周波数分解能ｆ_ｔｇと所望の時間窓長Ｔ_ｔｇを満たす。また本実施形態でも、近似計算を行わないため、計算するピーク周波数の精度劣化はない。

　本実施形態では、べき乗部と第二デジタルＢＰＦ部を多段構成にしたため、第一実施形態と比べると計算量が増える。しかし、より高い乗数ｎまで、所望の周波数分解能ｆ_ｔｇと所望の時間窓長Ｔ_ｔｇを満たして精度劣化無くピーク周波数を計算できるという長所は、これらのマイナス面を補って余りあるものがある。

　以下、図６に示したドップラー計測器２に本実施例のピーク周波数検出装置６を適用した場合について説明する。なお、下記の数値データ以外の必要な数値データは、第一実施形態の数値例と同じとする。後の説明を簡略化するため、第一実施形態同様、ｎ＝１２、ｆ_ｓ＝５１０ｋＨｚ、Ｎ＝４０９６、ｎ＝ｍ_１×ｍ_２とする。またｍ_１＝４、ｍ_２＝３とする。

　ｙ＝ｓｉｎ（２π１７８００ｔ）をサンプリング周波数５１０ｋＨｚのＡ／Ｄ変換器２８への入力信号と仮定して擬似的なデジタルデータ列作成した上で、第一デジタルＢＰＦ１１に入力し、第一デジタルＢＰＦ１１の出力を、Ａ（１）：ａ_０，ａ_１，ａ_２，…とする。

　べき乗部（１）で、Ａ（１）の各要素を４乗し、べき乗部（１）通過後のデジタルデータ列を、Ｂ（４）：ｂ_０，ｂ_１，ｂ_２，…
とする。すなわち、ｂ_ｉ＝（ａ_ｉ）^４（ｉ＝０，１，２，３，４，…）となる。

　このＢ（４）をデジタルＢＰＦ（１）に通し、デジタルＢＰＦ（１）デジタル通過後のデジタルデータ列をＣ（４）：ｃ_０，ｃ_１，ｃ_２，…とする。デジタルＢＰＦ（１）は、次数８次、カットオフ周波数６３．２ｋＨｚ（４×ｆ_ｃｌ）、７２．８ｋＨｚ（４×ｆ_ｃｈ）のバタワース型のＩＩＲフィルタとする。

　次に、べき乗部（２）で、Ｃ（４）の各要素を３乗し、べき乗部（２）通過後のデジタルデータ列をＤ（１２）：ｄ_０，ｄ_１，ｄ_２，…とする。すなわち、ｄ_ｉ＝（ｃ_i）^３（ｉ＝０，１，２，３，４，…）となる。

　このＤ（１２）をデジタルＢＰＦ（２）に通し、デジタルＢＰＦ（２）通過後のデジタルデータ列をＥ（１２）：ｅ_０，ｅ_１，ｅ_２，…とする。デジタルＢＰＦ（２）は、次数８次、カットオフ周波数１８９.６ｋＨｚ（４×３×ｆ_ｃｌ）、２１８.４ｋＨｚ（４×３×ｆ_ｃｈ）のバタワース型のＩＩＲフィルタとする。

　このＥ（１２）について、サンプル周波数ｆ_ｓ＝５１０ｋＨｚ、サンプル数Ｎ＝４０９６のＦＦＴで周波数を求めると、（１２×ｆ'）≒２１３６６２．１Ｈｚが得られる。求めるｆ'は、ｆ'＝（１２×ｆ'）／１２≒１７８０５．２Ｈｚとなる。ｆ'は、ｆ_ｔｇを満たして計算されたものである。実際に、ここでの理論値ｆ＝１７８００Ｈｚに対する誤差εは、ε＝ｆ'－f＝１７８０５．２－１７８００＝５．２Ｈｚであり、±（ｆ_ｔｇ）／２＝±１２／２＝±６Ｈｚに収まっている。

　また、時間窓長も第一実施形態同様、Ｎ／ｆ_ｓ＝４０９６／５１０ｋＨｚ＝８．０ｍｓ＜Ｔ_ｔｇ＝１０ｍｓであり、所望の時間窓長Ｔ_ｔｇを満たしている。

　６．第六実施形態
　ここまで説明したピーク周波数検出装置１、３－６に、パラメータ設定部を追加してもよい。図１５は、第一実施形態のピーク周波数検出装置１にパラメータ設定部２０を追加した構成を示すブロック図である。パラメータ設定部２０は、プロセッサ、メモリ、入出力機構を備えるコンピュータであって、図示しないキーボードやマウスやタッチパネルディスプレイ等の図示しない操作部を用いたユーザの入力に応じてピーク周波数検出装置１、３－６にパラメータの値を設定する。パラメータは、これまで説明したとおり、次のようなものである。
受信デジタルデータ列のサンプリング周波数ｆ_ｉｓ
デジタルデータ列のサンプリング周波数ｆ_ｓ
所望の周波数分解能ｆ_ｔｇ
所望の時間窓長Ｔ_ｔｇ
ｎ乗部の乗数ｎ（ｎは２以上の整数）
第一デジタルＢＰＦの周波数帯域　ほぼｆ_ｃｌ～ｆ_ｃｈ（ｆ_ｃｌ＜ｆ_ｃｈ）
第二デジタルＢＰＦの周波数帯域　ほぼｎ×ｆ_ｃｌ～ｎ×ｆ_ｃｈ
ＦＦＴのサンプリング周波数ｆ_ｓ
ＦＦＴのサンプリング数Ｎ

　パラメータ設定部２０は、上記全パラメータの数値をユーザが選択可能な入力値と対応付けて予めメモリに記憶しておき、ユーザの入力値に応じて設定してもよいし、一部のパラメータの数値をユーザが選択可能な入力値と対応付けて予めメモリに記憶しておき、ユーザの入力に応じて設定し、残りのパラメータは、入力された数値に応じてプロセッサが計算で求めても良い。

　本実施形態では、周波数に対応した音階番号のみをユーザに入力させ、入力された音階番号に応じて全パラメータの数値をメモリから取得して設定する調律補助装置として使用されるピーク周波数検出装置について説明する。すなわち、パラメータ設定部２０は、音階番号１～８８のいずれか一つの入力に対応するパラメータの数値をメモリから取得して設定する。図１６は、音階番号を入力する状態におけるタッチパネルディスプレイの画面構成図の一例である。全パラメータの数値は、音階番号Ｐに応じて、前もって全て決められて不揮発性メモリに記憶されている。音階番号は平均律ピアノの鍵に対応し、音階番号１は最低音の周波数、音階番号８８は最高音の周波数に対応する。すなわち、音階番号Ｐに対応する周波数ｆ_ｐは次式［数２１］によって定まる値が予め不揮発性メモリに記憶されている。

　そして例えば、Ｐ＝４９（ｆ_ｐ＝４４２Ｈｚ）に対する各パラメータの数値は、以下のように定められている。

　なおｆ_ｔｇは、ｆ_ｐの２セント分にあたる周波数である。セント値は、周知の通り、２音間の周波数比を対数表現で表した値であって、１００セントが平均律１２音階の半音に相当する。Ｔ_ｔｇは、時間窓長（１／ｆ_ｔｇ）の１／１０である。ｆ_ｃｌは、ｆ_ｐの５０セント下の周波数である。ｆ_ｃｈは、ｆ_ｐの５０セント上の周波数である。Ｎは、ｆ_ｓ／（ｎ×ｆ_ｔｇ）＝２４０００／（１６×０．５１０９）≒２９３６～ｆ_ｓ×Ｔ_ｔｇ＝２４０００×０．１９５７≒４６９７間の２のべき乗の整数でなければならない。

　以上のように前もって取り決めておけば、パラメータ設定部２０は、音階番号Ｐの入力に応じて全パラメータの数値を設定することができる。このように設定されたパラメータを用いることにより、受信デジタルデータ列のピーク周波数が計算できる。なお、ピーク周波数検出装置１に入力するデジタルデータ列は、図示しないマイク、Ａ／Ｄ変換器を通じてピーク周波数検出装置１にリアルタイムで順次入力されても良いし、メモリに格納されていても良い。

　Ｐ＝４９（ｆ_ｐ＝４４２Ｈｚ）とし、ｙ＝ｓｉｎ（２π４４２ｔ）をサンプリング周波数２４ｋＨｚでサンプリングした擬似的なデジタルデータ列を作成した上で、実際に周波数解析を行なった結果は、１６×ｆ'≒７０７２．２６５６Ｈｚが得られる。求めるｆ'は、ｆ'＝（１６×ｆ'）／１６≒４４２．０１６６Ｈｚとなる。求めたピーク周波数の周波数分解能は、所望の周波数分解能ｆ_ｔｇを満たしている。実際に、ここでの理論値f＝４４２Ｈｚに対する誤差εは、ε＝ｆ'－ f＝４４２．０１６６－４４２＝０．０１６６Ｈｚであり、±（ｆ_ｔｇ／２）＝±０．５１０９／２≒±０．２５５Ｈｚに収まっている。また、時間窓長は、Ｎ／ｆ_ｓ＝４０９６／２４ｋＨｚ＝０．１７１ｓであり、所望の時間窓長Ｔ_ｔｇ＝０．１９５ｓを満たしている。なお、求めたピーク周波数の出力形式は、単位Ｈｚでも良いし、ｆp＝４４２Ｈｚに対するセント値（＋０．０６５セント）でも良い。

　７．効果
　以上説明した本発明の実施形態によると、ＦＦＴによる周波数解析の際に制約となる周波数分解能（ｆ_０）と時間窓長（Ｔ_０）との間にあるｆ_０＝１／Ｔ_０と言う相反関係の問題を回避し、前記所望の周波数分解能ｆ_ｔｇと前記所望の時間窓長Ｔ_ｔｇでもって、信号波のピーク周波数の検出を可能にすることができる。受信デジタルデータ列のサンプリング周波数は、２×ｆ_ｃｈ以上であればピーク周波数計算が可能となる。そして、近似計算やカーブフィッティング、平均化などの処理を必要としないため、ピーク周波数の計算精度の劣化は無い。

　また、入力信号に余分な周波数帯域の成分（直流分、低周波分、高周波分）が含まれていても、Ｓ／Ｎ比が悪くても、殆ど振幅が無いような場合でも、第一デジタルＢＰＦ、第二デジタルＢＰＦにより、ピーク周波数を計算するのに必要な帯域の周波数成分を抽出するため、問題なくピーク周波数を計算することができる。したがって、ピーク周波数検出装置の前段のハードウェアの要求仕様を下げることができ、小型化やコストダウンをはかることができる。具体的には例えば、第一デジタルＢＰＦのカットオフ周波数ｆ_ｃｈに対し、入力信号のサンプリング周波数が２×ｆ_ｃｈ以上であればピーク周波数が計算できるため、前段のアナログフィルタやＡ／Ｄ変換器などのハードウェアの要求仕様を下げることができる。また、べき乗部に設定する乗数ｎを大きめに設定すれば、デジタルバンドパスフィルタの次数を上げて急峻にし、外来ノイズに強くすることができる。

　８．他の実施形態
　尚、本発明の技術的範囲は、上述した実施例に限定されるものではなく、本発明の要旨を逸脱しない範囲内において種々変更を加え得ることは勿論である。

　例えば、ＩＩＲ、バタワース型、８次のデジタルＢＰＦを例示したが、他の形式のものも好適に使用できる。例えば、ＦＩＲ（Finite impulse response）、チェビシェフ型などが好適に使用できる例として挙げられる。次数も８次に限らない。また、カットオフ周波数も、状況に応じて、抽出範囲を広くあるいは狭くすれば良い。また入力信号が、０を中心に振れる、高調波成分の殆ど無い信号（つまり、単一周波数のｓｉｎカーブに近い信号）であるなら、第一デジタルＢＰＦをローパスフィルタにしても良い。

　また反射エコーにおけるドップラー周波数の検出や調律補助装置への適用を例に挙げて本発明を説明したが、本発明の適用範囲はこれに限られるものではない。ＦＦＴにより信号波のピーク周波数の検出を行う用途に、本発明は広く一般的に適用可能である。

　また、上記実施の形態に係る各機能部は、１または複数のＬＳＩ（ＬａｒｇｅＳｃａｌｅＩｎｔｅｇｒａｔｉｏｎ）により実現されてもよく、また、複数の機能部が１のＬＳＩにより実現されてもよい。また、集積化の手法としてはＬＳＩに限定されるものではなく、専用回路または汎用プロセッサで実現することとしてもよい。さらには、ＬＳＩ製造後にプログラムすることが可能なＦＰＧＡ（Ｆｉｅｌｄ　Ｐｒｏｇｒａｍｍａｂｌｅ　Ｇａｔｅ　Ａｒｒａｙ）や、ＬＳＩ内部の回路セルの接続や設定を再構成可能なリコンフィギュラブル・プロセッサ（Ｒｅｃｏｎｆｉｇｕｒａｂｌｅ　Ｐｒｏｃｅｓｓｏｒ）を利用してもよい。また、半導体技術の進歩又は派生する別技術により、ＬＳＩに置き換わる集積回路化の技術が登場すれば、当然、その技術を用いて機能ブロックの集積化を行ってもよい。

１…ピーク周波数検出装置、２…ドップラー計測器、３…ピーク周波数検出装置、４…ピーク周波数検出装置、５…ピーク周波数検出装置、６…ピーク周波数検出装置、１１…第一デジタルＢＰＦ部、１２…ｎ乗部、１３…第二デジタルＢＰＦ部、１４…ＦＦＴ部、１５…１／ｎ倍部、１６…間引き部、１７…間引き部、１８…補間部、１９…多重べき乗部、２０…パラメータ設定部、２１…送受信器、２２…送受信切替回路、２３…送信回路、２４…受信増幅器、２５…変調器、２６…局部発振回路、２７…アナログフィルタ、２８…Ａ／Ｄ変換器

Claims

　所定の周波数帯域（ｆ_ｃｌ～ｆ_ｃｈ）においてパワースペクトルが最大となるピーク周波数を検出するピーク周波数検出装置であって、
　デジタルデータ列の各要素をｎ乗（ｎは２以上の整数）するｎ乗部と、
　サンプリング周波数ｆ_ｓと周波数分解能ｆ_ｔｇと時間窓長Ｔ_ｔｇに応じて決まるＮ（Ｎは２のべき乗の整数）個のｎ乗されたサンプリング周波数ｆ_ｓのデジタルデータ列に対して高速フーリエ変換を行って得られるパワースペクトルの最大値に対応する周波数を仮想ピーク周波数として導出するＦＦＴ部と、
　前記仮想ピーク周波数を１／ｎ倍した値をデジタルデータ列のピーク周波数として出力する１／ｎ倍部と、
　を備え、
ｎ≧１／（ｆ_ｔｇ×Ｔ_ｔｇ）
ｆ_ｓ／（ｎ×ｆ_ｔｇ）≦Ｎ≦ｆ_ｓ×Ｔ_ｔｇ
ｆ_ｓ＞２×ｎ×ｆ_ｃｈ
　を満たすピーク周波数検出装置。
　前記所定の周波数帯域に含まれる周波数のデジタルデータ列を抽出する第一デジタルバンドパスフィルタを備え、
　前記ｎ乗部には、前記第一デジタルバンドパスフィルタの出力が入力される、
　請求項１に記載のピーク周波数検出装置。
　サンプリング周波数ｆ_ｉｓのデジタルデータ列を１／ｒ（ｒは２以上の整数）に間引いてサンプリング周波数をｆ_ｓにする間引き部を備え、
　前記第一デジタルバンドパスフィルタには前記間引き部の出力が入力される、
ことを特徴とする請求項２に記載のピーク周波数検出装置。
　前記デジタルデータ列をｇ倍（ｇは２以上の整数）に補間してサンプリング周波数をｆ_ｓにする補間部を備え、
　前記第一デジタルバンドパスフィルタには前記補間部の出力が入力される、
ことを特徴とする請求項２記載のピーク周波数検出装置。
　ｎ乗されたＮ個の前記デジタルデータ列から第二の周波数帯域に含まれるデジタルデータ列を抽出する第二デジタルバンドパスフィルタを備え、
　前記ＦＦＴ部には、前記第二デジタルバンドパスフィルタによって抽出されたデジタルデータ列が入力され、
　前記第二の周波数帯域は、ほぼｎ×ｆ_ｃｌ～ｎ×ｆ_ｃｈである、
　請求項１に記載のピーク周波数検出装置。
　前記第二デジタルバンドパスフィルタで抽出されたデジタルデータ列を１／ｒ（ｒは２以上の整数）に間引いてサンプリング周波数をｆ_ｓとする間引き部を備え、
　前記ＦＦＴ部には、前記間引き部の出力が入力される、
　請求項５に記載のピーク周波数検出装置。
　前記特定の周波数帯域に含まれる周波数のデジタルデータ列を抽出する第一デジタルバンドパスフィルタと、
　前記ｎ乗部の出力から第二の周波数帯域に含まれるデジタルデータ列を抽出する第二デジタルバンドパスフィルタとを備え、
　前記ｎ乗部には、前記第一デジタルバンドパスフィルタの出力が入力され、
　前記ＦＦＴ部には、前記第二デジタルバンドパスフィルタの出力が入力され、
　前記第二の周波数帯域は、ほぼｎ×ｆ_ｃｌ～ｎ×ｆ_ｃｈである、
　請求項１記載のピーク周波数検出装置。
　前記ｎ乗部に代えて、入力されるデジタルデータ列をｍ_ｊ乗（ｍ_ｊは２以上の整数）するべき乗部（ｊ）（ｊ＝１，２，…，ｋ）と前記べき乗部（ｊ）の出力から特定の周波数帯域ｆ_ｃｌ（ｊ）～ｆ_ｃｈ（ｊ）の信号を抽出するデジタルバンドパスフィルタ（ｊ）とを備えるべき乗ブロック（ｊ）をｋ（ｋは２以上の整数）段備える多重べき乗部を備え、
ｎ＝ｍ_１×ｍ_２×…×ｍ_ｋ
ｆ_ｃｌ（ｊ）≒（ｍ_１×ｍ_２×…×ｍ_ｊ）×ｆ_ｃｌ
ｆ_ｃｈ（ｊ）≒ｍ_１×ｍ_２×…×ｍ_ｊ）×ｆ_ｃｈ
である請求項１～４のいずれか一項に記載のピーク周波数検出装置。
　ユーザの指示を受け付ける操作部と、
　前記指示に応じたｎ、ｆ_ｓ、Ｎの少なくともいずれかを設定するパラメータ設定部と、
　を備える請求項１から８のいずれか一項に記載のピーク周波数検出装置。
　所定の周波数帯域（ｆ_ｃｌ～ｆ_ｃｈ）においてパワースペクトルが最大となるピーク周波数を検出するピーク周波数検出方法であって、
　デジタルデータ列の各要素をｎ乗（ｎは２以上の整数）し、
　サンプリング周波数ｆ_ｓと周波数分解能ｆ_ｔｇと時間窓長Ｔ_ｔｇに応じて決まるＮ（Ｎは２のべき乗の整数）個のｎ乗されたサンプリング周波数ｆ_ｓのデジタルデータ列に対して高速フーリエ変換を行って得られるパワースペクトルの最大値に対応する周波数を仮想ピーク周波数として導出し、
　前記仮想ピーク周波数を１／ｎ倍した値をデジタルデータ列のピーク周波数として出力する、
　ことを含み、
ｎ≧１／（ｆ_ｔｇ×Ｔ_ｔｇ）
ｆ_ｓ／（ｎ×ｆ_ｔｇ）≦Ｎ≦ｆ_ｓ×Ｔ_ｔｇ
ｆ_ｓ＞２×ｎ×ｆ_ｃｈ
　を満たすピーク周波数検出方法。
　所定の周波数帯域（ｆ_ｃｌ～ｆ_ｃｈ）においてパワースペクトルが最大となるピーク周波数を検出するピーク周波数検出プログラムであって、
　デジタルデータ列の各要素をｎ乗（ｎは２以上の整数）するｎ乗部と、
　サンプリング周波数ｆ_ｓと周波数分解能ｆ_ｔｇと時間窓長Ｔ_ｔｇに応じて決まるＮ（Ｎは２のべき乗の整数）個のｎ乗されたサンプリング周波数ｆ_ｓのデジタルデータ列に対して高速フーリエ変換を行って得られるパワースペクトルの最大値に対応する周波数を仮想ピーク周波数として導出するＦＦＴ部と、
　前記仮想ピーク周波数を１／ｎ倍した値をデジタルデータ列のピーク周波数として出力する１／ｎ倍部と、
　してコンピュータを機能させ、
ｎ≧１／（ｆ_ｔｇ×Ｔ_ｔｇ）
ｆ_ｓ／（ｎ×ｆ_ｔｇ）≦Ｎ≦ｆ_ｓ×Ｔ_ｔｇ
ｆ_ｓ＞２×ｎ×ｆ_ｃｈ
　を満たすピーク周波数検出プログラム。