フーリエ変換でのサンプル数を削減した、短時間信号のピークパワースペクトルを検出する方法及び装置
【課題】測定対象とする時間範囲の分解能を高めつつ周波数分解能を高め、しかもフーリエ変換におけるサンプリング周波数を低くしてサンプル数を少なくし、演算量を減らす。
【解決手段】受信信号を変調器5によって周波数変換し、アナログフィルタ7を介して周波数変換後の信号をA/D変換器61によりオーバーサンプリングでデジタル信号に変換し、デジタルBPF62により所望の周波数帯域幅Δfpに相当する成分を抽出し、抽出された信号をダウンサンプリング部63により最低周波数帯域にダウンサンプリングする。その後、信号切り出しゲート11により、ダウンサンプリング部63からのデジタルデータ列を所望の時間幅内で切り出し、周波数分解能を満たすようにゼロ付加部12によってゼロデータを付加し、ゼロデータが付加されたデジタルデータ列に対しFFT処理部13によって高速フーリエ変換を行う。
【解決手段】受信信号を変調器5によって周波数変換し、アナログフィルタ7を介して周波数変換後の信号をA/D変換器61によりオーバーサンプリングでデジタル信号に変換し、デジタルBPF62により所望の周波数帯域幅Δfpに相当する成分を抽出し、抽出された信号をダウンサンプリング部63により最低周波数帯域にダウンサンプリングする。その後、信号切り出しゲート11により、ダウンサンプリング部63からのデジタルデータ列を所望の時間幅内で切り出し、周波数分解能を満たすようにゼロ付加部12によってゼロデータを付加し、ゼロデータが付加されたデジタルデータ列に対しFFT処理部13によって高速フーリエ変換を行う。
【発明の詳細な説明】
【技術分野】
【0001】
本発明は、時間に応じて変化する信号の中から所望の時間幅内の部分を切り出し、その切り出された部分に対してフーリエ変換を適用することによってピークパワースペクトルを検出する方法及び装置に関し、特に、所望の周波数分解能力を維持しながらフーリエ変換におけるサンプル数を減らして演算処理時間を短縮でき、かつ検出されるスペクトルの周波数精度を向上させることができる方法及び装置に関する。
【背景技術】
【0002】
所定の周波数帯域幅を有して時々刻々と変化する信号の中から所望の時間幅の部分のみを切り出し、その切り出された部分内での信号のピークパワースペクトルを検出することは、後述するように、例えばドップラー計測などで広く行われている。ここでいう所望の時間幅は、信号全体の継続時間に比べた場合には十分に短く、例えば、その開始時刻が与えられているかあるいは開始時刻を自動検出できるものである。また、検出されるピークパワースペクトルは、切り出された部分内での信号の周波数スペクトルの中で一番強いパワースペクトルを有するものである。パワースペクトルの検出には、一般的にフーリエ変換が用いられ、計算量の観点から、中でも、FFT(高速フーリエ変換:Fast Fourier Transform)が使用される。
【0003】
ところで、ある信号を媒質中に発射(送信)し、その媒質中から反射した信号(エコー信号)を受信したとき、送信した時刻をスタート時刻とすると、スタート時刻から反射によるエコー信号を受信するまでの時間は、媒質中でのその信号の伝搬距離すなわち対象物までの距離に依存する。したがって、信号の伝搬速度(音波であれば音速)を考慮して、スタート時刻から所望の時間経過した時点で受信されるエコー信号は、その時間経過に応じた距離の位置にある対象物からのエコー信号となる。スタート時刻からの時間(すなわち距離)に応じてエコー信号の強度を映像として表示するようにした応用機器が、例えば、レーダーや魚群探知機である。
【0004】
信号を反射する対象物が相対的に移動しているとすると、反射エコーでは、ドップラー効果により、送信した周波数からの周波数変位が観測される。この周波数差、すなわち媒質に対して送信したときの周波数と反射エコーで観測される周波数との差をドップラー周波数という。ドップラー計測では、一定の周波数の信号を媒質中に発射し、媒質中の対象物によって反射された信号(反射エコー)を受信し、その周波数解析を行って、ドップラー周波数を決定するものである。ドップラー計測を行うものとして、例えば、ドップラーレーダー、ドップラーソナー、ADCP(音波ドップラー流速プロファイラー;acoustic Doppler current profiler)、潮流計、音波レーダー、超音波ドップラー速度計などがあり、これらは、ドップラー周波数に基づいて、対象物との相対速度を決定する。
【0005】
送信回路によって一定の送信周波数ftxの信号を生成して、この信号を送受波器から媒質に発射する。媒質中を対象物(反射体)が移動しているとして、ドップラー周波数fdopは、対象物の移動速度をV(m/s)、媒質中での信号の伝搬速度をC(m/s)とすると、媒質中の信号の伝搬速度が対象物の速度よりも十分に速い場合(C≫V)には、次の近似式で表される。
【0006】
【数1】
【0007】
ここでθは、対象物から送受波器を見たときの方向と対象物の移動方向とがなす角度であり、同じ方向であれば(対象物が送受波器にむかってまっすぐ進んでいるときには)、θ=0となる。したがって、ドップラー周波数fdopは、対象物が近付いているときには正(+)となり、遠ざかるときには負(−)となる。
【0008】
エコー信号の周波数解析を行ってドップラー周波数fdopを求めるが、この周波数解析には、一般的に、フーリエ変換の手法が使われている。フーリエ変換として、パワースペクトル解析を行う場合にはDFT(離散フーリエ変換;Discrete Fourier Transform)が用いられ、あるいは信号の標本数(サンプル数)が2のべき乗個であるときには、高速演算処理が可能なFFT(高速フーリエ変換:Fast Fourier Transform)が通常用いられている。DFTもFFTも複素乗算演算を反復するものであるが、サンプル数をNとすると、DFTでの複素乗算数回数はN×N回(すなわちN2回)であり、N=2nとすると、FFTでの複素乗算回数は、2n×log2 2n回(すなわちN×log2 N回)であり、Nが大きくなるほど、FFTはDFTに比べて高速演算ができる。
【0009】
ドップラー周波数を含んだ受信信号からFFTあるいはDFTなどの手法によりフーリエ変換で周波数を求めるためには、アナログ信号である受信信号をデジタル信号に変換することが必要である。サンプリング定理から、このアナログ/デジタル(A/D)変換でのサンプリング周波数は、入力信号の2倍以上の周波数でなければならない。言い換えると、A/D変換の対象となる入力信号はサンプリング周波数の1/2以下でなければならない。
【0010】
ところで、ドップラー計測では、対象物の位置と速度の両方に関する情報を取得したい場合が多い。例えば、海洋において深さごとの潮流の方向と速度とを求めて流速プロファイルを得るためには、海中に対して斜め方向に音波パルスを放射して反射エコーを取得し、周波数解析を行う。この際、音波パルスの放射から所定の時間の経過後にある時間幅で反射エコーを切り出し、切り出された部分に対して周波数解析を行うことになる。音波パルスの発射からの所定の時間によって、どの深度の流速を計測しようかが定まることになり、反射エコーを切り出す時の時間幅によって、どの程度の深度範囲(厚さ)にわたる平均として流速が算出されるかが定まることになる。例えば、切り出す時の時間幅が大きければ、海中における厚みの大きい層での平均的な流速が求められることになる。これに対し、時間幅が小さければ、薄い層での流速を求めることが可能になって、結果として、細かな流速プロファイルが得られることになる。
【0011】
しかしながら、ドップラー周波数を求めるための周波数解析としてのフーリエ変換の問題点として、周波数分解能と時間(距離)分解能とが相反する関係にある、ということが挙げられる。ここでいう周波数分解能は、ドップラー周波数を求める際の分解能であり、これは速度の分解能に比例する。時間分解能とは、反射エコーを切り出す時にどこまで時間幅を小さくできるかということであり、これは、ドップラー計測における距離あるいは位置の分解能に比例する。時間分解能が粗ければ、上述した流速プロファイルの場合であれば、深度方向の厚みが大きな層でのおおざっぱな平均的な流速しか得られないことになる。
【0012】
離散フーリエ変換(DFT)は、次の式で表される。
【0013】
【数2】
【0014】
DFTまたはFFTでのサンプリング周波数をfft、フーリエ変換の基本波の周期をT0とすると、周期の逆数は、フーリエ変換での基本周波数f0と呼ばれ、
f0=1/T0 :式1−4
と表される。この基本周波数f0は、DFTやFFTでの最小分解周波数すなわち周波数の分解能となる。サンプルされた信号点の総数すなわちサンプル数をNとすると、
T0=N/fft :式1−5
が当然に成立し、周波数分解能である最小分解周波数f0は、
f0=1/T0=fft/N :式1−6
と表されることになる。逆に、所望の周波数分解能f0が与えられたとして、必要なサンプル数Nは、
N=fft/f0 :式1−7
と表されることになる。フーリエ変換では、基本波の周期T0と等しい時間幅の信号から等間隔でN個のサンプルを取り出すから、T0が信号の切り出し時間(切り出しの時間幅)ということになる。
【0015】
FFTで得られるパワースペクトルの例が図1に示されている。FFTのパワースペクトルは、サンプル数Nに対し、0から(N/2)−1までの一定の周波数間隔で周波数ごとの値14,15が与えられるスペクトルとして表される。ここで周波数間隔は上記のf0と一致し、f0=1/T0である。
【0016】
式1−4に示すように、周波数分解能f0は、信号の切り出し時間幅T0の逆数であるから、周波数分解能を上げる(f0を低くする)ためには、信号の切り出し時間幅T0を長くしなければならない。その結果、周波数分解能を2倍細かくすると、時間分解能すなわち位置(距離)の分解能は2倍粗くなって劣化する。
【0017】
対象物との相対的な速度によって、ドップラー周波数fdopは変化する。このドップラー周波数ffopが変化する範囲をΔfpとおくことにすると、近づく方向と遠ざかる方向での移動が同等に起こるものと仮定すれば(この過程は、潮流計などの用途において妥当な仮定である)、ドップラー信号の周波数範囲ΔfPの中心周波数fc、つまりドップラー周波数fdopがゼロのときの受信信号の周波数は、送信周波数ftxと一致する(fc=ftx)。
【0018】
以下の条件に従い、超音波のパルス状の信号を海中に送信し、移動物からの反射エコーのドップラー周波数を測定する場合を考える。
【0019】
海水中の音波伝搬速度C:C=1500(m/s)
送信周波数ftx:ftx=120kHz
検出最大速度(水平方向)V:V=15(m/s)
また、音波は水平方向から斜め方向(θ=60°)に送受信するものとする。式1−1により、
【0020】
【数3】
【0021】
となり、検出最大速度が15m/sであるとすれば、ドップラー信号として観測される周波数範囲Δfpは、120±1.2(kHz)となる。
【0022】
測定速度範囲すなわち検出最大速度は、通常、予め決められているので、ドップラー信号についても、その周波数範囲が予め定められて測定する。
【0023】
ドップラー計測による速度精度を0.1m/sとすると、それに必要な周波数分解能f0は12Hzとなる。ドップラー計測において検出対象とする位置の範囲(位置分解能)を7.5mとすると、信号の切り出し時間幅は、海水中で音波が1mを往復する時間である2/1500(m/s)を乗ずると、7.5×2/1500=0.01となって、10msとなる。
【0024】
しかしながら、周波数分解能f0である12Hzを式1−4に代入すると、T0=1/f0=1/12≒83.3msとなり、83.3msは、海水中で音波が1mを往復する時間2/1500(m/s)で割ることにより距離に換算すると、83.3(ms)×1500(m/s)/2≒62.5(m)となる。
【0025】
すなわち、速度精度を0.1m/sとするために周波数分解能f0を12Hzとすると、信号の切り出し時間幅として83msが必要であり、これは距離に換算すれば63mに相当し、位置(距離)に関する分解能を7.5mとする要求を満足することできない。逆に、位置に関する分解能が7.5mである信号でFFTを行うとすると、式1−4より周波数分解能f0が100Hzとなり、速度精度から定められた周波数分解能12Hzの要求を満足できない。
【0026】
このようにフーリエ変換には、周波数分解能と、測定対象に関する時間分解能とが両立しない、という基本的な問題があり、これは、ドップラー計測の場合には、周波数分解能あるいは速度分解能と、距離あるいは位置の分解能とが両立しないことを意味する。
【0027】
図2は、一般的なドップラー計測機器の構成の一例として、受信信号の中から所望の時間幅の部分のみを切り出し、その切り出された部分内での信号のピークパワースペクトルを検出する従来のピークパワースペクトル検出装置の構成を示している。ここでは、信号が音響信号であり、信号が伝搬する媒質が水(あるいは海水)であるとする。音響信号の送受波器1は、送受信切替え回路2を介し、送信回路3の出力と受信増幅器4の入力に接続している。送信回路3は、送信周波数ftxの信号を生成するものである。受信増幅器4の出力には、受信信号を中間周波数信号に周波数変換する変調器5が設けられており、変調器5には、局部発振回路6から局部発振周波数flocの信号が供給されている。変調器5の出力すなわち中間周波数信号は、高次アナログフィルタ60を介してアナログ/デジタル変換器(A/D)10に供給されてデジタル信号に変換され、このデジタル信号は、所定のサンプル数(すなわち時間幅)で信号を切り出す信号切り出しゲート11に供給される。信号切り出しゲート11で切り出された信号は、次にFFT処理部13に送られてFFTによる周波数解析の対象となる。FFT処理部13から、計算されたパワースペクトルが出力される。
【0028】
この装置では、水中に音波パルスを発射する時には、送信周波数ftxの信号を送信回路3から送受信切替え回路2を介して送受波器1に供給する、その後、直ちに送受信切替え回路2を受信側とし、媒質(水中)からの反射エコーを送受波器1で受信して電気信号である受信信号に変換し、送受信切替え回路2を介して受信信号を受信増幅器4に供給する。その結果、受信信号は受信増幅器4において信号増幅され、変調器5に対して供給されることになる。
【0029】
図3は、受信信号や中間周波数信号などの周波数間の関係を示している。受信増幅器4から出力される信号自体は、図3において符号52で示されるような周波数帯域を有するものとする。受信信号の帯域は、符号52で示す周波数帯域の一部であり、検出対象となるドップラー信号は、符号51で示される周波数範囲Δfpの中にある。周波数範囲Δfpの中心、すなわち、ドップラー周波数がゼロの時の周波数fc(47)は、送信周波数ftxと同じである(fc=ftx)。局部発振回路6は、周波数fcよりも中間周波数fmidだけ高い局部発振周波数floc(45)の局部発振信号を生成し(すなわちfloc>fc)、変調器5に供給する。その結果、変調器5では、局部発振信号を受信信号に乗算することにより、受信信号を中間周波数fmidの中間周波数信号に変換する。
【0030】
受信信号の周波数帯域のうち、トップラー計測のために着目する周波数範囲Δfpの上限を上限周波数fmax(49)とし、下限を下限周波数fmin(48)とすると、フーリエ変換すべき周波数範囲Δfp(51)は、Δfp=fmax−fminで表される。
【0031】
変調器5における中間周波数fmidの信号の生成は、下記の式で表されるものであり、2つの変換周波数から、高次アナログフィルタ60により、低い方の周波数を得る。
【0032】
【数4】
【0033】
ドップラー信号の周波数範囲Δfpは、中間周波数fmidを中心に、
【0034】
【数5】
【0035】
の範囲に変換されることとなり、検出対象の周波数範囲Δfpは、周波数変換しても保たれる。上記の式は、floc>fcの時の処理を示しているが、floc<fcの時も、複号(符号(±))が逆になるだけで、同様に処理がなされる。
【0036】
図3において、右側部分が受信信号と局部発振周波数flocの信号との関係を表しており、ここでfmid(46)は、fmid=floc−fcの関係を満たしている。一方、図3の左側部分では、中間周波数信号に変換後のスペクトルを示しており、ここでは中間周波数信号の中心周波数はfmid(50)となっている。また、上述の周波数帯域52も周波数帯域54に変換され、検出対象の周波数範囲Δfp(51)は、中心周波数がfmidである周波数範囲Δfp(53)に変換されている。すなわち帯域53は、変調器5によって周波数変換して低い方の周波数成分だけを取り出すことによって中間周波数信号とされ、高次アナログフィルタ60を通過させた後の受信信号におけるΔfpの帯域を表している。
【0037】
高次アナログフィルタ60は、変調器5から出力される2つの変換周波数のうちの低い方を選択し、さらに、A/D変換でのエリアシングを防ぐためのものであるので、このフィルタには、バンドパスフィルタ(BPF;band-pass filter)かローパスフィルタ(LPF;low-pass filter)が使用される。以下では、説明のため、通過帯域におけるリップルが少ないバターワース型LPFを高次アナログフィルタ60に使用するものとする。A/D変換でのエリアシングを考えると、A/D変換器10におけるA/D変換のデータ幅が8ビットであるときには、A/D変換器10のサンプリング周波数f2adの1/2の周波数すなわちf2ad/2での信号レベルが、フィルタにおける通過帯域での信号レベルを基準として−48dB以下であれば、完全にエリアシングは起こらない。同様に、16ビット幅のA/D変換の場合には、信号レベルが−96dB以下であればエリアシングは完全に起らない。したがって、高次アナログフィルタ60は、中間周波数信号における着目する周波数帯域(周波数fmidを中心とする範囲Δfp)をその通過帯域内に含むとともに、周波数f2ad/2では上述した減衰量(8ビット幅のA/D変換であれば−48dB、16ビット幅であれば−96dB)を有するものである必要がある。逆に言えば、減衰量がこの値となるまでfft/2を高くし、その分、サンプリング周波数f2adも高くする必要がある。
【0038】
以下では、A/D変換のデータ幅が8ビットであるものとして説明を進める。図4は、LPFの周波数特性の一例を示している。以下、図4を用いて、受信信号の周波数帯域を説明する。
【0039】
受信信号の周波数fcが120kHz、ドップラー周波数の最大値が1200Hz、フーリエ変換に要求される周波数分解能が12Hzであるとする。式2−3に示すように、受信信号を中間周波数信号に変換したとしても、ドップラー周波数の周波数範囲Δfpそのものは変化しない。その一方で、受信信号と局部発振周波数の信号がクロストークしたり混変調を起したりしないように、局部発振周波数flocを選ぶ必要がある。ここでは、受信信号からその周波数で1割以上離れた局部発振周波数flocを選ぶこととし、17KHz離すこととする。局部発振周波数flocはfloc=fc+17kHz=137KHzとなり、中間周波数fmidは17KHzとなる。このとき、アナログLPFの通過帯域(全帯域)Wは、
【0040】
【数6】
【0041】
とする必要がある。Δfp/2=1200Hzであり、これに適切なα(減衰域)とfmidとを加算すると、LPFは、図4に示すように、20KHzまでの信号を完全に通過させるものであればよい。図4において、符号74は2次のLPFの特性を示しており、符号73は4次のLPFの特性を示しており、符号72は8次のLPFの特性を示している。アナログLPFでは、構成部品の精度や温度係数などから、8次のLPFを構成することはコスト面からも実用的ではないので、4次のLPFを高次アナログフィルタ26に用いるものとする。
【0042】
上述したように、8ビット幅A/D変換では、サンプリング周波数の1/2の周波数での信号レベルを通過帯域に比べて−48dBとすれば、エリアシングの問題は生じない。図5は、中間周波数に変換された信号とA/D変換のサンプリング周波数の関係を示している。図5においては、サンプリング周波数は、f2ad(78)で示されている。ここでは4次のLPFを使うこととしたので、−48dBの周波数は80kHzとなり、サンプリング周波数f2adの半分の周波数f2ad/2(77)も80kHzに設定される。したがって、A/D変換は、160kHzのサンプリング周波数f2adで行なわれることになる。
【0043】
その結果、式1−6から、分解能の条件を満たすサンプル数Nは、N=f2ad/f0=160000/16=10000となる。FFTのデータ点数は2のべき乗個とすべきであることから、FFTではサンプル数を16384(=214)とし、このときの周波数分解能f0は、f0=160000/16384=9.77[Hz]となる。このとき、FFTで実行される複素乗算の回数は、
FFT複素乗算回数: 214×log2 214=229376
となる。FFTの代わりにDFTを用いて周波数解析を行うものとすれば、DFTにおいてはデータ点数が2のべき乗である必要はないので、FFTの場合と同じ条件でDFT計算を行ったときの複素乗算の回数は、
DFT複素乗算回数: 10000×10000=100000000)
となる。DFTでは、必要とする複素乗算回数が1×108回と極めて膨大となり、本発明のような周波数解析の用途には実用的ではない。FFTを用いる場合であっても、膨大な演算が必要となり、FFT演算器に高い処理能力を要求するとともに長大な処理時間を必要とし、その結果、強力なDSP(デジタル信号プロセッサ)やソフト演算能力を備えたCPUが必要となり、周波数検出装置のコストを上昇させる。
【0044】
図6は、潮流における流速プロファイル測定を行うために、観測船からθ=60°の方向(水平方向から60°下を向いた方向)に対して音波パルスを発射して反射エコーを受信したときの、受波波形の一例を示している。図において横軸は時刻を表している。ここでは反射エコー16は、中間周波数信号に変換されたのちの信号として示されている。時刻0の直後の振幅の大きな信号17は、送信した音波パルス自体を送受波器が受信したことによるものである。また、振幅の大きい波形18は、海底によって反射されたエコー(海底エコー)を示している。海底エコーにも、観測船の船速に応じ、±Δfp/2の周波数範囲内でのドップラーシフトが含まれる。
【0045】
ここで120±1.2(kHz)の帯域において周波数分解能12Hzでピークパワースペクトルの検出を行う場合を考えると、周波数変換の容易さ及び高次アナログフィルタの実現の容易さを考慮すれば、上述した考察により、中間周波数fmid(すなわち周波数変換された信号の中心周波数)は17kHzとなり、アナログ/デジタル変換でのサンプリング周波数及びFFTでのサンプリング周波数fftは160kHzとなり、フーリエ変換でのサンプル数Nは、N=214=16384となる。このとき、周波数分解能は9.766Hzとなって、上述した速度精度の0.1m/sに対応する周波数分解能12Hzよりも細かく、要求事項を満足している。
【0046】
図6において、符号20は、音波パルスの送信と同時に取り込みを開始するとして、FFTに必要な16384個の信号サンプルを取り込むのに要する時間T0を表している。式1−5からT0はT0=102.4msとなる。このT0の期間内には、海底エコーに起因する波形18の他に、送信した音波パルス自体の回り込みによる信号17が含まれている。振幅強度としては、信号17の方が波形18よりも大きいから、時間幅T0内の信号サンプルを取り込んでFFTを行った結果においても、最もパワーが大きな送信信号自体がピークパワースペクトルとして検出されることとなり、海底エコーに含まれるドップラーシフトした周波数をピークパワースペクトルとして検出することはできない。これは、102.4msという時間幅T0を対象としてFFTを行うことは、実質的には、送受波器の位置から海底までをまとめてFFTの対象としているからである。つまり、所望する時間範囲以外に強い信号があるFFT入力信号では、その所望する時間範囲内でのピークパワーとなるスペクトルあるいは周波数を求めることはできない。
【0047】
図6において、音波パルスの送信時刻を始点として時間t1(22)だけシフトさせ、ピーク18だけを含む時間範囲T1(21)の信号だけを切り出してFFTを行うことも考えられる。T1=6.4msであるとすると、周波数分解能f0は、f0=1/0.0064≒156Hzであり、これは、FFTステップが1つずれれば156Hz異なるという、粗い精度である。この周波数分解能は、要求されている分解能12Hzと比べて極めて粗く、要求に応えることができない。
【0048】
従来のドップラー計測では、所望の時間間隔内でピークパワーとなる周波数を求める場合、まず、距離(時間)分解能を優先させて、受信信号から必要な範囲T1を切り出し、周波数分解能が粗いままでFFTによる周波数解析を実行する。解析結果において最も高い値となったサンプリング点、次に大きな値となったサンプリング点、また次に大きな値となったサンプリング点などから、スペクトルにおけるピークの形状を仮定したカーブフィッティングなどの手法によって補間を行うことにより、真のピークパワー周波数を推定する方法などを用いていた。補間によるので、真のピークパワー周波数は、FFTにおける周波数分解能よりも細かい分解能で算出できる。
【0049】
しかしながら、このような従来のやり方では、補間によってピークパワー周波数を求めるための回路や推定アルゴリズムが複雑になり、その結果、FFT以外の処理にも演算パワーを費やす必要が生じ、演算速度の低下やコストアップの要因となっていた。また、FFTの後で周波数推定を行うため、測定精度が悪化しやすいという問題もある。
【0050】
結局、ドップラー計測における従来のピークパワースペクトル検出の手法では、ドップラー周波数すなわちドップラーシフトの帯域幅は小さくても、A/D変換におけるエリアシングを防ぐためのアナログフィルタを容易に実現し、周波数変換を確実に行うことの必要性から、サンプリング周波数を高く設定する必要があって、その分、サンプル数が大きくなり、フーリエ変換の複素乗算回数が膨大になり、演算処理が困難になることがあった。また膨大な複素演算を行うために超高速な演算処理器を必要とし、多大なコストが掛かる問題があった。
【0051】
また、フーリエ変換では、周波数分解能と測定する範囲の距離(時間)分解能は両立しない、という基本的な問題がある。従来のピークパワースペクトル検出方法では、これを解決するために、距離(時間)分解能を優先させて必要な範囲を切り出して粗い周波数分解能のFFTを行い、解析結果において最も高い値となったサンプリング点、次に大きな値となったサンプリング点、また次に大きな値となったサンプリング点などから、スペクトルにおけるピークの形状を仮定したカーブフィッティングなどの手法によって補間を行うことにより、真のピークパワー周波数を推定していた。補間によるので、真のピークパワー周波数は、FFTにおける周波数分解能よりも細かい分解能で算出できる。
【0052】
しかしながら、このような従来のやり方では、補間によってピークパワー周波数を求めるための回路や推定アルゴリズムが複雑になり、その結果、FFT以外の処理にも演算パワーを費やす必要が生じ、演算速度の低下やコストアップの要因となっていた。また、FFTの後で周波数推定を行うため、測定精度が悪化しやすいという問題もある。
【0053】
特許文献1には、FFTにおける演算ポイント数を減らしつつ、周波数分解能を高めるようにした信号処理方法が示されている。
【先行技術文献】
【特許文献】
【0054】
【特許文献1】WO2006/043511
【発明の概要】
【発明が解決しようとする課題】
【0055】
反射エコーを受信して得られた受信信号に対してドップラー周波数を求めるためなどにフーリエ変換を実行する場合、サンプリング周波数が高くなってフーリエ変換での複素乗算回数が膨大となり、演算処理が困難になることがあり、また、周波数分解能を向上させることと測定対象とする範囲の時間(あるいは距離)についての分解能を向上させることとが両立しない、というフーリエ変換の基本的な問題に起因する問題が生じる。
【0056】
本発明の目的は、フーリエ変換における上述した問題を解決し、測定対象とする時間範囲の分解能を高めつつ、周波数分解能を高め、しかもフーリエ変換におけるサンプリング周波数を低くしてサンプル数を少なくし、これによってフーリエ変換での複素乗算回数を減らすことができる検出方法及び装置を提供することにある。
【課題を解決するための手段】
【0057】
本発明の検出方法は、時間変化する受信信号における第1の周波数を含む所定の周波数帯域幅内で、かつ所望の時間幅の範囲内で、受信信号のピークパワー周波数を検出する方法であって、第1の周波数とは異なる中間周波数の信号に受信信号を周波数変換する段階と、アナログフィルタを適用して、周波数変換された信号から高域成分を除去する段階と、アナログフィルタを適用した後の信号に対して第2の周波数でサンプリングしてA/D変換する段階と、A/D変換で得られたデジタル信号に対してデジタルバンドパスフィルタを適用し、所定の周波数帯域幅に相当する信号のみを抽出する段階と、第2の周波数の2のべき乗分の1の周波数である第3の周波数をサンプリング周波数として、所望の時間幅内において、デジタルバンドパスフィルタから出力される信号をダウンサンプリングして抽出するダウンサンプリング段階と、必要とされる周波数分解能を満たすために高速フーリエ変換において必要なサンプル数を基準サンプル数として、この基準サンプル数に達するように、ダウンサンプリング段階から出力されるデジタルデータ列にゼロデータを付加する付加段階と、ゼロデータが付加されたデジタルデータ列に対して高速フーリエ変換を行うFFT段階と、を有し、ダウンサンプリング段階において、ダウンサンプリングによって抽出される信号の周波数帯域が周波数ゼロから所定の周波数帯域幅に相当する帯域幅内に配置され、高速フーリエ変換の結果において最大値を示す周波数をピークパワー周波数として検出する。
【0058】
本発明の検出装置は、時間変化する受信信号における第1の周波数を含む所定の周波数帯域幅内で、かつ所望の時間幅の範囲内で、受信信号のピークパワー周波数を検出する装置であって、受信信号を第1の周波数とは異なる中間周波数の信号に周波数変換する周波数変換手段と、周波数変換された信号から高域成分を除去するアナログフィルタと、アナログフィルタから出力される信号を第2の周波数でサンプリングしてA/D変換するA/D変換器と、A/D変換器の出力に接続され、所定の周波数帯域幅に相当する信号のみを含むデジタル信号を出力するデジタルバンドパスフィルタと、第2の周波数の2のべき乗分の1の周波数である第3の周波数をサンプリング周波数として、所望の時間幅内において、デジタルバンドパスフィルタから出力される信号をダウンサンプリングして抽出するダウンサンプリング手段と、必要とされる周波数分解能を満たすために高速フーリエ変換において必要なサンプル数を基準サンプル数として、この基準サンプル数に達するように、ダウンサンプリング手段から出力されるデジタルデータ列にゼロデータを付加するゼロ付加手段と、ゼロデータが付加されたデジタルデータ列に対して高速フーリエ変換を行うFFT手段と、を有し、ダウンサンプリング手段から出力されるデジタルデータ列における信号の周波数帯域が周波数ゼロから所定の周波数帯域幅に相当する帯域幅内に配置され、高速フーリエ変換の結果において最大値を示す周波数をピークパワー周波数として検出する。
【発明の効果】
【0059】
本発明によれば、FFT後の補間や推定などを行うことなく、所望の周波数分解能でFFTからピークパワー周波数を直接求めることができ、また、FFTにおけるサンプル数を大幅に減らして演算負荷を小さくすることができるので、複雑な回路やアルゴリズムを必要とせずに、シンプルな構成で、低コストで、精度良く、かつ高速に、ピークパワー周波数を求めることができるようになる。
【図面の簡単な説明】
【0060】
【図1】FFTによる周波数スペクトル解析結果を示すグラフである。
【図2】従来のピークパワースペクトル測定装置の構成の一例を示すブロック図である。
【図3】受信信号の中間周波数への変換などを説明するスペクトル図である。
【図4】中間周波数への周波数変換を行う場合に用いられるアナログLPFの特性の一例を示すグラフである。
【図5】A/D変換でのサンプリング周波数と中間周波数との関係を示すスペクトル図である。
【図6】斜め方向に海中に音響信号を送波したときのその反射エコーの受信波形の一例を示す波形図である。
【図7】本発明の実施の一形態のピークパワースペクトル検出装置の構成を示すブロック図である。
【図8】図7に示した装置における受信信号の中間周波数信号への変換を説明するスペクトル図である。
【図9】図7に示す装置におけるFFTの対象となる信号を示すスペクトル図である。
【図10】2のべき乗でダウンサンプリングした時の信号を説明するスペクトル図である。
【図11】デジタルBPF(バンドパスフィルタ)により不要な周波数成分を削除した後の信号を示すのスペクトル図である。
【図12】中間周波数とA/D変換器のサンプリング周波数との関係を示すスペクトル図である。
【図13】FFTへの入力信号を説明する波形図である。
【発明を実施するための形態】
【0061】
図7に示した本発明の実施の一形態のピークパワースペクトル検出装置は、フーリエ変換でのサンプル数を最小にして複素乗算回数を少なくすることでFFT演算時間の短縮を図り、かつフーリエ変換の基本的な問題である、周波数分解能と測定する範囲の距離(時間)分解能は両立しない、という問題を解決したものである。まず、図7に示す装置において、どのようにしてフーリエ変換でのサンプル数を少なくしたかについて説明する。
【0062】
媒質に対してある決まった周波数の音波信号や電波信号を発射(送信)して媒質中の対象物(ターゲット)で反射させ、反射してきた信号(反射エコー)を受信する場合を考える。ドップラー効果によって反射エコーの周波数は元の音波信号あるいは電波信号の周波数からずれるが、そのずれの大きさすなわちドップラー周波数の周波数帯域幅は、対象物との間の相対的な移動速度が媒質中での信号の伝搬速度に比べて十分に小さい場合には、元の信号の周波数に比べて十分に小さくなる。例えば、θ=60°とし、ftx=120kHzの超音波信号を水中に反射し、最大測定速度Vを15m/sとしたとき、上記の式1−1から、ドップラー周波数fdopの最大値は1200Hzとなる。対象物の移動方向は近づく方向の場合と遠ざかる方向の場合があるから、ドップラー周波数の周波数帯域幅Δfpは2400Hzとなるが、これは120kHzである送信周波数ftxに比べて十分に小さい。周波数帯域幅Δfpに隣接する減衰域を考慮しないとすれば、受信信号からこのΔfpの範囲内の周波数成分のみを取り出して0から2400Hzの範囲内に配置し、これに対してフーリエ変換を行えば、十分に小さなサンプリング周波数で十分な周波数分解能が得られるはずであり、必要な周波数分解能を得るために最小のサンプル数でフーリエ変換を実行できるはずである。従来技術においては、きわめて低い周波数への周波数変換の困難さや、急峻な特性で高域成分を除去するためのアナログフィルタの実現の困難さのために、最小のサンプル数でフーリエ変換を行って周波数解析を行うことができなかった。
【0063】
ここで注意すべきことは、式1−6に示したように、フーリエ変換による周波数分解能f0自体は、フーリエ変換の基本波の周期T0すなわち図6に示すように信号の切り出しの長さT0の逆数で表され、サンプリング周波数fftには依存しないことである。周波数分解能f0を一定としてサンプリング周波数fftが高くなれば、式1−7に基づいてサンプル数Nも増加して演算量が増大するが、増加したサンプル数N自体は、周波数解析によって検出可能な周波数範囲をより高周波側に広げることには寄与するが、周波数分解能の向上には寄与しない。サンプリング周波数が高いということは、ドップラー周波数の考え得る変化範囲を超えて、無駄に広い周波数範囲に対して実質的に周波数解析を行っていることを意味する。
【0064】
そこで図7に示すピークパワースペクトル検出装置は、オーバーサンプリングによってA/D変換を行った後、デジタルバンドパスフィルタ(D−BPF)を通して必要とする周波数帯域を取り出し、その後、ダウンサンプリングを行うことで、受信信号における周波数解析の対象となる周波数帯域を周波数ゼロから最低周波数範囲に実質的に移動させている。これにより、周波数解析の対象となる周波数帯域幅の2倍に相当するサンプリング周波数でフーリエ変換を行うことが可能となり、最小のサンプル数でフーリエ変換を行って周波数解析を行うことが可能になる。つまり、本実施形態のピークパワースペクトル検出装置によれば、所望する周波数分解能によって、周波数解析の対象としたい周波数帯域に対するパワースペクトル解析を行うことができ、かつ、サンプル数が最小であるので、最速のフーリエ変換を行うことができる。
【0065】
次に、図7に示す装置において、どのようにして周波数分解能と距離(時間)分解能とを両立させたかについて説明する。
【0066】
式1−6に示すように、フーリエ変換による周波数分解能f0自体は、フーリエ変換の基本波の周期T0すなわち図6に示すように信号の切り出しの長さT0の逆数で表され、サンプル数には依存しない。一方、距離(時間)分解能は、信号の切り出しの長さによって決定する。そこで本実施形態では、オーバーサンプリングによるA/D変換とデジタルBPFによる周波数帯域抽出とを行ってダウンサンプリングを行う際に、要求される距離分解能に対応する時間幅Td(ただし、Td<T0)においてのみダウンサンプリングを行ってサンプルを抽出してデジタルデータ列とし、その後、周波数分解能を満たすようにそのデジタルデータ列に値がゼロであるサンプルを付加し、ゼロデータを付加したことにより時間幅がT0となったデジタルデータ列に対してFFTを行うようにしている。
【0067】
以下、図7に示したピークパワースペクトル検出装置について説明する。
【0068】
図7に示したピークパワースペクトル検出装置は、図2に示す装置において高次アナログフィルタ60及びA/D変換器10と取り除き、その代わりに、オーバーサンプリングにより変調器5の出力をデジタル信号に変換するA/D変換器61と、A/D変換器61から出力されるデジタル信号から、後述するように中間周波数fmidを中心とする所定の周波数帯域Δfpの成分のみを抽出するデジタルBPF62と、デジタルBPF62から出力されるデジタル信号のダウンサンプリングを行うダウンサンプリング部63とを設け、ダウンサンプリング部63から出力されるデジタルデータ列が信号切り出しゲート11に供給されるようにしたものである。さらに、信号切り出しゲート11とFFT処理部13との間に、ゼロ付加部12を挿入したものである。信号切り出しゲート11は、ダウンサンプリング部63からのデジタルデータ列から、測定対象の時間範囲Tdに対応して部分のみを切り出し、ゼロ付加部12は、必要とされる周波数分解能を満たすためにFFTにおいて必要なサンプル数を基準サンプル数として、この基準サンプル数に達するように、信号切り出しゲート11からのデジタルデータ列にゼロデータを付加する。ゼロ付加部12での処理の詳細については後述する。
【0069】
最初に、送受波器1からダウンサンプリング部63までの動作について説明する。
【0070】
このピークパワースペクトル検出装置では、水中に超音波パルスを発射する時には、送信周波数ftxの信号を送信回路3から送受信切替え回路2を介して送受波器1に供給する、その後、直ちに送受信切替え回路2を受信側とし、媒質(水中)からの反射エコーを送受波器1で受信して電気信号である受信信号に変換し、送受信切替え回路2を介して受信信号を受信増幅器4に供給する。その結果、受信信号は受信増幅器4において信号増幅され、変調器5において中間周波数信号に変換される。
【0071】
図8は、受信信号と、局部発振回路6から変調器5に供給される局部発振周波数floc(45)の信号と周波数変換後の中間周波数fmid(46)の信号との関係を示している。受信増幅器2から出力される信号自体は、符号52で示されるような周波数帯域を有するものとする。受信信号の帯域51は、符号52で示す周波数帯域の一部であり、受信信号の周波数帯域の中心は、測定しようとするドップラー周波数の中心となる。これは、ドップラー周波数がゼロの時の周波数fc(47)であって、送信周波数ftxと同じである(fc=ftx)。変調器5は、受信信号と局部発振周波数flocの信号とを乗算して、受信信号を中間周波数fmidの信号に周波数変換する。ここでは、上述の図3の場合と同様に、局部発振周波数flocは受信信号よりも高い周波数に設定される、すなわちfloc>fc+Δfp/2であるものとする。受信信号の上限周波数をfmax(49)とし、下限周波数をfmin(48)とすると、フーリエ変換を実行すべき周波数範囲Δfp(51)は、Δfp=fmax−fminで表される。本実施形態のピークパワースペクトル検出装置では、フーリエ変換における最小サンプル数を実現するために、図9に示すように、ドップラー周波数の変化範囲Δfpを最終的に最小の周波数範囲に収めるようにする。言い換えれば、フーリエ変換の対象とする最高周波数をfbsとして、fbs=Δfpとなるように、周波数ゼロから周波数fbsまでの範囲内に周波数解析の範囲Δfpを収めるようにする。言い換えれば、フーリエ変換の対象とする最高周波数をfbsとして、fbs=Δfpとなるように、周波数ゼロから周波数fbsまでの範囲内に周波数解析の範囲Δfpを収めるようにする。以下、Δfpの範囲をこのように配置することについて説明する。
【0072】
図9のスペクトル図において、所望する周波数解析の範囲はΔfp(71)であり、その最高周波数はfbs(69)であるから、フーリエ変換のサンプリング周波数fftは、これの2倍でよい。なお、範囲Δfp(71)の中心周波数fmfが符号70で表されている。(サンプリング周波数fftとして、Δfpの2倍よりは少し高めの周波数を使用してもよく、周波数範囲Δfpは、A/D変換の時にエリアシングが発生しない範囲にあればよい。
【0073】
まず、中間周波数fmidについて説明する。図10は、2のべき乗でダウンサンプリングした時の信号を説明する図であって、受信信号を中間周波数fmidの信号に周波数変換したときのΔfpの帯域53と、目的とする最小周波数帯に投影したときのΔfpの帯域71とを示している。帯域71の最高周波数69をfbsとすると、fmidを次の式に示すように定めればよい。
【0074】
【数7】
【0075】
fbs(69)をΔfpの帯域幅と同じに設定するとと、
【0076】
【数8】
【0077】
が得られる。このようにして中間周波数fmidを求めたら、次に、局部発振周波数flocを決定する。図8の左側部分に示す中間周波数fmid(50)は、右側部分のfmid(46)と等しく、これはfmid=floc−fcであるから、
floc=fc+fmid=ftx+3.5×Δfp :式3−3
とflocを定めればよい。
【0078】
次に、A/D変換器61でのオーバーサンプリングによるサンプリング周波数fadcについて説明する。図11は、受信信号を中間周波数fmidの信号に周波数変換し、A/D変換器61によりデジタル信号に変換し、その後、デジタルBPF62により不要な周波数成分を削除してΔfpの範囲のみを抽出した信号53と、A/D変換器61のサンプリング周波数fadc(64)との関係を示している。ここで、中間周波数fmidの信号に変換した後のΔfpの範囲の上限周波数66をfmmとおくと、後述するダウンサンプリングを考えて、図10や式3−2に示すように、
fmm=4×fbs :式3−4
とすればよいことが分かる。fmmは、A/D変換器61への入力信号の最高周波数となるから、mを2のべき乗(すなわちnを1以上の整数としてm=2n)として、m倍のオーバーサンプリングを行うこととすれば、式3−1から、サンプリング周波数fadcは、
fadc=fmm×m=4×m×fbs :式3−5
で表されることになる。
【0079】
次に、アナログフィルタ7に要求される特性について説明する。図12は、中間周波数fmid(50)とオーバーサンプリングによるサンプリング周波数fadc(64)との関係を示している。
【0080】
A/D変換器61に対して供給される入力信号は、サンプリング定理に従い、アナログフィルタ7によって予め高域成分が減衰されていなければならない。上述したように、fmm(66)が信号の最高周波数であるから、これをアナログフィルタ7での通過帯域の上限とし、サンプリング周波数fadc(64)の半分の周波数fadc/2(65)において信号が十分に減衰(例えば、8ビット幅A/D変換であれば−48dB)していればよい。式3−5よりfmm(66)はfadc/mでもあるから、結局、信号の最高周波数であるfmmからみてそのm/2倍の周波数において十分な減衰が確保できればよいことになる。m倍のオーバーサンプリングを行うことにより、アナログフィルタにおいて必要とされるオクターブあたりの減衰量を1/mに緩和することが可能になって、アナログフィルタの設計が容易になる。m=8とした場合には、4次のバタワース型LPFに相当する減衰特性を有するフィルタをアナログフィルタ7に使用することが可能になり、アナログフィルタ7を容易に構成することが可能になる。
【0081】
次に、デジタルBPF62について説明する。アナログフィルタ7を介してA/D変換器61に入力してオーバーサンプリングされデジタル信号に変換された信号は、図12において符号54で示すように、解析対象とする周波数範囲Δfp(53)よりも広い周波数帯域を有している。そこでデジタルBPF62は、符号54で示す周波数範囲の信号から不要な信号成分を取り除き、図11に示すように所望の帯域Δfp(53)のみを取り出してダウンサンプリング部63に供給する。デジタルBPF62としては、BPF特性を有するFIR(有限インパルス応答;finite impulse response)デジタルフィルタを用いることができる。
【0082】
次に、ダウンサンプリング部63によるダウンサンプリングについて説明する。デジタルBPF62によって処理することによって、中間周波数fmidに周波数変換された受信信号は、所望の帯域Δfp(53)の成分のみを有する信号となっている。ダウンサンプリング部63は、この信号を、mによってダウンサンプリングする。言い換えれば、m倍のオーバーサンプリングであるfadcでサンプリングされた信号を、サンプル数が2m分の1になるように、ダウンサンプリングを行う。図10及び図11の符号66は、帯域Δfpの最高周波数fmmを表しているが、式3−5から(fadc/m)=fmmであるから、1/mにダウンサンプリングした周波数fadc/mは、A/D変換時にエリアシングを発生させない最高周波数fmm(66)と等しいことになる。
【0083】
ところで、符号53で表される、中間周波数fmidに変換されたのちの帯域Δfpは、全体として、ダウンサンプリングの周波数fadc/m(66)の半分の周波数fadc/2m(67)よりは高い周波数領域にあり、ダウンサンプリングによるエリアシング領域にあることになる。その結果、ダウンサンプリングを行うことにより、符号53で表される帯域Δfpは、周波数fadc/2m(67)を中心として鏡像対称で折り返され、図12において符号71で表される領域に移ることになる。fmm(66)はこの折り返しにより、符号68で示す周波数に移ることになる。ここでは、Δfp=fbs=fadc/4mとしているので、符号68で示す周波数は0Hzとなる。したがって、ドップラー周波数測定のための周波数解析の対象となる帯域Δfp(71)は、最終的に、0Hzからfbsまでの周波数範囲に移ったことになる。
【0084】
ところで、図9及び図10において、帯域Δfp(71)の最高周波数は符号69で表されるfbsである。デジタルBPF62によるデジタルフィルタリングを行っているので、ダウンサンプリング部63からの出力には、周波数fbsを超える成分は現れない。したがって、fbsの2倍である周波数fadc/2m(67)をフーリエ変換のサンプリング周波数fftとすることができる。
【0085】
【数9】
【0086】
次に、信号切り出しゲート11からFFT処理部13までの動作を説明する。
【0087】
ダウンサンプリング部63からは、周波数fmf(70)を中心周波数とするように周波数が移された受信信号における瞬時値を表すデジタル値からなるデジタルデータ列が、クロック周波数をfadc/mとして逐次出力されている。そこで信号切り出しゲート11は、所望する時間範囲で、例えば、海中における深さ方向での潮流の流速分布を求める場合であれば、所望の深さに対応した時間範囲で、ダウンサンプリング部63から出力されるデジタルデータ列を切り出す。この切り出されたデータのサンプル数を切り出しサンプル数とする。ゼロ付加部12は、要求されている周波数分解能に対応するFFTのサンプル数を基準サンプル数として求め、基準サンプル数に達するように、基準サンプル数と切り出しサンプル数との差に相当する数のゼロデータ(値がゼロであるサンプル)を切り出されたデジタルデータ列に付加する。そして、ゼロデータが付加されてサンプル数が基準サンプル数とされたデジタルデータ列に対し、FFT処理部13がFFTを実行し、ピークパワースペクトルを算出する。算出されたピークパワースペクトルは、信号切り出しゲート11において切り出された範囲に対するものであって、信号におけるその時間範囲でのピークパワー周波数を示している。また、ゼロデータが付加されたことにより、信号においてFFTの対象となる時間幅が実質的に拡大されており、この時間幅の逆数が周波数分解能に該当することから、周波数分解能が向上している。このように、本実施形態によれば、受信信号の中での所望の時間範囲(あるいは距離範囲)のピークパワースペクトルを、FFTの演算だけで、所望の周波数精度を満たして検出することができる。
【0088】
以下、数値例を挙げて、この実施形態における周波数分解能及びFFTにおけるサンプリング周波数について説明する。
【0089】
背景技術の欄で述べたものと同じ条件で、超音波のパルス状の信号を海中に送信し、移動物からの反射エコーのドップラー周波数を測定する場合を考える。すなわち、海水中の音波伝搬速度CをC=1500(m/s)とし、送信周波数ftxをftx=120kHzとし、検出最大速度(水平方向)VをV=15(m/s)とし、音波を水平方向から斜め方向(θ=60°)に送受信するものとする。
【0090】
このとき、ドップラー信号は120±1.2(kHz)の範囲であり、速度精度として0.1m/sを要求すれば周波数分解能f0は12Hzとなる。また、ドップラー計測において検出対象とする位置の範囲(位置分解能)を7.5mとする。
【0091】
観測対象とすべき周波数帯域Δfpは、Δfp≧2×1200=2400としなければならない。また、フーリエ変換での最大周波数fbsをfbs=Δfpとする。
【0092】
信号の減衰帯域(1300Hz)を考慮して、
Δfp=(1200+1300)×2=2500×2=5000
とする。1300Hzの減衰帯域の外側では信号はほとんどゼロと考えてよいので、5kHz(すなわちfbs)においては信号はほとんどゼロである。
【0093】
式3−1及び式3−3より、
中間周波数の中心: fmid=3.5×fbs=17500
局部発振周波数: floc=ftx+fmid=120000+17500=137500
となる。また、オーバーサンプリングの係数mをm=4とすると、A/D変換のサンプリング周波数fadcは、式3−5より、
A/D変換: fadc=4×m×fbs=16×5000=80000
フーリエ変換のサンプリング周波数fftは、
【0094】
【数10】
【0095】
となり、周波数分解能f0を満たすフーリエ変換でのサンプル数Nは、N=fft/f0≒833.3となる。834以上であって834に最も近い2のべき乗数は210=1024であるので、FFTでのサンプル数Nを1024とする。サンプル数を1024としたことにより、10000/1024≒9.77Hzの周波数分解能が達成される。
【0096】
海中に音波パルスを発射して受信信号を観測し、その受信信号を中間周波数信号に変換した後の波形を示す図6において、上述したように、符号16は反射エコーを表し、符号17は送信信号自体を表し、符号18は海底エコーを表している。海底エコー18の周波数を解析し、式1−1に基づき、ドップラー周波数から船の速度を求める場合を考える。ここでは、送信信号などの不要な信号が含まれないようにするために、海底エコー18の部分だけを切り出してFFTを行い、海底からの反射波に関してドップラー周波数を求める。
【0097】
図13は、図6に示す信号波形から海底エコーだけを切り出した信号19を示している。ただし、図6はFFTでのサンプリング周波数を160kHzとしたものであるのに対し、本実施形態ではFFTでのサンプリング周波数を10kHzとしているので、図13では、時間長では図6と図13とが一致するように、サンプリング周波数の比に合わせて、横軸もサンプル数表示では1/16にしている。ここでは、海底エコー19の部分として、音波パルスの送信から時間t2(26)の経過後に時間幅T2(24)を切り出すものとする。ここで切り出された信号は、サンプル数としては50個である。
【0098】
時間幅T2の50個のサンプルのみを用いてフーリエ変換を行った場合の周波数分解能f0は、f0=10000Hz/50=200Hzとなり、要求されている周波数分解能である12Hzを全く満たさない。
【0099】
そこでこの実施形態では、FFTについては周波数分解能f0を満たすようにサンプル数N=1024で実行するものとし、海底エコーの信号19の分のサンプル数は50なので、不足した974個のサンプルについては、ゼロ値からなるサンプルを挿入する。図13においてT0(23)は、1024個のサンプルに対応する時間であり、FFT解析の対象となる時間幅である。
【0100】
ここでT0=T2+T3の関係が成り立っており、測定したい範囲の信号19を切り出してFFT解析に必要なサンプル数Nを満たすように付加されるゼロ値信号の継続時間が、T3に対応することになる。時間幅T3の挿入位置は、時間幅T2すなわち切り出された信号19の前でも後でもよい。あるいは、信号19を挟むように前後に分割してゼロ値信号を挿入する時間幅を配置してFFT解析の時間幅T0(27)を設定してもよい。
【0101】
ここで示す例では、音波パルスの送信から信号の切り出しを開始するまでの経過時間t2(図13での符号26)と時間幅T2(図13での符号24)を予め定めているが、海底エコーのような強度の大きな反射エコーを検出対象とする場合であれば、切り出しまでの経過時間t2も切り出しの時間幅T2も自動検出して決定することも可能である。深度ごとの潮流の流速を求める場合などには、自動検出を行うことができないのでt2及びT2について予め設定しておくことになる。
【0102】
ここでゼロ付加部12によるゼロ値データの付加について説明する。
【0103】
一般にFFTではその対象とするサンプル数は2のべき乗である必要がある。その一方で、FFTによる解析を行いたい対象データにおけるサンプル数が2のべき乗となっていないことはよくあることである。そこで、サンプル数が2のべき乗となるように、ゼロデータのサンプルを付加することがゼロパディングなどと称して行われている。その場合、とにかくFFTを行えるだけのサンプル数となればよいのであるから、対象データのサンプル数をNsとして、2L-1<Ns≦2Lを満たす自然数Lが存在する時に、2L−Ns個のゼロ値データを付加することになる。それに対して本実施形態では、単にFFTの実行が可能になるようにゼロ値データを付加するのではなく、FFTにおけるサンプリング周波数fftを考慮して所望の周波数分解能が得られるように、ゼロ値データを付加する。したがって、付加されるゼロ値データの数も、単にFFTを実行可能になるようにする場合に比べ、はるかに多くなることがある。上述した例では、海底エコーを切り出してサンプル数が50になったとして、単にFFTを実行するためだけであれば26=64であるので14個のゼロ値データを付加すればよいだけであるが、所望の周波数分解能を満たすようにするために、974個もののゼロ値データを付加しているのである。
【0104】
またこの実施形態では、FFTにおけるサンプル数は1024であるので、複素乗算の演算回数は10240回で済む。図2〜図5に示した従来の例では、本実施形態におけるものと同じ周波数分解能を得るために、FFTにおいて229376回の複素乗算を行わなければならないので、本実施形態によれば、同じ周波数分解能を得るために必要な複素乗算の回数が大幅に減少することがわかる。また、図2〜図5に示した従来の例では、海底エコーの部分だけを切り出してFFTを行ったときの周波数分解能は156Hzであったのに対し、本実施形態では9.77Hzの周波数分解能が維持される。
【0105】
このように本実施形態によれば、周波数分解能と位置分解能とを両立させることができ、しかも、フーリエ変換で必要なサンプル数を大幅に減らして演算量を大幅に削減することが可能になる。
【0106】
以上、本発明の好ましい実施形態について、反射エコーにおけるドップラー周波数の検出を例に挙げて説明したが、本発明の適用範囲はこれに限られるものではない。検出対象とすべき周波数範囲を想定できる信号から短い時間範囲の部分を切り出してFFTによる解析を行う場合に、本発明は広く一般的に適用可能である。
【符号の説明】
【0107】
1 送受波器
2 送受信切替え回路
3 送信回路
4 受信増幅器
5 変調器
6 局部発振回路
7 アナログフィルタ
10,61 A/D変換器
11 信号切り出しゲート
12 ゼロ付加部
13 FFT処理部
60 高次アナログフィルタ
62 デジタルBPF(バンドパスフィルタ)
63 ダウンサンプリング部
【技術分野】
【0001】
本発明は、時間に応じて変化する信号の中から所望の時間幅内の部分を切り出し、その切り出された部分に対してフーリエ変換を適用することによってピークパワースペクトルを検出する方法及び装置に関し、特に、所望の周波数分解能力を維持しながらフーリエ変換におけるサンプル数を減らして演算処理時間を短縮でき、かつ検出されるスペクトルの周波数精度を向上させることができる方法及び装置に関する。
【背景技術】
【0002】
所定の周波数帯域幅を有して時々刻々と変化する信号の中から所望の時間幅の部分のみを切り出し、その切り出された部分内での信号のピークパワースペクトルを検出することは、後述するように、例えばドップラー計測などで広く行われている。ここでいう所望の時間幅は、信号全体の継続時間に比べた場合には十分に短く、例えば、その開始時刻が与えられているかあるいは開始時刻を自動検出できるものである。また、検出されるピークパワースペクトルは、切り出された部分内での信号の周波数スペクトルの中で一番強いパワースペクトルを有するものである。パワースペクトルの検出には、一般的にフーリエ変換が用いられ、計算量の観点から、中でも、FFT(高速フーリエ変換:Fast Fourier Transform)が使用される。
【0003】
ところで、ある信号を媒質中に発射(送信)し、その媒質中から反射した信号(エコー信号)を受信したとき、送信した時刻をスタート時刻とすると、スタート時刻から反射によるエコー信号を受信するまでの時間は、媒質中でのその信号の伝搬距離すなわち対象物までの距離に依存する。したがって、信号の伝搬速度(音波であれば音速)を考慮して、スタート時刻から所望の時間経過した時点で受信されるエコー信号は、その時間経過に応じた距離の位置にある対象物からのエコー信号となる。スタート時刻からの時間(すなわち距離)に応じてエコー信号の強度を映像として表示するようにした応用機器が、例えば、レーダーや魚群探知機である。
【0004】
信号を反射する対象物が相対的に移動しているとすると、反射エコーでは、ドップラー効果により、送信した周波数からの周波数変位が観測される。この周波数差、すなわち媒質に対して送信したときの周波数と反射エコーで観測される周波数との差をドップラー周波数という。ドップラー計測では、一定の周波数の信号を媒質中に発射し、媒質中の対象物によって反射された信号(反射エコー)を受信し、その周波数解析を行って、ドップラー周波数を決定するものである。ドップラー計測を行うものとして、例えば、ドップラーレーダー、ドップラーソナー、ADCP(音波ドップラー流速プロファイラー;acoustic Doppler current profiler)、潮流計、音波レーダー、超音波ドップラー速度計などがあり、これらは、ドップラー周波数に基づいて、対象物との相対速度を決定する。
【0005】
送信回路によって一定の送信周波数ftxの信号を生成して、この信号を送受波器から媒質に発射する。媒質中を対象物(反射体)が移動しているとして、ドップラー周波数fdopは、対象物の移動速度をV(m/s)、媒質中での信号の伝搬速度をC(m/s)とすると、媒質中の信号の伝搬速度が対象物の速度よりも十分に速い場合(C≫V)には、次の近似式で表される。
【0006】
【数1】
【0007】
ここでθは、対象物から送受波器を見たときの方向と対象物の移動方向とがなす角度であり、同じ方向であれば(対象物が送受波器にむかってまっすぐ進んでいるときには)、θ=0となる。したがって、ドップラー周波数fdopは、対象物が近付いているときには正(+)となり、遠ざかるときには負(−)となる。
【0008】
エコー信号の周波数解析を行ってドップラー周波数fdopを求めるが、この周波数解析には、一般的に、フーリエ変換の手法が使われている。フーリエ変換として、パワースペクトル解析を行う場合にはDFT(離散フーリエ変換;Discrete Fourier Transform)が用いられ、あるいは信号の標本数(サンプル数)が2のべき乗個であるときには、高速演算処理が可能なFFT(高速フーリエ変換:Fast Fourier Transform)が通常用いられている。DFTもFFTも複素乗算演算を反復するものであるが、サンプル数をNとすると、DFTでの複素乗算数回数はN×N回(すなわちN2回)であり、N=2nとすると、FFTでの複素乗算回数は、2n×log2 2n回(すなわちN×log2 N回)であり、Nが大きくなるほど、FFTはDFTに比べて高速演算ができる。
【0009】
ドップラー周波数を含んだ受信信号からFFTあるいはDFTなどの手法によりフーリエ変換で周波数を求めるためには、アナログ信号である受信信号をデジタル信号に変換することが必要である。サンプリング定理から、このアナログ/デジタル(A/D)変換でのサンプリング周波数は、入力信号の2倍以上の周波数でなければならない。言い換えると、A/D変換の対象となる入力信号はサンプリング周波数の1/2以下でなければならない。
【0010】
ところで、ドップラー計測では、対象物の位置と速度の両方に関する情報を取得したい場合が多い。例えば、海洋において深さごとの潮流の方向と速度とを求めて流速プロファイルを得るためには、海中に対して斜め方向に音波パルスを放射して反射エコーを取得し、周波数解析を行う。この際、音波パルスの放射から所定の時間の経過後にある時間幅で反射エコーを切り出し、切り出された部分に対して周波数解析を行うことになる。音波パルスの発射からの所定の時間によって、どの深度の流速を計測しようかが定まることになり、反射エコーを切り出す時の時間幅によって、どの程度の深度範囲(厚さ)にわたる平均として流速が算出されるかが定まることになる。例えば、切り出す時の時間幅が大きければ、海中における厚みの大きい層での平均的な流速が求められることになる。これに対し、時間幅が小さければ、薄い層での流速を求めることが可能になって、結果として、細かな流速プロファイルが得られることになる。
【0011】
しかしながら、ドップラー周波数を求めるための周波数解析としてのフーリエ変換の問題点として、周波数分解能と時間(距離)分解能とが相反する関係にある、ということが挙げられる。ここでいう周波数分解能は、ドップラー周波数を求める際の分解能であり、これは速度の分解能に比例する。時間分解能とは、反射エコーを切り出す時にどこまで時間幅を小さくできるかということであり、これは、ドップラー計測における距離あるいは位置の分解能に比例する。時間分解能が粗ければ、上述した流速プロファイルの場合であれば、深度方向の厚みが大きな層でのおおざっぱな平均的な流速しか得られないことになる。
【0012】
離散フーリエ変換(DFT)は、次の式で表される。
【0013】
【数2】
【0014】
DFTまたはFFTでのサンプリング周波数をfft、フーリエ変換の基本波の周期をT0とすると、周期の逆数は、フーリエ変換での基本周波数f0と呼ばれ、
f0=1/T0 :式1−4
と表される。この基本周波数f0は、DFTやFFTでの最小分解周波数すなわち周波数の分解能となる。サンプルされた信号点の総数すなわちサンプル数をNとすると、
T0=N/fft :式1−5
が当然に成立し、周波数分解能である最小分解周波数f0は、
f0=1/T0=fft/N :式1−6
と表されることになる。逆に、所望の周波数分解能f0が与えられたとして、必要なサンプル数Nは、
N=fft/f0 :式1−7
と表されることになる。フーリエ変換では、基本波の周期T0と等しい時間幅の信号から等間隔でN個のサンプルを取り出すから、T0が信号の切り出し時間(切り出しの時間幅)ということになる。
【0015】
FFTで得られるパワースペクトルの例が図1に示されている。FFTのパワースペクトルは、サンプル数Nに対し、0から(N/2)−1までの一定の周波数間隔で周波数ごとの値14,15が与えられるスペクトルとして表される。ここで周波数間隔は上記のf0と一致し、f0=1/T0である。
【0016】
式1−4に示すように、周波数分解能f0は、信号の切り出し時間幅T0の逆数であるから、周波数分解能を上げる(f0を低くする)ためには、信号の切り出し時間幅T0を長くしなければならない。その結果、周波数分解能を2倍細かくすると、時間分解能すなわち位置(距離)の分解能は2倍粗くなって劣化する。
【0017】
対象物との相対的な速度によって、ドップラー周波数fdopは変化する。このドップラー周波数ffopが変化する範囲をΔfpとおくことにすると、近づく方向と遠ざかる方向での移動が同等に起こるものと仮定すれば(この過程は、潮流計などの用途において妥当な仮定である)、ドップラー信号の周波数範囲ΔfPの中心周波数fc、つまりドップラー周波数fdopがゼロのときの受信信号の周波数は、送信周波数ftxと一致する(fc=ftx)。
【0018】
以下の条件に従い、超音波のパルス状の信号を海中に送信し、移動物からの反射エコーのドップラー周波数を測定する場合を考える。
【0019】
海水中の音波伝搬速度C:C=1500(m/s)
送信周波数ftx:ftx=120kHz
検出最大速度(水平方向)V:V=15(m/s)
また、音波は水平方向から斜め方向(θ=60°)に送受信するものとする。式1−1により、
【0020】
【数3】
【0021】
となり、検出最大速度が15m/sであるとすれば、ドップラー信号として観測される周波数範囲Δfpは、120±1.2(kHz)となる。
【0022】
測定速度範囲すなわち検出最大速度は、通常、予め決められているので、ドップラー信号についても、その周波数範囲が予め定められて測定する。
【0023】
ドップラー計測による速度精度を0.1m/sとすると、それに必要な周波数分解能f0は12Hzとなる。ドップラー計測において検出対象とする位置の範囲(位置分解能)を7.5mとすると、信号の切り出し時間幅は、海水中で音波が1mを往復する時間である2/1500(m/s)を乗ずると、7.5×2/1500=0.01となって、10msとなる。
【0024】
しかしながら、周波数分解能f0である12Hzを式1−4に代入すると、T0=1/f0=1/12≒83.3msとなり、83.3msは、海水中で音波が1mを往復する時間2/1500(m/s)で割ることにより距離に換算すると、83.3(ms)×1500(m/s)/2≒62.5(m)となる。
【0025】
すなわち、速度精度を0.1m/sとするために周波数分解能f0を12Hzとすると、信号の切り出し時間幅として83msが必要であり、これは距離に換算すれば63mに相当し、位置(距離)に関する分解能を7.5mとする要求を満足することできない。逆に、位置に関する分解能が7.5mである信号でFFTを行うとすると、式1−4より周波数分解能f0が100Hzとなり、速度精度から定められた周波数分解能12Hzの要求を満足できない。
【0026】
このようにフーリエ変換には、周波数分解能と、測定対象に関する時間分解能とが両立しない、という基本的な問題があり、これは、ドップラー計測の場合には、周波数分解能あるいは速度分解能と、距離あるいは位置の分解能とが両立しないことを意味する。
【0027】
図2は、一般的なドップラー計測機器の構成の一例として、受信信号の中から所望の時間幅の部分のみを切り出し、その切り出された部分内での信号のピークパワースペクトルを検出する従来のピークパワースペクトル検出装置の構成を示している。ここでは、信号が音響信号であり、信号が伝搬する媒質が水(あるいは海水)であるとする。音響信号の送受波器1は、送受信切替え回路2を介し、送信回路3の出力と受信増幅器4の入力に接続している。送信回路3は、送信周波数ftxの信号を生成するものである。受信増幅器4の出力には、受信信号を中間周波数信号に周波数変換する変調器5が設けられており、変調器5には、局部発振回路6から局部発振周波数flocの信号が供給されている。変調器5の出力すなわち中間周波数信号は、高次アナログフィルタ60を介してアナログ/デジタル変換器(A/D)10に供給されてデジタル信号に変換され、このデジタル信号は、所定のサンプル数(すなわち時間幅)で信号を切り出す信号切り出しゲート11に供給される。信号切り出しゲート11で切り出された信号は、次にFFT処理部13に送られてFFTによる周波数解析の対象となる。FFT処理部13から、計算されたパワースペクトルが出力される。
【0028】
この装置では、水中に音波パルスを発射する時には、送信周波数ftxの信号を送信回路3から送受信切替え回路2を介して送受波器1に供給する、その後、直ちに送受信切替え回路2を受信側とし、媒質(水中)からの反射エコーを送受波器1で受信して電気信号である受信信号に変換し、送受信切替え回路2を介して受信信号を受信増幅器4に供給する。その結果、受信信号は受信増幅器4において信号増幅され、変調器5に対して供給されることになる。
【0029】
図3は、受信信号や中間周波数信号などの周波数間の関係を示している。受信増幅器4から出力される信号自体は、図3において符号52で示されるような周波数帯域を有するものとする。受信信号の帯域は、符号52で示す周波数帯域の一部であり、検出対象となるドップラー信号は、符号51で示される周波数範囲Δfpの中にある。周波数範囲Δfpの中心、すなわち、ドップラー周波数がゼロの時の周波数fc(47)は、送信周波数ftxと同じである(fc=ftx)。局部発振回路6は、周波数fcよりも中間周波数fmidだけ高い局部発振周波数floc(45)の局部発振信号を生成し(すなわちfloc>fc)、変調器5に供給する。その結果、変調器5では、局部発振信号を受信信号に乗算することにより、受信信号を中間周波数fmidの中間周波数信号に変換する。
【0030】
受信信号の周波数帯域のうち、トップラー計測のために着目する周波数範囲Δfpの上限を上限周波数fmax(49)とし、下限を下限周波数fmin(48)とすると、フーリエ変換すべき周波数範囲Δfp(51)は、Δfp=fmax−fminで表される。
【0031】
変調器5における中間周波数fmidの信号の生成は、下記の式で表されるものであり、2つの変換周波数から、高次アナログフィルタ60により、低い方の周波数を得る。
【0032】
【数4】
【0033】
ドップラー信号の周波数範囲Δfpは、中間周波数fmidを中心に、
【0034】
【数5】
【0035】
の範囲に変換されることとなり、検出対象の周波数範囲Δfpは、周波数変換しても保たれる。上記の式は、floc>fcの時の処理を示しているが、floc<fcの時も、複号(符号(±))が逆になるだけで、同様に処理がなされる。
【0036】
図3において、右側部分が受信信号と局部発振周波数flocの信号との関係を表しており、ここでfmid(46)は、fmid=floc−fcの関係を満たしている。一方、図3の左側部分では、中間周波数信号に変換後のスペクトルを示しており、ここでは中間周波数信号の中心周波数はfmid(50)となっている。また、上述の周波数帯域52も周波数帯域54に変換され、検出対象の周波数範囲Δfp(51)は、中心周波数がfmidである周波数範囲Δfp(53)に変換されている。すなわち帯域53は、変調器5によって周波数変換して低い方の周波数成分だけを取り出すことによって中間周波数信号とされ、高次アナログフィルタ60を通過させた後の受信信号におけるΔfpの帯域を表している。
【0037】
高次アナログフィルタ60は、変調器5から出力される2つの変換周波数のうちの低い方を選択し、さらに、A/D変換でのエリアシングを防ぐためのものであるので、このフィルタには、バンドパスフィルタ(BPF;band-pass filter)かローパスフィルタ(LPF;low-pass filter)が使用される。以下では、説明のため、通過帯域におけるリップルが少ないバターワース型LPFを高次アナログフィルタ60に使用するものとする。A/D変換でのエリアシングを考えると、A/D変換器10におけるA/D変換のデータ幅が8ビットであるときには、A/D変換器10のサンプリング周波数f2adの1/2の周波数すなわちf2ad/2での信号レベルが、フィルタにおける通過帯域での信号レベルを基準として−48dB以下であれば、完全にエリアシングは起こらない。同様に、16ビット幅のA/D変換の場合には、信号レベルが−96dB以下であればエリアシングは完全に起らない。したがって、高次アナログフィルタ60は、中間周波数信号における着目する周波数帯域(周波数fmidを中心とする範囲Δfp)をその通過帯域内に含むとともに、周波数f2ad/2では上述した減衰量(8ビット幅のA/D変換であれば−48dB、16ビット幅であれば−96dB)を有するものである必要がある。逆に言えば、減衰量がこの値となるまでfft/2を高くし、その分、サンプリング周波数f2adも高くする必要がある。
【0038】
以下では、A/D変換のデータ幅が8ビットであるものとして説明を進める。図4は、LPFの周波数特性の一例を示している。以下、図4を用いて、受信信号の周波数帯域を説明する。
【0039】
受信信号の周波数fcが120kHz、ドップラー周波数の最大値が1200Hz、フーリエ変換に要求される周波数分解能が12Hzであるとする。式2−3に示すように、受信信号を中間周波数信号に変換したとしても、ドップラー周波数の周波数範囲Δfpそのものは変化しない。その一方で、受信信号と局部発振周波数の信号がクロストークしたり混変調を起したりしないように、局部発振周波数flocを選ぶ必要がある。ここでは、受信信号からその周波数で1割以上離れた局部発振周波数flocを選ぶこととし、17KHz離すこととする。局部発振周波数flocはfloc=fc+17kHz=137KHzとなり、中間周波数fmidは17KHzとなる。このとき、アナログLPFの通過帯域(全帯域)Wは、
【0040】
【数6】
【0041】
とする必要がある。Δfp/2=1200Hzであり、これに適切なα(減衰域)とfmidとを加算すると、LPFは、図4に示すように、20KHzまでの信号を完全に通過させるものであればよい。図4において、符号74は2次のLPFの特性を示しており、符号73は4次のLPFの特性を示しており、符号72は8次のLPFの特性を示している。アナログLPFでは、構成部品の精度や温度係数などから、8次のLPFを構成することはコスト面からも実用的ではないので、4次のLPFを高次アナログフィルタ26に用いるものとする。
【0042】
上述したように、8ビット幅A/D変換では、サンプリング周波数の1/2の周波数での信号レベルを通過帯域に比べて−48dBとすれば、エリアシングの問題は生じない。図5は、中間周波数に変換された信号とA/D変換のサンプリング周波数の関係を示している。図5においては、サンプリング周波数は、f2ad(78)で示されている。ここでは4次のLPFを使うこととしたので、−48dBの周波数は80kHzとなり、サンプリング周波数f2adの半分の周波数f2ad/2(77)も80kHzに設定される。したがって、A/D変換は、160kHzのサンプリング周波数f2adで行なわれることになる。
【0043】
その結果、式1−6から、分解能の条件を満たすサンプル数Nは、N=f2ad/f0=160000/16=10000となる。FFTのデータ点数は2のべき乗個とすべきであることから、FFTではサンプル数を16384(=214)とし、このときの周波数分解能f0は、f0=160000/16384=9.77[Hz]となる。このとき、FFTで実行される複素乗算の回数は、
FFT複素乗算回数: 214×log2 214=229376
となる。FFTの代わりにDFTを用いて周波数解析を行うものとすれば、DFTにおいてはデータ点数が2のべき乗である必要はないので、FFTの場合と同じ条件でDFT計算を行ったときの複素乗算の回数は、
DFT複素乗算回数: 10000×10000=100000000)
となる。DFTでは、必要とする複素乗算回数が1×108回と極めて膨大となり、本発明のような周波数解析の用途には実用的ではない。FFTを用いる場合であっても、膨大な演算が必要となり、FFT演算器に高い処理能力を要求するとともに長大な処理時間を必要とし、その結果、強力なDSP(デジタル信号プロセッサ)やソフト演算能力を備えたCPUが必要となり、周波数検出装置のコストを上昇させる。
【0044】
図6は、潮流における流速プロファイル測定を行うために、観測船からθ=60°の方向(水平方向から60°下を向いた方向)に対して音波パルスを発射して反射エコーを受信したときの、受波波形の一例を示している。図において横軸は時刻を表している。ここでは反射エコー16は、中間周波数信号に変換されたのちの信号として示されている。時刻0の直後の振幅の大きな信号17は、送信した音波パルス自体を送受波器が受信したことによるものである。また、振幅の大きい波形18は、海底によって反射されたエコー(海底エコー)を示している。海底エコーにも、観測船の船速に応じ、±Δfp/2の周波数範囲内でのドップラーシフトが含まれる。
【0045】
ここで120±1.2(kHz)の帯域において周波数分解能12Hzでピークパワースペクトルの検出を行う場合を考えると、周波数変換の容易さ及び高次アナログフィルタの実現の容易さを考慮すれば、上述した考察により、中間周波数fmid(すなわち周波数変換された信号の中心周波数)は17kHzとなり、アナログ/デジタル変換でのサンプリング周波数及びFFTでのサンプリング周波数fftは160kHzとなり、フーリエ変換でのサンプル数Nは、N=214=16384となる。このとき、周波数分解能は9.766Hzとなって、上述した速度精度の0.1m/sに対応する周波数分解能12Hzよりも細かく、要求事項を満足している。
【0046】
図6において、符号20は、音波パルスの送信と同時に取り込みを開始するとして、FFTに必要な16384個の信号サンプルを取り込むのに要する時間T0を表している。式1−5からT0はT0=102.4msとなる。このT0の期間内には、海底エコーに起因する波形18の他に、送信した音波パルス自体の回り込みによる信号17が含まれている。振幅強度としては、信号17の方が波形18よりも大きいから、時間幅T0内の信号サンプルを取り込んでFFTを行った結果においても、最もパワーが大きな送信信号自体がピークパワースペクトルとして検出されることとなり、海底エコーに含まれるドップラーシフトした周波数をピークパワースペクトルとして検出することはできない。これは、102.4msという時間幅T0を対象としてFFTを行うことは、実質的には、送受波器の位置から海底までをまとめてFFTの対象としているからである。つまり、所望する時間範囲以外に強い信号があるFFT入力信号では、その所望する時間範囲内でのピークパワーとなるスペクトルあるいは周波数を求めることはできない。
【0047】
図6において、音波パルスの送信時刻を始点として時間t1(22)だけシフトさせ、ピーク18だけを含む時間範囲T1(21)の信号だけを切り出してFFTを行うことも考えられる。T1=6.4msであるとすると、周波数分解能f0は、f0=1/0.0064≒156Hzであり、これは、FFTステップが1つずれれば156Hz異なるという、粗い精度である。この周波数分解能は、要求されている分解能12Hzと比べて極めて粗く、要求に応えることができない。
【0048】
従来のドップラー計測では、所望の時間間隔内でピークパワーとなる周波数を求める場合、まず、距離(時間)分解能を優先させて、受信信号から必要な範囲T1を切り出し、周波数分解能が粗いままでFFTによる周波数解析を実行する。解析結果において最も高い値となったサンプリング点、次に大きな値となったサンプリング点、また次に大きな値となったサンプリング点などから、スペクトルにおけるピークの形状を仮定したカーブフィッティングなどの手法によって補間を行うことにより、真のピークパワー周波数を推定する方法などを用いていた。補間によるので、真のピークパワー周波数は、FFTにおける周波数分解能よりも細かい分解能で算出できる。
【0049】
しかしながら、このような従来のやり方では、補間によってピークパワー周波数を求めるための回路や推定アルゴリズムが複雑になり、その結果、FFT以外の処理にも演算パワーを費やす必要が生じ、演算速度の低下やコストアップの要因となっていた。また、FFTの後で周波数推定を行うため、測定精度が悪化しやすいという問題もある。
【0050】
結局、ドップラー計測における従来のピークパワースペクトル検出の手法では、ドップラー周波数すなわちドップラーシフトの帯域幅は小さくても、A/D変換におけるエリアシングを防ぐためのアナログフィルタを容易に実現し、周波数変換を確実に行うことの必要性から、サンプリング周波数を高く設定する必要があって、その分、サンプル数が大きくなり、フーリエ変換の複素乗算回数が膨大になり、演算処理が困難になることがあった。また膨大な複素演算を行うために超高速な演算処理器を必要とし、多大なコストが掛かる問題があった。
【0051】
また、フーリエ変換では、周波数分解能と測定する範囲の距離(時間)分解能は両立しない、という基本的な問題がある。従来のピークパワースペクトル検出方法では、これを解決するために、距離(時間)分解能を優先させて必要な範囲を切り出して粗い周波数分解能のFFTを行い、解析結果において最も高い値となったサンプリング点、次に大きな値となったサンプリング点、また次に大きな値となったサンプリング点などから、スペクトルにおけるピークの形状を仮定したカーブフィッティングなどの手法によって補間を行うことにより、真のピークパワー周波数を推定していた。補間によるので、真のピークパワー周波数は、FFTにおける周波数分解能よりも細かい分解能で算出できる。
【0052】
しかしながら、このような従来のやり方では、補間によってピークパワー周波数を求めるための回路や推定アルゴリズムが複雑になり、その結果、FFT以外の処理にも演算パワーを費やす必要が生じ、演算速度の低下やコストアップの要因となっていた。また、FFTの後で周波数推定を行うため、測定精度が悪化しやすいという問題もある。
【0053】
特許文献1には、FFTにおける演算ポイント数を減らしつつ、周波数分解能を高めるようにした信号処理方法が示されている。
【先行技術文献】
【特許文献】
【0054】
【特許文献1】WO2006/043511
【発明の概要】
【発明が解決しようとする課題】
【0055】
反射エコーを受信して得られた受信信号に対してドップラー周波数を求めるためなどにフーリエ変換を実行する場合、サンプリング周波数が高くなってフーリエ変換での複素乗算回数が膨大となり、演算処理が困難になることがあり、また、周波数分解能を向上させることと測定対象とする範囲の時間(あるいは距離)についての分解能を向上させることとが両立しない、というフーリエ変換の基本的な問題に起因する問題が生じる。
【0056】
本発明の目的は、フーリエ変換における上述した問題を解決し、測定対象とする時間範囲の分解能を高めつつ、周波数分解能を高め、しかもフーリエ変換におけるサンプリング周波数を低くしてサンプル数を少なくし、これによってフーリエ変換での複素乗算回数を減らすことができる検出方法及び装置を提供することにある。
【課題を解決するための手段】
【0057】
本発明の検出方法は、時間変化する受信信号における第1の周波数を含む所定の周波数帯域幅内で、かつ所望の時間幅の範囲内で、受信信号のピークパワー周波数を検出する方法であって、第1の周波数とは異なる中間周波数の信号に受信信号を周波数変換する段階と、アナログフィルタを適用して、周波数変換された信号から高域成分を除去する段階と、アナログフィルタを適用した後の信号に対して第2の周波数でサンプリングしてA/D変換する段階と、A/D変換で得られたデジタル信号に対してデジタルバンドパスフィルタを適用し、所定の周波数帯域幅に相当する信号のみを抽出する段階と、第2の周波数の2のべき乗分の1の周波数である第3の周波数をサンプリング周波数として、所望の時間幅内において、デジタルバンドパスフィルタから出力される信号をダウンサンプリングして抽出するダウンサンプリング段階と、必要とされる周波数分解能を満たすために高速フーリエ変換において必要なサンプル数を基準サンプル数として、この基準サンプル数に達するように、ダウンサンプリング段階から出力されるデジタルデータ列にゼロデータを付加する付加段階と、ゼロデータが付加されたデジタルデータ列に対して高速フーリエ変換を行うFFT段階と、を有し、ダウンサンプリング段階において、ダウンサンプリングによって抽出される信号の周波数帯域が周波数ゼロから所定の周波数帯域幅に相当する帯域幅内に配置され、高速フーリエ変換の結果において最大値を示す周波数をピークパワー周波数として検出する。
【0058】
本発明の検出装置は、時間変化する受信信号における第1の周波数を含む所定の周波数帯域幅内で、かつ所望の時間幅の範囲内で、受信信号のピークパワー周波数を検出する装置であって、受信信号を第1の周波数とは異なる中間周波数の信号に周波数変換する周波数変換手段と、周波数変換された信号から高域成分を除去するアナログフィルタと、アナログフィルタから出力される信号を第2の周波数でサンプリングしてA/D変換するA/D変換器と、A/D変換器の出力に接続され、所定の周波数帯域幅に相当する信号のみを含むデジタル信号を出力するデジタルバンドパスフィルタと、第2の周波数の2のべき乗分の1の周波数である第3の周波数をサンプリング周波数として、所望の時間幅内において、デジタルバンドパスフィルタから出力される信号をダウンサンプリングして抽出するダウンサンプリング手段と、必要とされる周波数分解能を満たすために高速フーリエ変換において必要なサンプル数を基準サンプル数として、この基準サンプル数に達するように、ダウンサンプリング手段から出力されるデジタルデータ列にゼロデータを付加するゼロ付加手段と、ゼロデータが付加されたデジタルデータ列に対して高速フーリエ変換を行うFFT手段と、を有し、ダウンサンプリング手段から出力されるデジタルデータ列における信号の周波数帯域が周波数ゼロから所定の周波数帯域幅に相当する帯域幅内に配置され、高速フーリエ変換の結果において最大値を示す周波数をピークパワー周波数として検出する。
【発明の効果】
【0059】
本発明によれば、FFT後の補間や推定などを行うことなく、所望の周波数分解能でFFTからピークパワー周波数を直接求めることができ、また、FFTにおけるサンプル数を大幅に減らして演算負荷を小さくすることができるので、複雑な回路やアルゴリズムを必要とせずに、シンプルな構成で、低コストで、精度良く、かつ高速に、ピークパワー周波数を求めることができるようになる。
【図面の簡単な説明】
【0060】
【図1】FFTによる周波数スペクトル解析結果を示すグラフである。
【図2】従来のピークパワースペクトル測定装置の構成の一例を示すブロック図である。
【図3】受信信号の中間周波数への変換などを説明するスペクトル図である。
【図4】中間周波数への周波数変換を行う場合に用いられるアナログLPFの特性の一例を示すグラフである。
【図5】A/D変換でのサンプリング周波数と中間周波数との関係を示すスペクトル図である。
【図6】斜め方向に海中に音響信号を送波したときのその反射エコーの受信波形の一例を示す波形図である。
【図7】本発明の実施の一形態のピークパワースペクトル検出装置の構成を示すブロック図である。
【図8】図7に示した装置における受信信号の中間周波数信号への変換を説明するスペクトル図である。
【図9】図7に示す装置におけるFFTの対象となる信号を示すスペクトル図である。
【図10】2のべき乗でダウンサンプリングした時の信号を説明するスペクトル図である。
【図11】デジタルBPF(バンドパスフィルタ)により不要な周波数成分を削除した後の信号を示すのスペクトル図である。
【図12】中間周波数とA/D変換器のサンプリング周波数との関係を示すスペクトル図である。
【図13】FFTへの入力信号を説明する波形図である。
【発明を実施するための形態】
【0061】
図7に示した本発明の実施の一形態のピークパワースペクトル検出装置は、フーリエ変換でのサンプル数を最小にして複素乗算回数を少なくすることでFFT演算時間の短縮を図り、かつフーリエ変換の基本的な問題である、周波数分解能と測定する範囲の距離(時間)分解能は両立しない、という問題を解決したものである。まず、図7に示す装置において、どのようにしてフーリエ変換でのサンプル数を少なくしたかについて説明する。
【0062】
媒質に対してある決まった周波数の音波信号や電波信号を発射(送信)して媒質中の対象物(ターゲット)で反射させ、反射してきた信号(反射エコー)を受信する場合を考える。ドップラー効果によって反射エコーの周波数は元の音波信号あるいは電波信号の周波数からずれるが、そのずれの大きさすなわちドップラー周波数の周波数帯域幅は、対象物との間の相対的な移動速度が媒質中での信号の伝搬速度に比べて十分に小さい場合には、元の信号の周波数に比べて十分に小さくなる。例えば、θ=60°とし、ftx=120kHzの超音波信号を水中に反射し、最大測定速度Vを15m/sとしたとき、上記の式1−1から、ドップラー周波数fdopの最大値は1200Hzとなる。対象物の移動方向は近づく方向の場合と遠ざかる方向の場合があるから、ドップラー周波数の周波数帯域幅Δfpは2400Hzとなるが、これは120kHzである送信周波数ftxに比べて十分に小さい。周波数帯域幅Δfpに隣接する減衰域を考慮しないとすれば、受信信号からこのΔfpの範囲内の周波数成分のみを取り出して0から2400Hzの範囲内に配置し、これに対してフーリエ変換を行えば、十分に小さなサンプリング周波数で十分な周波数分解能が得られるはずであり、必要な周波数分解能を得るために最小のサンプル数でフーリエ変換を実行できるはずである。従来技術においては、きわめて低い周波数への周波数変換の困難さや、急峻な特性で高域成分を除去するためのアナログフィルタの実現の困難さのために、最小のサンプル数でフーリエ変換を行って周波数解析を行うことができなかった。
【0063】
ここで注意すべきことは、式1−6に示したように、フーリエ変換による周波数分解能f0自体は、フーリエ変換の基本波の周期T0すなわち図6に示すように信号の切り出しの長さT0の逆数で表され、サンプリング周波数fftには依存しないことである。周波数分解能f0を一定としてサンプリング周波数fftが高くなれば、式1−7に基づいてサンプル数Nも増加して演算量が増大するが、増加したサンプル数N自体は、周波数解析によって検出可能な周波数範囲をより高周波側に広げることには寄与するが、周波数分解能の向上には寄与しない。サンプリング周波数が高いということは、ドップラー周波数の考え得る変化範囲を超えて、無駄に広い周波数範囲に対して実質的に周波数解析を行っていることを意味する。
【0064】
そこで図7に示すピークパワースペクトル検出装置は、オーバーサンプリングによってA/D変換を行った後、デジタルバンドパスフィルタ(D−BPF)を通して必要とする周波数帯域を取り出し、その後、ダウンサンプリングを行うことで、受信信号における周波数解析の対象となる周波数帯域を周波数ゼロから最低周波数範囲に実質的に移動させている。これにより、周波数解析の対象となる周波数帯域幅の2倍に相当するサンプリング周波数でフーリエ変換を行うことが可能となり、最小のサンプル数でフーリエ変換を行って周波数解析を行うことが可能になる。つまり、本実施形態のピークパワースペクトル検出装置によれば、所望する周波数分解能によって、周波数解析の対象としたい周波数帯域に対するパワースペクトル解析を行うことができ、かつ、サンプル数が最小であるので、最速のフーリエ変換を行うことができる。
【0065】
次に、図7に示す装置において、どのようにして周波数分解能と距離(時間)分解能とを両立させたかについて説明する。
【0066】
式1−6に示すように、フーリエ変換による周波数分解能f0自体は、フーリエ変換の基本波の周期T0すなわち図6に示すように信号の切り出しの長さT0の逆数で表され、サンプル数には依存しない。一方、距離(時間)分解能は、信号の切り出しの長さによって決定する。そこで本実施形態では、オーバーサンプリングによるA/D変換とデジタルBPFによる周波数帯域抽出とを行ってダウンサンプリングを行う際に、要求される距離分解能に対応する時間幅Td(ただし、Td<T0)においてのみダウンサンプリングを行ってサンプルを抽出してデジタルデータ列とし、その後、周波数分解能を満たすようにそのデジタルデータ列に値がゼロであるサンプルを付加し、ゼロデータを付加したことにより時間幅がT0となったデジタルデータ列に対してFFTを行うようにしている。
【0067】
以下、図7に示したピークパワースペクトル検出装置について説明する。
【0068】
図7に示したピークパワースペクトル検出装置は、図2に示す装置において高次アナログフィルタ60及びA/D変換器10と取り除き、その代わりに、オーバーサンプリングにより変調器5の出力をデジタル信号に変換するA/D変換器61と、A/D変換器61から出力されるデジタル信号から、後述するように中間周波数fmidを中心とする所定の周波数帯域Δfpの成分のみを抽出するデジタルBPF62と、デジタルBPF62から出力されるデジタル信号のダウンサンプリングを行うダウンサンプリング部63とを設け、ダウンサンプリング部63から出力されるデジタルデータ列が信号切り出しゲート11に供給されるようにしたものである。さらに、信号切り出しゲート11とFFT処理部13との間に、ゼロ付加部12を挿入したものである。信号切り出しゲート11は、ダウンサンプリング部63からのデジタルデータ列から、測定対象の時間範囲Tdに対応して部分のみを切り出し、ゼロ付加部12は、必要とされる周波数分解能を満たすためにFFTにおいて必要なサンプル数を基準サンプル数として、この基準サンプル数に達するように、信号切り出しゲート11からのデジタルデータ列にゼロデータを付加する。ゼロ付加部12での処理の詳細については後述する。
【0069】
最初に、送受波器1からダウンサンプリング部63までの動作について説明する。
【0070】
このピークパワースペクトル検出装置では、水中に超音波パルスを発射する時には、送信周波数ftxの信号を送信回路3から送受信切替え回路2を介して送受波器1に供給する、その後、直ちに送受信切替え回路2を受信側とし、媒質(水中)からの反射エコーを送受波器1で受信して電気信号である受信信号に変換し、送受信切替え回路2を介して受信信号を受信増幅器4に供給する。その結果、受信信号は受信増幅器4において信号増幅され、変調器5において中間周波数信号に変換される。
【0071】
図8は、受信信号と、局部発振回路6から変調器5に供給される局部発振周波数floc(45)の信号と周波数変換後の中間周波数fmid(46)の信号との関係を示している。受信増幅器2から出力される信号自体は、符号52で示されるような周波数帯域を有するものとする。受信信号の帯域51は、符号52で示す周波数帯域の一部であり、受信信号の周波数帯域の中心は、測定しようとするドップラー周波数の中心となる。これは、ドップラー周波数がゼロの時の周波数fc(47)であって、送信周波数ftxと同じである(fc=ftx)。変調器5は、受信信号と局部発振周波数flocの信号とを乗算して、受信信号を中間周波数fmidの信号に周波数変換する。ここでは、上述の図3の場合と同様に、局部発振周波数flocは受信信号よりも高い周波数に設定される、すなわちfloc>fc+Δfp/2であるものとする。受信信号の上限周波数をfmax(49)とし、下限周波数をfmin(48)とすると、フーリエ変換を実行すべき周波数範囲Δfp(51)は、Δfp=fmax−fminで表される。本実施形態のピークパワースペクトル検出装置では、フーリエ変換における最小サンプル数を実現するために、図9に示すように、ドップラー周波数の変化範囲Δfpを最終的に最小の周波数範囲に収めるようにする。言い換えれば、フーリエ変換の対象とする最高周波数をfbsとして、fbs=Δfpとなるように、周波数ゼロから周波数fbsまでの範囲内に周波数解析の範囲Δfpを収めるようにする。言い換えれば、フーリエ変換の対象とする最高周波数をfbsとして、fbs=Δfpとなるように、周波数ゼロから周波数fbsまでの範囲内に周波数解析の範囲Δfpを収めるようにする。以下、Δfpの範囲をこのように配置することについて説明する。
【0072】
図9のスペクトル図において、所望する周波数解析の範囲はΔfp(71)であり、その最高周波数はfbs(69)であるから、フーリエ変換のサンプリング周波数fftは、これの2倍でよい。なお、範囲Δfp(71)の中心周波数fmfが符号70で表されている。(サンプリング周波数fftとして、Δfpの2倍よりは少し高めの周波数を使用してもよく、周波数範囲Δfpは、A/D変換の時にエリアシングが発生しない範囲にあればよい。
【0073】
まず、中間周波数fmidについて説明する。図10は、2のべき乗でダウンサンプリングした時の信号を説明する図であって、受信信号を中間周波数fmidの信号に周波数変換したときのΔfpの帯域53と、目的とする最小周波数帯に投影したときのΔfpの帯域71とを示している。帯域71の最高周波数69をfbsとすると、fmidを次の式に示すように定めればよい。
【0074】
【数7】
【0075】
fbs(69)をΔfpの帯域幅と同じに設定するとと、
【0076】
【数8】
【0077】
が得られる。このようにして中間周波数fmidを求めたら、次に、局部発振周波数flocを決定する。図8の左側部分に示す中間周波数fmid(50)は、右側部分のfmid(46)と等しく、これはfmid=floc−fcであるから、
floc=fc+fmid=ftx+3.5×Δfp :式3−3
とflocを定めればよい。
【0078】
次に、A/D変換器61でのオーバーサンプリングによるサンプリング周波数fadcについて説明する。図11は、受信信号を中間周波数fmidの信号に周波数変換し、A/D変換器61によりデジタル信号に変換し、その後、デジタルBPF62により不要な周波数成分を削除してΔfpの範囲のみを抽出した信号53と、A/D変換器61のサンプリング周波数fadc(64)との関係を示している。ここで、中間周波数fmidの信号に変換した後のΔfpの範囲の上限周波数66をfmmとおくと、後述するダウンサンプリングを考えて、図10や式3−2に示すように、
fmm=4×fbs :式3−4
とすればよいことが分かる。fmmは、A/D変換器61への入力信号の最高周波数となるから、mを2のべき乗(すなわちnを1以上の整数としてm=2n)として、m倍のオーバーサンプリングを行うこととすれば、式3−1から、サンプリング周波数fadcは、
fadc=fmm×m=4×m×fbs :式3−5
で表されることになる。
【0079】
次に、アナログフィルタ7に要求される特性について説明する。図12は、中間周波数fmid(50)とオーバーサンプリングによるサンプリング周波数fadc(64)との関係を示している。
【0080】
A/D変換器61に対して供給される入力信号は、サンプリング定理に従い、アナログフィルタ7によって予め高域成分が減衰されていなければならない。上述したように、fmm(66)が信号の最高周波数であるから、これをアナログフィルタ7での通過帯域の上限とし、サンプリング周波数fadc(64)の半分の周波数fadc/2(65)において信号が十分に減衰(例えば、8ビット幅A/D変換であれば−48dB)していればよい。式3−5よりfmm(66)はfadc/mでもあるから、結局、信号の最高周波数であるfmmからみてそのm/2倍の周波数において十分な減衰が確保できればよいことになる。m倍のオーバーサンプリングを行うことにより、アナログフィルタにおいて必要とされるオクターブあたりの減衰量を1/mに緩和することが可能になって、アナログフィルタの設計が容易になる。m=8とした場合には、4次のバタワース型LPFに相当する減衰特性を有するフィルタをアナログフィルタ7に使用することが可能になり、アナログフィルタ7を容易に構成することが可能になる。
【0081】
次に、デジタルBPF62について説明する。アナログフィルタ7を介してA/D変換器61に入力してオーバーサンプリングされデジタル信号に変換された信号は、図12において符号54で示すように、解析対象とする周波数範囲Δfp(53)よりも広い周波数帯域を有している。そこでデジタルBPF62は、符号54で示す周波数範囲の信号から不要な信号成分を取り除き、図11に示すように所望の帯域Δfp(53)のみを取り出してダウンサンプリング部63に供給する。デジタルBPF62としては、BPF特性を有するFIR(有限インパルス応答;finite impulse response)デジタルフィルタを用いることができる。
【0082】
次に、ダウンサンプリング部63によるダウンサンプリングについて説明する。デジタルBPF62によって処理することによって、中間周波数fmidに周波数変換された受信信号は、所望の帯域Δfp(53)の成分のみを有する信号となっている。ダウンサンプリング部63は、この信号を、mによってダウンサンプリングする。言い換えれば、m倍のオーバーサンプリングであるfadcでサンプリングされた信号を、サンプル数が2m分の1になるように、ダウンサンプリングを行う。図10及び図11の符号66は、帯域Δfpの最高周波数fmmを表しているが、式3−5から(fadc/m)=fmmであるから、1/mにダウンサンプリングした周波数fadc/mは、A/D変換時にエリアシングを発生させない最高周波数fmm(66)と等しいことになる。
【0083】
ところで、符号53で表される、中間周波数fmidに変換されたのちの帯域Δfpは、全体として、ダウンサンプリングの周波数fadc/m(66)の半分の周波数fadc/2m(67)よりは高い周波数領域にあり、ダウンサンプリングによるエリアシング領域にあることになる。その結果、ダウンサンプリングを行うことにより、符号53で表される帯域Δfpは、周波数fadc/2m(67)を中心として鏡像対称で折り返され、図12において符号71で表される領域に移ることになる。fmm(66)はこの折り返しにより、符号68で示す周波数に移ることになる。ここでは、Δfp=fbs=fadc/4mとしているので、符号68で示す周波数は0Hzとなる。したがって、ドップラー周波数測定のための周波数解析の対象となる帯域Δfp(71)は、最終的に、0Hzからfbsまでの周波数範囲に移ったことになる。
【0084】
ところで、図9及び図10において、帯域Δfp(71)の最高周波数は符号69で表されるfbsである。デジタルBPF62によるデジタルフィルタリングを行っているので、ダウンサンプリング部63からの出力には、周波数fbsを超える成分は現れない。したがって、fbsの2倍である周波数fadc/2m(67)をフーリエ変換のサンプリング周波数fftとすることができる。
【0085】
【数9】
【0086】
次に、信号切り出しゲート11からFFT処理部13までの動作を説明する。
【0087】
ダウンサンプリング部63からは、周波数fmf(70)を中心周波数とするように周波数が移された受信信号における瞬時値を表すデジタル値からなるデジタルデータ列が、クロック周波数をfadc/mとして逐次出力されている。そこで信号切り出しゲート11は、所望する時間範囲で、例えば、海中における深さ方向での潮流の流速分布を求める場合であれば、所望の深さに対応した時間範囲で、ダウンサンプリング部63から出力されるデジタルデータ列を切り出す。この切り出されたデータのサンプル数を切り出しサンプル数とする。ゼロ付加部12は、要求されている周波数分解能に対応するFFTのサンプル数を基準サンプル数として求め、基準サンプル数に達するように、基準サンプル数と切り出しサンプル数との差に相当する数のゼロデータ(値がゼロであるサンプル)を切り出されたデジタルデータ列に付加する。そして、ゼロデータが付加されてサンプル数が基準サンプル数とされたデジタルデータ列に対し、FFT処理部13がFFTを実行し、ピークパワースペクトルを算出する。算出されたピークパワースペクトルは、信号切り出しゲート11において切り出された範囲に対するものであって、信号におけるその時間範囲でのピークパワー周波数を示している。また、ゼロデータが付加されたことにより、信号においてFFTの対象となる時間幅が実質的に拡大されており、この時間幅の逆数が周波数分解能に該当することから、周波数分解能が向上している。このように、本実施形態によれば、受信信号の中での所望の時間範囲(あるいは距離範囲)のピークパワースペクトルを、FFTの演算だけで、所望の周波数精度を満たして検出することができる。
【0088】
以下、数値例を挙げて、この実施形態における周波数分解能及びFFTにおけるサンプリング周波数について説明する。
【0089】
背景技術の欄で述べたものと同じ条件で、超音波のパルス状の信号を海中に送信し、移動物からの反射エコーのドップラー周波数を測定する場合を考える。すなわち、海水中の音波伝搬速度CをC=1500(m/s)とし、送信周波数ftxをftx=120kHzとし、検出最大速度(水平方向)VをV=15(m/s)とし、音波を水平方向から斜め方向(θ=60°)に送受信するものとする。
【0090】
このとき、ドップラー信号は120±1.2(kHz)の範囲であり、速度精度として0.1m/sを要求すれば周波数分解能f0は12Hzとなる。また、ドップラー計測において検出対象とする位置の範囲(位置分解能)を7.5mとする。
【0091】
観測対象とすべき周波数帯域Δfpは、Δfp≧2×1200=2400としなければならない。また、フーリエ変換での最大周波数fbsをfbs=Δfpとする。
【0092】
信号の減衰帯域(1300Hz)を考慮して、
Δfp=(1200+1300)×2=2500×2=5000
とする。1300Hzの減衰帯域の外側では信号はほとんどゼロと考えてよいので、5kHz(すなわちfbs)においては信号はほとんどゼロである。
【0093】
式3−1及び式3−3より、
中間周波数の中心: fmid=3.5×fbs=17500
局部発振周波数: floc=ftx+fmid=120000+17500=137500
となる。また、オーバーサンプリングの係数mをm=4とすると、A/D変換のサンプリング周波数fadcは、式3−5より、
A/D変換: fadc=4×m×fbs=16×5000=80000
フーリエ変換のサンプリング周波数fftは、
【0094】
【数10】
【0095】
となり、周波数分解能f0を満たすフーリエ変換でのサンプル数Nは、N=fft/f0≒833.3となる。834以上であって834に最も近い2のべき乗数は210=1024であるので、FFTでのサンプル数Nを1024とする。サンプル数を1024としたことにより、10000/1024≒9.77Hzの周波数分解能が達成される。
【0096】
海中に音波パルスを発射して受信信号を観測し、その受信信号を中間周波数信号に変換した後の波形を示す図6において、上述したように、符号16は反射エコーを表し、符号17は送信信号自体を表し、符号18は海底エコーを表している。海底エコー18の周波数を解析し、式1−1に基づき、ドップラー周波数から船の速度を求める場合を考える。ここでは、送信信号などの不要な信号が含まれないようにするために、海底エコー18の部分だけを切り出してFFTを行い、海底からの反射波に関してドップラー周波数を求める。
【0097】
図13は、図6に示す信号波形から海底エコーだけを切り出した信号19を示している。ただし、図6はFFTでのサンプリング周波数を160kHzとしたものであるのに対し、本実施形態ではFFTでのサンプリング周波数を10kHzとしているので、図13では、時間長では図6と図13とが一致するように、サンプリング周波数の比に合わせて、横軸もサンプル数表示では1/16にしている。ここでは、海底エコー19の部分として、音波パルスの送信から時間t2(26)の経過後に時間幅T2(24)を切り出すものとする。ここで切り出された信号は、サンプル数としては50個である。
【0098】
時間幅T2の50個のサンプルのみを用いてフーリエ変換を行った場合の周波数分解能f0は、f0=10000Hz/50=200Hzとなり、要求されている周波数分解能である12Hzを全く満たさない。
【0099】
そこでこの実施形態では、FFTについては周波数分解能f0を満たすようにサンプル数N=1024で実行するものとし、海底エコーの信号19の分のサンプル数は50なので、不足した974個のサンプルについては、ゼロ値からなるサンプルを挿入する。図13においてT0(23)は、1024個のサンプルに対応する時間であり、FFT解析の対象となる時間幅である。
【0100】
ここでT0=T2+T3の関係が成り立っており、測定したい範囲の信号19を切り出してFFT解析に必要なサンプル数Nを満たすように付加されるゼロ値信号の継続時間が、T3に対応することになる。時間幅T3の挿入位置は、時間幅T2すなわち切り出された信号19の前でも後でもよい。あるいは、信号19を挟むように前後に分割してゼロ値信号を挿入する時間幅を配置してFFT解析の時間幅T0(27)を設定してもよい。
【0101】
ここで示す例では、音波パルスの送信から信号の切り出しを開始するまでの経過時間t2(図13での符号26)と時間幅T2(図13での符号24)を予め定めているが、海底エコーのような強度の大きな反射エコーを検出対象とする場合であれば、切り出しまでの経過時間t2も切り出しの時間幅T2も自動検出して決定することも可能である。深度ごとの潮流の流速を求める場合などには、自動検出を行うことができないのでt2及びT2について予め設定しておくことになる。
【0102】
ここでゼロ付加部12によるゼロ値データの付加について説明する。
【0103】
一般にFFTではその対象とするサンプル数は2のべき乗である必要がある。その一方で、FFTによる解析を行いたい対象データにおけるサンプル数が2のべき乗となっていないことはよくあることである。そこで、サンプル数が2のべき乗となるように、ゼロデータのサンプルを付加することがゼロパディングなどと称して行われている。その場合、とにかくFFTを行えるだけのサンプル数となればよいのであるから、対象データのサンプル数をNsとして、2L-1<Ns≦2Lを満たす自然数Lが存在する時に、2L−Ns個のゼロ値データを付加することになる。それに対して本実施形態では、単にFFTの実行が可能になるようにゼロ値データを付加するのではなく、FFTにおけるサンプリング周波数fftを考慮して所望の周波数分解能が得られるように、ゼロ値データを付加する。したがって、付加されるゼロ値データの数も、単にFFTを実行可能になるようにする場合に比べ、はるかに多くなることがある。上述した例では、海底エコーを切り出してサンプル数が50になったとして、単にFFTを実行するためだけであれば26=64であるので14個のゼロ値データを付加すればよいだけであるが、所望の周波数分解能を満たすようにするために、974個もののゼロ値データを付加しているのである。
【0104】
またこの実施形態では、FFTにおけるサンプル数は1024であるので、複素乗算の演算回数は10240回で済む。図2〜図5に示した従来の例では、本実施形態におけるものと同じ周波数分解能を得るために、FFTにおいて229376回の複素乗算を行わなければならないので、本実施形態によれば、同じ周波数分解能を得るために必要な複素乗算の回数が大幅に減少することがわかる。また、図2〜図5に示した従来の例では、海底エコーの部分だけを切り出してFFTを行ったときの周波数分解能は156Hzであったのに対し、本実施形態では9.77Hzの周波数分解能が維持される。
【0105】
このように本実施形態によれば、周波数分解能と位置分解能とを両立させることができ、しかも、フーリエ変換で必要なサンプル数を大幅に減らして演算量を大幅に削減することが可能になる。
【0106】
以上、本発明の好ましい実施形態について、反射エコーにおけるドップラー周波数の検出を例に挙げて説明したが、本発明の適用範囲はこれに限られるものではない。検出対象とすべき周波数範囲を想定できる信号から短い時間範囲の部分を切り出してFFTによる解析を行う場合に、本発明は広く一般的に適用可能である。
【符号の説明】
【0107】
1 送受波器
2 送受信切替え回路
3 送信回路
4 受信増幅器
5 変調器
6 局部発振回路
7 アナログフィルタ
10,61 A/D変換器
11 信号切り出しゲート
12 ゼロ付加部
13 FFT処理部
60 高次アナログフィルタ
62 デジタルBPF(バンドパスフィルタ)
63 ダウンサンプリング部
【特許請求の範囲】
【請求項1】
時間変化する受信信号における第1の周波数を含む所定の周波数帯域幅内で、かつ所望の時間幅の範囲内で、前記受信信号のピークパワー周波数を検出する方法であって、
前記第1の周波数とは異なる中間周波数の信号に前記受信信号を周波数変換する段階と、
アナログフィルタを適用して、前記周波数変換された信号から高域成分を除去する段階と、
前記アナログフィルタを適用した後の信号に対して第2の周波数でサンプリングしてA/D変換する段階と、
A/D変換で得られたデジタル信号に対してデジタルバンドパスフィルタを適用し、前記所定の周波数帯域幅に相当する信号のみを抽出する段階と、
前記第2の周波数の2のべき乗分の1の周波数である第3の周波数をサンプリング周波数として、前記所望の時間幅内において、前記デジタルバンドパスフィルタから出力される信号をダウンサンプリングして抽出するダウンサンプリング段階と、
必要とされる周波数分解能を満たすために高速フーリエ変換において必要なサンプル数を基準サンプル数として、該基準サンプル数に達するように、前記ダウンサンプリング段階から出力されるデジタルデータ列にゼロデータを付加する付加段階と、
前記ゼロデータが付加された前記デジタルデータ列に対して高速フーリエ変換を行うFFT段階と、
を有し、
前記ダウンサンプリング段階において、前記ダウンサンプリングによって抽出される信号の周波数帯域が周波数ゼロから前記所定の周波数帯域幅に相当する帯域幅内に配置され、
前記高速フーリエ変換の結果において最大値を示す周波数を前記ピークパワー周波数として検出する方法。
【請求項2】
前記所定の周波数帯域幅の周波数幅をΔfとし、第2及び第3の周波数をそれぞれf2及びf3とし、前記中間周波数をfmidとし、nを1以上の整数としてm=2nとし、
fmid=3.5×Δf,
f2=4×m×Δf,
f3=f2/(2×m)
とする、請求項1に記載の方法。
【請求項3】
前記基準サンプル数は、前記サンプリング段階で得られる前記デジタルデータ列のサンプル数の2倍以上である、請求項1または2に記載の方法。
【請求項4】
時間変化する受信信号における第1の周波数を含む所定の周波数帯域幅内で、かつ所望の時間幅の範囲内で、前記受信信号のピークパワー周波数を検出する装置であって、
前記受信信号を前記第1の周波数とは異なる中間周波数の信号に周波数変換する周波数変換手段と、
前記周波数変換された信号から高域成分を除去するアナログフィルタと、
前記アナログフィルタから出力される信号を第2の周波数でサンプリングしてA/D変換するA/D変換器と、
A/D変換器の出力に接続され、前記所定の周波数帯域幅に相当する信号のみを含むデジタル信号を出力するデジタルバンドパスフィルタと、
前記第2の周波数の2のべき乗分の1の周波数である第3の周波数をサンプリング周波数として、前記所望の時間幅内において、前記デジタルバンドパスフィルタから出力される信号をダウンサンプリングして抽出するダウンサンプリング手段と、
必要とされる周波数分解能を満たすために高速フーリエ変換において必要なサンプル数を基準サンプル数として、該基準サンプル数に達するように、前記ダウンサンプリング手段から出力されるデジタルデータ列にゼロデータを付加するゼロ付加手段と、
前記ゼロデータが付加された前記デジタルデータ列に対して高速フーリエ変換を行うFFT手段と、
を有し、
前記ダウンサンプリング手段から出力される前記デジタルデータ列における信号の周波数帯域が周波数ゼロから前記所定の周波数帯域幅に相当する帯域幅内に配置され、
前記高速フーリエ変換の結果において最大値を示す周波数を前記ピークパワー周波数として検出する、装置。
【請求項5】
前記所定の周波数帯域幅の周波数幅をΔfとし、第2及び第3の周波数をそれぞれf2及びf3とし、前記中間周波数をfmidとして、nを1以上の整数としてm=2nとし、
fmid=3.5×Δf,
f2=4×m×Δf,
f3=f2/(2×m)
とする、請求項4に記載の装置。
【請求項6】
前記基準サンプル数は、前記サンプリング手段で得られる前記デジタルデータ列のサンプル数の2倍以上である、請求項5に記載の装置。
【請求項1】
時間変化する受信信号における第1の周波数を含む所定の周波数帯域幅内で、かつ所望の時間幅の範囲内で、前記受信信号のピークパワー周波数を検出する方法であって、
前記第1の周波数とは異なる中間周波数の信号に前記受信信号を周波数変換する段階と、
アナログフィルタを適用して、前記周波数変換された信号から高域成分を除去する段階と、
前記アナログフィルタを適用した後の信号に対して第2の周波数でサンプリングしてA/D変換する段階と、
A/D変換で得られたデジタル信号に対してデジタルバンドパスフィルタを適用し、前記所定の周波数帯域幅に相当する信号のみを抽出する段階と、
前記第2の周波数の2のべき乗分の1の周波数である第3の周波数をサンプリング周波数として、前記所望の時間幅内において、前記デジタルバンドパスフィルタから出力される信号をダウンサンプリングして抽出するダウンサンプリング段階と、
必要とされる周波数分解能を満たすために高速フーリエ変換において必要なサンプル数を基準サンプル数として、該基準サンプル数に達するように、前記ダウンサンプリング段階から出力されるデジタルデータ列にゼロデータを付加する付加段階と、
前記ゼロデータが付加された前記デジタルデータ列に対して高速フーリエ変換を行うFFT段階と、
を有し、
前記ダウンサンプリング段階において、前記ダウンサンプリングによって抽出される信号の周波数帯域が周波数ゼロから前記所定の周波数帯域幅に相当する帯域幅内に配置され、
前記高速フーリエ変換の結果において最大値を示す周波数を前記ピークパワー周波数として検出する方法。
【請求項2】
前記所定の周波数帯域幅の周波数幅をΔfとし、第2及び第3の周波数をそれぞれf2及びf3とし、前記中間周波数をfmidとし、nを1以上の整数としてm=2nとし、
fmid=3.5×Δf,
f2=4×m×Δf,
f3=f2/(2×m)
とする、請求項1に記載の方法。
【請求項3】
前記基準サンプル数は、前記サンプリング段階で得られる前記デジタルデータ列のサンプル数の2倍以上である、請求項1または2に記載の方法。
【請求項4】
時間変化する受信信号における第1の周波数を含む所定の周波数帯域幅内で、かつ所望の時間幅の範囲内で、前記受信信号のピークパワー周波数を検出する装置であって、
前記受信信号を前記第1の周波数とは異なる中間周波数の信号に周波数変換する周波数変換手段と、
前記周波数変換された信号から高域成分を除去するアナログフィルタと、
前記アナログフィルタから出力される信号を第2の周波数でサンプリングしてA/D変換するA/D変換器と、
A/D変換器の出力に接続され、前記所定の周波数帯域幅に相当する信号のみを含むデジタル信号を出力するデジタルバンドパスフィルタと、
前記第2の周波数の2のべき乗分の1の周波数である第3の周波数をサンプリング周波数として、前記所望の時間幅内において、前記デジタルバンドパスフィルタから出力される信号をダウンサンプリングして抽出するダウンサンプリング手段と、
必要とされる周波数分解能を満たすために高速フーリエ変換において必要なサンプル数を基準サンプル数として、該基準サンプル数に達するように、前記ダウンサンプリング手段から出力されるデジタルデータ列にゼロデータを付加するゼロ付加手段と、
前記ゼロデータが付加された前記デジタルデータ列に対して高速フーリエ変換を行うFFT手段と、
を有し、
前記ダウンサンプリング手段から出力される前記デジタルデータ列における信号の周波数帯域が周波数ゼロから前記所定の周波数帯域幅に相当する帯域幅内に配置され、
前記高速フーリエ変換の結果において最大値を示す周波数を前記ピークパワー周波数として検出する、装置。
【請求項5】
前記所定の周波数帯域幅の周波数幅をΔfとし、第2及び第3の周波数をそれぞれf2及びf3とし、前記中間周波数をfmidとして、nを1以上の整数としてm=2nとし、
fmid=3.5×Δf,
f2=4×m×Δf,
f3=f2/(2×m)
とする、請求項4に記載の装置。
【請求項6】
前記基準サンプル数は、前記サンプリング手段で得られる前記デジタルデータ列のサンプル数の2倍以上である、請求項5に記載の装置。
【図1】
【図2】
【図3】
【図4】
【図5】
【図6】
【図7】
【図8】
【図9】
【図10】
【図11】
【図12】
【図13】
【図2】
【図3】
【図4】
【図5】
【図6】
【図7】
【図8】
【図9】
【図10】
【図11】
【図12】
【図13】
【公開番号】特開2012−247302(P2012−247302A)
【公開日】平成24年12月13日(2012.12.13)
【国際特許分類】
【出願番号】特願2011−119217(P2011−119217)
【出願日】平成23年5月27日(2011.5.27)
【出願人】(303057044)株式会社ソニック (17)
【Fターム(参考)】
【公開日】平成24年12月13日(2012.12.13)
【国際特許分類】
【出願日】平成23年5月27日(2011.5.27)
【出願人】(303057044)株式会社ソニック (17)
【Fターム(参考)】
[ Back to top ]