JP4077092B2 - ドップラ周波数測定方法およびドップラソナー - Google Patents
ドップラ周波数測定方法およびドップラソナー Download PDFInfo
- Publication number
- JP4077092B2 JP4077092B2 JP33713998A JP33713998A JP4077092B2 JP 4077092 B2 JP4077092 B2 JP 4077092B2 JP 33713998 A JP33713998 A JP 33713998A JP 33713998 A JP33713998 A JP 33713998A JP 4077092 B2 JP4077092 B2 JP 4077092B2
- Authority
- JP
- Japan
- Prior art keywords
- frequency
- doppler
- doppler frequency
- prf
- signal
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Lifetime
Links
Images
Description
【発明の属する技術分野】
この発明は、パルス信号を用いて高精度のドップラ周波数測定を可能にしたドップラ周波数測定方法およびドップラソナーに関する。
【0002】
【発明が解決しようとする課題】
船舶の種類によっては非常に高精度の速度計測が必要なものがあり、また、通常の船舶でも港内とりわけ接岸時には極めて低速の速度制御が必要となる。しかし、通常のパルスドップラ装置では、単発のパルスによる計測であるため1ノット程度の速度分解能しか実現できない。これは、速度分解能(周波数分解能)はパルス幅に比例するが、超音波パルスビームは、パルスの先頭が対象物に反射してそのエコー信号が返ってくるまでの時間以下のパルス幅でなければビームを送信することができないからである。
【0003】
したがって、対象物に接近するほど低速航行になり、より高い速度分解能が必要になるが、このような場面ほどエコーの遅れ時間が短くなるため送信できるパルス幅が短くなり速度分解能が低下するという問題点があった。
【0004】
また、複数回の測定結果を平均して測定精度を上げる技術も実用化されているが、この方式で1ヘルツ以下の周波数分解能を実現するためには、極めて長時間の積算平均処理が必要であり、上記岸壁に接近している時などに適用できるものではなかった。
【0005】
一方、送信用トランスデューサと受信用トランスデューサを別々に設け、送信用トランスデューサは専ら連続波の超音波ビームを送信し、受信用トランスデューサでは専らそのエコー信号を受信するようにして、長い時間幅のFFTを可能にすることにより高い周波数分解能を実現することも考えられるが、連続波ビームを送信するためには送信用トランスデューサに常時駆動電力を印加する必要があるため、トランスデューサの発熱や駆動回路の負荷など負担が大きく、パルス波に比べてビーム出力を抑えざるを得ない。このため、連続波ビームは、対象物に対して50メートル程度以内に接近してからでなければ使用することができないという問題点があった。
【0006】
この発明は、パルスビーム方式で高精度にドップラ周波数を測定することができるドップラ周波数測定方法およびドップラソナーを提供することを目的とする。
【0007】
【課題を解決するための手段】
請求項1の発明は、一定間隔で送信された複数のパルス信号のエコー信号を受信してサンプリングデータとし、このサンプリングデータを用いて離散フーリエ変換を行ってベースバンドのスペクトルのピーク周波数を求める手順と、
パルス信号のエコー信号に基づいて概略のドップラ周波数を求め、この概略のドップラ周波数に基づいて前記ベースバンドのスペクトルが真のスペクトルであるかまたは写像であるかを求め、さらに、写像である場合には真のスペクトルとの周波数差を求める手順と、
前記ピーク周波数を前記周波数差で補正することによって精ドップラ周波数を求める手順と、
を有することを特徴とする。
【0008】
請求項2の発明は、
(1) 周期信号であるパルス信号を送信して対象物で反射したエコー信号を受信し、該エコー信号をフーリエ変換して粗ドップラ周波数を求める
(2) 周期信号であるパルス信号を送信して対象物で反射したエコー信号を受信し、該エコー信号から前記パルス信号の周波数成分を除去してドップラ周波数の波形成分のみをサンプリング周波数PRFの複素数データ(以下「PRFサンプリングデータ」という)として取り出す
(3) 手順(2)をパルス送信周波数PRFで繰り返し実行する
(4) 複数のPRFサンプリングデータを用いてフーリエ変換し、ベースバンド(0Hz〜PRF)におけるスペクトルのピーク周波数fdを求める
(5) 手順(1)で求めた粗ドップラ周波数に基づき、ドップラ周波数のスペクトルの属する第n次帯域(PRF×n〜PRF×(n+1)…n:整数)を割り出す
(6) 手順(4)で求めたfd、および、手順(5)で求めたnを用い、fd+PRF×nの演算によって精ドップラ周波数を求める
の手順を有することを特徴とする。
【0009】
請求項3の発明は、請求項2の発明において、手順(1)を複数回繰り返し、求められた複数の粗ドップラ周波数を平均した値を手順(5)の粗ドップラ周波数として用いることを特徴とする。
【0010】
請求項4の発明は、請求項2,3の発明において、手順(1)および手順(2)のうち、周期信号であるパルス信号を送信して対象物で反射したエコー信号を受信する手順を共通にしたことを特徴とする。
【0011】
請求項5の発明は、請求項2〜4の発明において、サンプリング周波数fsがエコー信号の注目領域の中心周波数fcに対してfs=4fc/(4k+1)またはfs=4fc/(4k+3)、(k=0または正の整数)となるようなサンプリング周波数fsを用いて該エコー信号をサンプリングし、このサンプリングデータ列に対して+1,−j,−1,+jまたは+1,+j,−1,−jを順次乗算することによって、手順(2)のエコー信号から前記パルス信号の周波数成分を除去する手順を処理することを特徴とする。
【0015】
≪発明の概要≫
図面を参照してこの発明について説明する。
パルス信号を送信して、そのエコー信号を受信する。一般のパルスドップラ装置の場合このエコー信号をFFTしてピーク周波数を検出し、このピーク周波数と送信した超音波ビームの周波数との周波数差をドップラ周波数として割り出す。この方式で検出したドップラ周波数の精度(周波数分解能)は、数十ヘルツ〜数百ヘルツであり、これを複数回繰り返して平均化することにより精度を1桁程度あげることができる。この従来の単一パルスモードによるドップラ周波数の測定結果をこの発明では概略のドップラ周波数である「粗ドップラ周波数」として用いる。なお、上記処理は、エコー信号を中間周波数にダウンコンバートした信号で行われる場合が多い。
【0016】
そして、この発明では、エコー信号から元のパルス信号の周波数成分を除去する(該周波数だけ周波数シフトする)。パルス状のエコー信号は、極めて低周波で波長の長いドップラ周波数の波形に対しては、図1に示すようにその波形上の1点のサンプリングデータ程度の時間幅である。そこで、これをPRFサンプリングデータとして用いる。実際には、パルス幅の中央値を用いたり、一定区間を平均化するなどして1つのサンプリングデータを求める。
【0017】
パルス信号を一定周波数PRF(pulse repetition frequency)で繰り返し送信して複数のPRFサンプリングデータを獲得し、これを用いて離散フーリエ解析を行う。たとえば、送信周期306ms(PRF=3.27Hz)のエコー信号のPRFサンプリングデータを32点求めてフーリエ解析を行った場合には、ほぼ0.1Hzの周波数分解能を得ることができる。このようにパルスドップラでありながら、極めて高い周波数分解能を得ることができる(図1(A)参照)。
【0018】
ただし、離散フーリエ解析のサンプリング周波数はパルスの繰り返し周波数PRFであり、これは数ヘルツ程度の低いものである(上記の例では3.27Hzである)。このため、真のドップラ周波数がこれよりも高い場合があり、この場合、0ヘルツ(DC)〜PRFのベースバンド(n=0)に現れるスペクトルは、より高い帯域からの写像である(図1(B)参照)。そこで、前記粗ドップラ周波数を用いて真のドップラ周波数がどの帯域(PRF×n〜PRF×(n+1)…n:整数)からの写像であるかを割り出す。図2に信号の帯域Bとサンプリング周波数fsに応じて発生する信号帯域の写像パターンを示す。上記粗ドップラ周波数によって求められた写像帯域の基底周波数n×PRFをベースバンドスペクトルのピーク周波数fdに加算することによって精ドップラ周波数を求める。なお、請求項1,2のピーク周波数は、単純なピークであってもよく、ピーク領域を加重平均または積分してピークの中心値を求めてもよい。
【0019】
【発明の実施の形態】
図面を参照してこの発明の実施形態であるドップラソナーについて説明する。図3は同ドップラソナーの構成図、図4は同ドップラソナーの受信部のブロック図である。このドップラソナーは、ドップラ周波数の測定モードとして、単一パルスモード、連続パルスモードの2種類のモードを備えている。
【0020】
図3において、ドップラソナーは、トランスデューサ1、切換器2、送信部3,受信部4および制御部5からなっている。送信部3は超音波ビームを形成するための超音波信号を発生する駆動回路を有し、トランスデューサ1から超音波ビームを送信する回路である。受信部4は、前記トランスデューサ1が受信したエコー信号を受信してそのドップラ周波数を測定する。切換器2はトランスデューサ1の接続を送信部3側または受信部4側に切り換える回路である。制御部5は、切換器2、送信部3、受信部4を制御する。
【0021】
トランスデューサ1から送信された超音波ビームは対象物で反射し、そのエコーが再度トランスデューサ1に受信される。装置(このドップラソナーを搭載した移動体)が対象物に対して相対的に接近していると(または遠ざかっていると)、エコーの周波数は、送信周波数に対して上に(または下に)ドップラシフトしている。このドップラ効果によるシフト量であるドップラ周波数が両者の相対速度に比例する。水中では音波の伝搬速度は約1500m/秒であるため、送信周波数が400kHzであれば、1Hzのドップラシフトが約0.01ノット(約0.005m/秒)の相対速度に対応する。
【0022】
また、パルスビームを送信したのちそのエコーを受信するまでの時間差が装置と対象物との距離に比例する。水中であれば、約0.133秒の時間遅れが100mの距離に対応する。
【0023】
図4の受信回路において、トランスデューサ1は、送信した超音波ビームのエコー信号を受信して電気信号に変換する。このエコー信号は、送信周波数とほぼ同じ周波数(たとえば438kHz)であるが、対象物との相対速度に応じて0〜±1kHz程度のドップラシフトを受けている。このエコー信号は、高周波アンプ11に入力される。高周波アンプ11はこの信号を後段の回路で処理可能なレベルまで増幅してバンドパスフィルタ13は不要周波数成分を除去し、帯域制限を加える。
【0024】
アンプ14はこの中間周波に変換されたエコー信号を適当なレベルまで増幅したのちA/D変換器15に入力する。A/D変換器15は、この信号を24kHzのサンプリング周波数でサンプリング(アンダーサンプリング)し、量子化・符号化することによってディジタルデータに変換する。438kHzの信号を24kHzでアンダーサンプリングすると、438kHz±6kHzの周波数領域のスペクトルが6kHz±6kHz(0(DC)〜12kHz)および18kHz±6kHz(12kHz〜24kHz)に現れる。したがって、438kHz付近に生じていたエコー信号のスペクトルは、6kHz付近および18kHz付近に現れる。
【0025】
なお、このように438kHzの信号を24kHzでサンプリングすることにより、15ノット程度までの相対速度を測定するのであれば、必要な周波数帯域幅は3kHz程度であるため、上記±6kHzの周波数領域がサンプリングされることで十分なレンジをカバーすることができる。なお、エリアシングを生じさせないため、前記バンドパスフィルタによって帯域を十分に絞り込んでおけばよい。
【0026】
A/D変換器15のアンダーサンプリングによって、ディジタル変換されたエコー信号はDSP16に入力される。超音波パルスビームの送信・エコーの受信は一定間隔(PRF)で繰り返し行われるため、DSP16へのディジタルエコー信号の入力も定期的に行われる。
【0027】
DSP16は、入力されたエコー信号を用いて単一パルスモードのドップラ周波数測定を行い、これを粗ドップラ周波数として、正確なドップラ周波数である精ドップラ周波数の属する帯域の特定に用いる。そして、直前までのパルスモード送受信データを用いて精ドップラ周波数を割り出す。
【0028】
以下、図5のフローチャートおよび図6のベースバンドおよびn次写像バンドの図を参照してDSP16の処理動作を説明する
図5において、A/D変換器15からエコー信号が入力されると(s1)、この反射パルス信号を用いて従来と同様の単一パルスによるドップラ周波数の計測を行う(s2)。すなわち、この24kHzでサンプリングされた信号でFFTを実行し、周波数分解能数十Hz程度の粗ドップラ周波数を割り出し、これを記憶する。
【0029】
つぎにこのエコー信号を6kHz下方に周波数シフトする。ただし、サンプリングにおいてスペクトル反転が発生する場合には、18kHz下方に周波数シフト、または6kHz上方にシフトする。すなわち、超音波パルスの送信周波数成分を除去し、ドップラ周波数成分のみが現れるようにデモデュレーションする(s4)。そして、ベースバンドのスペクトル以外を除去するようにフィルタリングする(s5)ことで、データを複素数化する(s6)。デモジュレーションされたエコー信号は図1に示したようにドップラ周波数波形の一部であるため、この信号の中央部の一部区間を平均するなどしてPRFサンプリングデータを割り出す(s7)。s1で受け取るエコー信号はパルス繰り返し周波数PRF毎に入力されるため、s1〜s7の処理はエコー信号が入力される毎に繰り返し実行される。
【0030】
所定数(たとえば32)のエコー信号が入力されるごとに(s8)、s9以下の精ドップラ周波数の算出動作に進む。または、FFTに必要な点数の一部を更新しながら精ドップラ周波数の算出動作に進む。
【0031】
s9以下では以下の動作を実行する。まず精度がPRF以下になるように粗ドップラ周波数を平均化する。すなわち、1回の粗ドップラ周波数の測定では真のドップラ周波数を含む帯域(周波数PRFのn倍の周波数にある周波数帯域幅の周波数領域)を動揺などの外乱雑音によって絞り込むことができないため、精度がPRF以下になるように蓄積してきた粗ドップラ周波数のうち最新のものから必要な回数分を平均する。そして蓄積したPRFサンプリングデータのうち最新の所定点数(たとえば128)を用いてFFTを実行し(s10)、ベースバンドに変換してピーク周波数fdを検出する(s11)。このピーク周波数fdは、単純なピークであってもよく、ピーク領域を加重平均または積分してピークの中心値を求めてもよい。ただし、実際のドップラ周波数はベースバンドよりも高い周波数または低い周波数である可能性があるため、前記s9で算出された平均の粗ドップラ周波数に基づいて真のドップラ周波数がどの写像帯域(第n写像帯域:n×PRF〜(n+1)×PRF)にあるかを割り出す(s12)。図6の場合にはn=5の写像帯域を示しているとする。そして、前記ピーク周波数fdにn×PRFを加算して精ドップラ周波数を算出する(s13)。
【0032】
なお、上記では粗ドップラ周波数の検出において平均処理を用いた例を示したが、粗ドップラ周波数の検出において、FFTによって求めた周波数範囲をさらに狭い範囲に納めてより周波数精度をあげるために、MEM(最大エントロピー法やAR法(自己回帰法)などのスペクトル推定手法をFFT結果に組み合わせる。これによって、PRFが低いときに精度の高い周波数範囲を決定することができる。
【0033】
MEM等のスペクトル推定法では、少ないデータ点数においても、鋭いピーク周波数を検出することができるという特徴を有するが、正しいピーク以外に誤ピークも多く発生させるため、スペクトル推定法で検出された複数のピークのうち、FFTによって求められたブロードなピーク周波数の近傍のものを選択することによって正確に鋭いピークを選択することが可能になる。
【0034】
ただし、粗ドップラ周波数を用いてスペクトル折返し写像の帯域番号nを求める場合、ドップラ周波数(ピーク周波数)fdがPRFの整数倍の場合には、この周波数がどの写像に属するものかは上記方法では求めることができなくなる。特にfd=PRFの場合にはドップラ周波数成分が離散周波数0の上に折り返されるためDC成分と区別がつかなくなる。このため、このような場合には、再サンプリングした元のデータを直接調べるなどの別の方法で求める。
【0035】
なお、PRFは常に一定ではなく、対象物との距離に応じて決定される。すなわち、対象物との距離が長いとエコーの遅延時間が長いためPRFは小さくなり(その変わりパルス幅が長くなる)、対象物との距離が短いとエコーの遅延時間が短いためPRFを大きくすることができる。
【0036】
このように単一のトランスデューサ1を用いたドップラソナーで1Hz以下の周波数分解能を得ることができるため、従来より船舶等に搭載されているドップラソナーの処理アルゴリズムを変更するのみで計測精度を著しく向上させることも可能である。たとえば、岸壁に接岸する場合などは、微速航行となりドップラ周波数も小さくなり、岸壁との距離も短くなるためPRFを10Hz(1秒間に10回パルスビームを送受信)程度にあげることができるため、極めて有効にこの方式を適用することができる。
【0037】
図7〜図10は、上記方式でドップラ周波数を算出するシミュレーション結果を示す図である。
図7,図8は、0〜PRF(Hz)の基本範囲(n=0)に真のドップラ周波数が存在する場合(図1(A)参照)をシミュレートしたものである。このシミュレーションは基本周波数438kHz、ドップラ周波数=1.2Hzのパルス信号をPRF=5.8594Hzで繰り返し受信し、24kHzでサプリングしたとして実行した。図7(A)は連続した送信パルスに対するエコー信号を示す図、図7(B)は同図(A)に示したもののうち先頭のエコー信号を時間軸方向に拡大した図である。また、図7(C)は同図(B)の波形をFFT解析した結果を示すスペクトル波形であり、6kHz付近と18kHz付近にエコー信号のスペクトルのピークが存在している。図7(D)は同図(C)のエコー信号をデモジュレーションした結果であり、6kHz付近にあったスペクトルのピークが0Hz付近に移動している。こののち、12kHz近傍の周波数成分をフィルタリングによって除去し、複素数の信号を形成する。
【0038】
図8(A)は、ローパスのフィルタリングの後に出力される複素ドップラ信号の実数部を示している。このシミュレーションでは図8(A)の信号の前後に現れる不連続点におけるリンギングを避けるため、同図(B)のようにその中央の一部分のみを用いて平均し、1つの受信波信号に対するサンプル値を求めている。図8(C)は図7(A)に示す複数の受信波信号に対する実数部のサンプル値を時系列に表示したものである。すなわち、復調されたエコー信号をサンプリング周波数PRFによって再サンプリングした信号波形を表している。そして、図8(D)が同図(C)の信号をFFT解析した結果を示し、このパワースペクトルのピーク周波数から求めたドップラ周波数は、1.2016Hzであり、0.016Hzの誤差でドップラ周波数を求めることができた。なお、図8(C)のサンプル点数は8点であるがゼロデータを末尾に追加して(補間はしていない)1024点でFFTを行った。
【0039】
図9は、真のドップラ周波数がPRFよりも高い場合(図1(B)参照)をシミュレートしたものである。このシミュレーションは基本周波数438kHz、ドップラ周波数=15.0Hzのパルス信号をPRF=2.9297Hzで繰り返し受信し、24kHzでサプリングしたとして実行した。同図(A)は受信した信号を直接FFTした結果を示す図である。複数のパワースペクトルから平均粗ドップラ周波数を計算すると14.65Hzが得られる。これによると、ドップラ周波数fdに対して、
|fd|≧n・PRF
を満足する最大整数nは5になる。ただし、より精度を上げるためには、粗ドップラ周波数の計算はデモジュレーション、フィルタリング後の複素ドップラ信号から求めることも可能である。同図(B)は図8(C)に対応する図であり、サンプリング周波数PRFによって再サンプリングしたドップラ信号の実数部を表している。同図(C)は(B)のFFT解析結果を示している。
【0040】
(C)のパワースペクトルのピーク周波数0.3519Hzにn・PRF、すなわち、n=5であるため、2.9297×5=14.6485Hzを加算すると、精ドップラ周波数15.0004Hzが求まる。(D)は(C)の信号に折り返し分に相当するゼロデータを挿入し折り返しスペクトルを生成している。右端のスペクトルピークが所望のドップラ周波数であり、左端のスペクトルが(C)に示している再サンプリング後に得られるスペクトルに相当している。
【0041】
以上は、ドップラ周波数がプラスの場合、すなわち、装置が対象物に接近している場合について検討したが、ここではドップラ周波数がマイナスの場合について検討する。すなわち、図10(A)に示すようにデモジュレーション・フィルタリング後の複素ドップラ信号周波数fdが、たとえば−2fsよりも低い場合を想定する。この信号をfs(=PRF)で再サンプリングすると、同図(B)の離散周波数スペクトルが得られる。このディジタル信号をFFTすると、(C)の離散的周波数スペクトルが得られる。この場合には|fd|≧n・PRFを満足する最大整数nは2になるが、マイナス周波数の場合には−(n+1)を用いるため、属する写像帯域は第−3帯域である。したがって、(C)の離散周波数kdから求まる連続時間での換算周波数に−3・PRFを加算する。
【0042】
図11は上記ドップラ周波数がマイナスの場合のシミュレーションを示す図である。このシミュレーションではドップラ周波数=−6.6Hz、PRF=2.9297Hzの場合を示している。同図(A)は1回の送信に対して受信した信号を直接FFTした結果を示す図である。同図(B)は図8(C)に対応する図であり、サンプリング周波数PRFによって再サンプリングしたドップラ信号の実数部を表している。同図(C)は(B)のFFT解析結果を示している。
【0043】
同図(A)のFFT結果に基づき、複数のパワースペクトルから平均粗ドップラ周波数を計算する。また、(C)のパワースペクトルから求まるピーク周波数は2.1887Hzであり、これに粗ドップラ周波数より求めたn=2に基づき、−(n+1)・PRFを加算する、すなわち、2.9297×(−3)=−8.7891Hzを加算すると、精ドップラ周波数−6.6004Hzが求まる。(D)はCの信号に折り返し分に相当するゼロデータを挿入し折り返しスペクトルを生成している。ただし、ドップラ周波数が負値であるため、(D)のスペクトルが全体的にn・FPF=3PRFだけ左にシフトしたスペクトルになっている。
【0044】
以上のシミュレーションで明らかなように、粗ドップラ周波数が得られる平均処理以降では精ドップラ周波数の計算に数秒間の測定データを用いることで可能である。つまり、数秒前からの測定データを基に計算するため、計測ドップラ周波数の時間分解能が飛躍的に向上することが判る。
【0045】
なお、図9(D)に示しているように、直接補間によっても所望ピーク周波数を決定できる。このため、補間法を用いる方法も可能である。また、このシミュレーションではカウンタ回路を想定しているため、サンプリング周波数24kHzの基準周波数から波形する周波数PRFを用いているが、PRFは任意である。
【0046】
なお、DSP16のデモジュレーション処理は、サンプリングデータ例に対して+1,−j,−1,+jを順次乗算するのみでよい。以下、この簡略化された演算処理方式について説明する。なお、この処理を可能にするために、A/D変換器15のサンプリング周波数を24kHzに設定し、438kHz付近にスペクトルが展開しているエコー信号をアンダーサンプリングした。
【0047】
デモジュレーションは、エコー信号のスペクトルが展開している注目領域の中心周波数をゼロ周波数(DC)にシフトする処理、すなわち、周波数軸に対してスペクトルを並行移動する処理である。ここで、A/D変換器15から入力されるエコー信号のサンプリングデータ列をx(n)と表す。このサンプリングデータ列x(n)は、上記のようにサンプリング周波数fs(=24kHz)でサンプリングされ、中心周波数fc(=6kHz)の離散時間信号となったものである。これに対して、以下のような離散複数指数関数列c(n)を乗算する。
【0048】
【数1】
このc(n)をx(n)の各項に乗算することによって中心周波数fcが0(DC)になるように周波数スペクトルをシフトすることができる。すなわち、データ列x(n)のDFT変換から求まる周波数スペクトルが、
【数2】
であるのに対し、データ数列x(n)に離散複素指数関数c(n)を乗算したデータ列の周波数スペクトルXshift(k)が、
【数3】
となることから、この乗算によりデータ列の周波数スペクトルX(k)が周波数軸に沿ってシフトされていることが分かる。すなわち、
【数4】
によってスペクトルの注目領域の中心周波数fcを周波数ゼロとするように、スペクトル全体を周波数軸に沿ってシフトすることができる。
【0049】
また、前記c(n)の指数部(−jΩc n)のnを、自然数Mを加算することによって(n+M)に置き換えた場合、すなわち,離散複素指数関数をM個シフトしてデータ数列に乗算した場合でも、
【数5】
で明らかなように、周波数パワースペクトルはこのずれに影響されることなく同様にシフトされる。
【0050】
ここで、サンプリング周波数fsと注目領域の中心周波数fcは、上述したように
fs=4fc/(4k+1) ,k=0,1,2,…
となるような関係に設定されている。このため、上記、Ωc は、
Ωc =2π(fc/fs)=2πk+π/2
となり、前記離散複素指数関数c(n)は、
【数6】
となる。したがって、任意の整数値nに対して、
【数7】
となり、+1,−j,−1,+jの4種類の値のみを取ることが分かる。
【0051】
したがって、x(n)に対するc(n)の乗算は、x(n)に対して+1,−j,−1,+jを順次乗算するのみでよい。また、この順序が維持されていれば、〔数5〕で証明されたように、x(n)の各データと+1,−j,−1,+jとの対応は任意である。
【0052】
上記の例では6kHz付近に展開しているエコー信号のスペクトルを用いているが、18kHz付近に展開しているスペクトルを用いることも可能である。この場合、注目領域の中心周波数fc=18kHzとなり、
fs=4fc/(4k+3) ,k=0,1,2,…
となるため、
Ωc =2π(fc/fs)=3π/2
となり、Ωc =π/2の場合と逆回りの+1,+j,−1,−jをx(n)に乗算することで18kHzを0Hzにシフトすることができる。
【0053】
しかし、実際に乗算を行う必要はなく、乗算した場合に合わせて正負符号制御および実数虚数制御をするだけでよい。すなわち、c(n)がマイナス符号の場合には符号反転計算のみを行い、c(n)が実数の場合はx(n)の値を全て実数部として処理し、c(n)が虚数の場合はx(n)の値を全て虚数部として処理すればよい。
【0054】
このように、
fs=4fc/(4k+1) ,k=0,1,2,…
となるようなサンプリング周波数fsでサンプリングすることにより、サンプリングデータのサンプリング番号に基づいて符号制御および実数部,虚数部に割り振るのみの処理で周波数スペクトルのシフトを行うことができ、上記指数関数を実際に乗算して演算する必要がなくなるため、処理を大幅に簡略化することができる。
【0055】
そして、図12のようなレジスタの転記を行うことによって、x(n)の周波数スペクトルをシフトすると同時に、実数値としてサンプリングされたデータ列x(n)を複素数に変換することができる。ここで、A/D変換器15から入力されたデータ列x(n)を記憶するA/Dデータバッファをx(n)で表し、複素数バッファをX(n)で表す。
【0056】
この図において、x(0)はそのままX(0)の実数部に転記され、X(0)の虚数部には0が書き込まれる。x(1)は正負の符号を反転されたのちX(1)の虚数部に転記され、X(1)の実数部には0が書き込まれる。x(2)は正負の符号を反転されたのちX(2)の実数部に転記され、X(2)の虚数部には0が書き込まれる。x(3)はそのままX(3)の虚数部に転記され、X(3)の実数部には0が書き込まれる。このように、離散複素指数関数の演算結果を複素単位乗数データ列(+1、−j、−1、+jの任意の値から開始する数列)の値にしたがって順次符号反転および転記を繰り返すのみでこの周波数シフトを行うことができ、指数関数を実際に乗算して演算する必要がなくなり、処理を大幅に簡略化することができる。
【0057】
さらに、前記複素数バッファX(n)の実数部Real(n)、Imaginary(n)のうち一方は必ず0であるため、上記規則に基づいて0になる側が分かっていれば0を記憶するバッファを省略してバッファの記憶領域を実質的に半分にすることも可能である。
【0058】
なお、図4に示す受信部はアンダーサンプリングを利用するためのサンプリング周波数fsと注目領域の中心周波数fcとの関係をfs=4fc/(4k+1)としているが、送信する超音波パルスの周波数を任意に設定できない場合には、一旦中間周波数にダウンコンバートし、上記関係式を満足する周波数としてサンプリングすることも可能である。また、アンダーサンプリングを用いない従来方式のパルスドップラ法においても同様に処理可能である。
【0059】
図13に受信したエコー信号を一旦中間周波信号にダウンコンバートし、この信号を通常のサンプリングによってA/D変換する構成の受信部を示しておく。この図において図4の受信部と同様の回路は同一番号を付す。この例ではトランスデューサ1が送信する超音波パルスの周波数は440kHzである。このエコー信号を一旦12.5kHzの中間周波信号に変換するため、ミキサ12で442.5kHzの信号とミキシングし、バンドパスフィルタ13で高域の信号を除去する。このバンドパスフィルタ13はローパスフィルタでもよい。こののち、この12.5kHzを中心に展開しているエコー信号を50kHzのサンプリング周波数でサンプリングする。これにより、サンプリング周波数fsと注目領域の中心周波数fcとの関係をfs=4fsにすることができる。
【0060】
ここで、パルス波でドップラ周波数(速度)を測定するパルスモードと連続波でドップラ周波数を測定する連続波モードの両モードによる測定が可能なドップラソナーにこの発明を適用した場合について説明する。連続波モードは、対象物との距離が数十メートル以下になったとき連続波を送受信して0.1Hz程度の極めて高精度の周波数分解能を実現するモードである。ただし、連続波はトランスデューサや駆動回路への負担が大きいため高出力にすることができず、使用可能範囲は上記のように数十メートルの範囲に限定される。
【0061】
図14は、パルスモード・連続波モードの切り換えが可能なドップラソナーの概略ブロック図である。このドップラソナーは、送信用トランスデューサ21、受信用トランスデューサ22、送信部23、受信部24および制御部25からなっている。このようにこのドップラソナーは連続波モードの動作を可能にするため、送信用トランスデューサ21と受信用トランスデューサ22を別々に備えている。ただし、受信部24のハードウェア構成は図3に示したものと同じでよい。
【0062】
図15は同受信部24のDSPの動作を示すフローチャートである。まず、超音波パルスを送受信し(s21)、距離を測定する(s22)。これはエコーの伝搬時間遅れによって測定する。この距離に基づいてパルスモード(s24以下)または連続波モード(s40以下)を選択する(s23)。
【0063】
測定された距離が50メートル以上でパルスモードが選択されたとき、まず単一パルスモードで速度を測定する(s24〜27)。すなわち、パルスを送受信し(s24)、受信信号をFFTして(s26)、そのピークから粗ドップラ周波数fdを求める(s27)。このとき同時に距離データもアップデートしておく(s25)。これは精度は悪いものであるがどのような場合でも問題なく速度が得られる。なお、ここではデータは複素数として処理する。
【0064】
つぎにfdがPRFの何倍かを調べる。すなわち、|fd|≧n×PRFを満足する最大整数nを求める(s28)。受信したエコー信号からこのエコー信号の代表値を再サンプリングデータとして求める(s29)。たとえば平均値、エコー信号の中心点の値等を再サンプリングデータとする。この処理は継続的に実行されており、再サンプリングデータは過去から蓄積されているものとする。
【0065】
蓄積されている再サンプリングデータからFFTに必要な点数(32,64,128点等)だけ集め、FFT処理前に窓関数を乗算する(s30)。窓関数の乗算は省略してもよい。そしてFFT演算を行い(s31)、ドップラ周波数を求めるために、周波数パワースペクトルの最大ピークを求める(s32)。このとき補間法によって分解能をあげることもできる。補間法としては、田部井・上田法(電子情報通信学会論文誌A,vol.J70−A,pp.798−805,1987)を用いればよい。このピーク周波数にn・PRFを加算して所望ドップラ周波数を求め(s33)、この値を出力する(s34)。fdが負の場合には−(n+1)・PRFを加算する。
【0066】
距離が50メートル以下で連続波モードが選択された場合には、送信用トランスデューサから連続波を送信し、このエコー信号を受信してこの信号を連続パルスモード時と同様にデモジュレーションするとともに、少ないサンプル点数で周波数分解能を上げるためにデシメーションを行い(s42)、1024点程度のサンプル点数でs30以下のFFT処理を実行する。このように少ないサンプル点数を用いても10秒程度の時間幅のエコー信号をFFTすることができるため、連続波モードでは0.1Hz程度の極めて高精度の周波数分解能を得ることができる。そして、連続波モードの場合でも一時的にパルスモードに復帰して(s40→s43)、距離を測定する。この距離が連続波モードを維持できる範囲内の距離であれば、このパルスエコー信号を連続波のエコー信号と結合してs41にもどる。距離が連続波モードを維持できないほど離れた場合にはパルスモードに移行するためs24に進む。
【0067】
このようにこのドップラソナーでは、距離に応じて連続波モードとパルスモード(連続パルスモード)とを切り換えることができるため、対象物との距離に応じた最も精度の高い速度測定(ドップラ周波数測定)をすることができる。なお、連続波モードとパルスモードとを切り換える距離にはヒステリシスを持たせてもよい。
【0068】
なお、上記実施形態では超音波を用いたドップラソナーについて説明したが、この発明は超音波以外でも電磁波など波動信号のドップラ周波数を測定する場合に広く適用できるものである。さらに、測定されたドップラ周波数に基づいて割り出される変量も速度に限定されない。
【0069】
【発明の効果】
この発明によれば、一定周期で送信されるパルス信号に起因するエコー信号のそれぞれを1つの(再)サンプリングデータとして用いることにより、短時間でのフーリエ解析が可能になり、高い周波数分解能でドップラ周波数を求めることができる。
【0070】
この場合に、エコー信号(パルス信号)の繰り返し周波数PRFよりもドップラ周波数が高い場合には、上記フーリエ解析で得られるドップラ周波数は基本帯域への写像周波数になるが、単一パルスをフーリエ解析して求めた粗ドップラ周波数により真のドップラ周波数の属する帯域が求められるため、この帯域の基底周波数(n×PRF)を前記写像周波数に加算することによって正確なドップラ周波数を求めることができる。
【0071】
このようにこの発明によれば、ドップラ周波数の帯域にかかわらず高い周波数分解能を得ることができるうえ、必要な周波数分解能を得るための時間幅やドップラ周波数の属する帯域を特定するための粗ドップラ周波数平均回数は大きくないため、短時間での測定が可能であり、時間分解能も向上することができる。
【図面の簡単な説明】
【図1】この発明のドップラ周波数の割り出しの方式を説明する図である。
【図2】この発明のドップラ周波数の割り出しで生じる写像を説明する図である。
【図3】この発明の実施形態であるドップラソナーの概略構成図である。
【図4】同ドップラソナーの受信部の構成を示す図である。
【図5】同ドップラソナーのDSPの動作を示すフローチャートである。
【図6】同DSPでPRFで再サンプリングした場合に形成される写像を示す例である。
【図7】この発明の方式によるドップラ周波数計測のシミュレーションを示す図である。
【図8】この発明の方式によるドップラ周波数計測のシミュレーションを示す図である。
【図9】この発明の方式によるドップラ周波数計測のシミュレーションを示す図である。
【図10】この発明の方式によるドップラ周波数計測のシミュレーションを示す図である。
【図11】この発明の方式によるドップラ周波数計測のシミュレーションを示す図である。
【図12】前記DSPにおけるデモジュレーション処理を説明する図である。
【図13】この発明の他の実施形態であるドップラソナーの概略ブロック図である。
【図14】この発明の他の実施形態であるドップラソナーの概略ブロック図である。
【図15】同ドップラソナーのDSPの動作を示すフローチャートである。
【符号の説明】
1…トランスデューサ、2…切換器、3…送信部、4…受信部、5…制御部、15…A/D変換器、16…DSP
21…送信用トランスデューサ、22…受信用トランスデューサ、23…送信部、24…受信部、25…制御部
Claims (5)
- 一定間隔で送信された複数のパルス信号のエコー信号を受信してサンプリングデータとし、このサンプリングデータを用いて離散フーリエ変換を行ってベースバンドのスペクトルのピーク周波数を求める手順と、
パルス信号のエコー信号に基づいて概略のドップラ周波数を求め、この概略のドップラ周波数に基づいて前記ベースバンドのスペクトルが真のスペクトルであるかまたは写像であるかを求め、さらに、写像である場合には真のスペクトルとの周波数差を求める手順と、
前記ピーク周波数を前記周波数差で補正することによって精ドップラ周波数を求める手順と、
を有するドップラ周波数測定方法。 - (1) 周期信号であるパルス信号を送信して対象物で反射したエコー信号を受信し、該エコー信号をフーリエ変換して粗ドップラ周波数を求める
(2) 周期信号であるパルス信号を送信して対象物で反射したエコー信号を受信し、該エコー信号から前記パルス信号の周波数成分を除去してドップラ周波数の波形成分のみをサンプリング周波数PRFの複素数データ(以下「PRFサンプリングデータ」という)として取り出す
(3) 手順(2)をパルス送信周波数PRFで繰り返し実行する
(4) 複数のPRFサンプリングデータを用いてフーリエ変換し、ベースバンド(0Hz〜PRF)におけるスペクトルのピーク周波数fdを求める
(5) 手順(1)で求めた粗ドップラ周波数に基づき、ドップラ周波数のスペクトルの属する第n次帯域(PRF×n〜PRF×(n+1)…n:整数)を割り出す
(6) 手順(4)で求めたfd、および、手順(5)で求めたnを用い、fd+PRF×nの演算によって精ドップラ周波数を求める
の手順を有するドップラ周波数測定方法。 - 手順(1)を複数回繰り返し、求められた複数の粗ドップラ周波数を平均した値を手順(5)の粗ドップラ周波数として用いる請求項2に記載のドップラ周波数測定方法。
- 手順(1)および手順(2)のうち、周期信号であるパルス信号を送信して対象物で反射したエコー信号を受信する手順を共通にした請求項2または請求項3に記載のドップラ周波数測定方法。
- サンプリング周波数fsがエコー信号の注目領域の中心周波数fcに対して
fs=4fc/(4k+1) ,k=0,1,2,…
または
fs=4fc/(4k+3) ,k=0,1,2,…
となるようなサンプリング周波数fsを用いて該エコー信号をサンプリングし、このサンプリングデータ列に対して+1,−j,−1,+jまたは+1,+j,−1,−jを順次乗算することによって、手順(2)のエコー信号から前記パルス信号の周波数成分を除去する手順を処理する請求項2、請求項3または請求項4に記載のドップラ周波数測定方法。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP33713998A JP4077092B2 (ja) | 1998-11-27 | 1998-11-27 | ドップラ周波数測定方法およびドップラソナー |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP33713998A JP4077092B2 (ja) | 1998-11-27 | 1998-11-27 | ドップラ周波数測定方法およびドップラソナー |
Publications (2)
Publication Number | Publication Date |
---|---|
JP2000162317A JP2000162317A (ja) | 2000-06-16 |
JP4077092B2 true JP4077092B2 (ja) | 2008-04-16 |
Family
ID=18305818
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
JP33713998A Expired - Lifetime JP4077092B2 (ja) | 1998-11-27 | 1998-11-27 | ドップラ周波数測定方法およびドップラソナー |
Country Status (1)
Country | Link |
---|---|
JP (1) | JP4077092B2 (ja) |
Families Citing this family (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP4667613B2 (ja) * | 2001-02-07 | 2011-04-13 | 古野電気株式会社 | 信号周波数検出方法 |
JP5341288B2 (ja) * | 2001-04-27 | 2013-11-13 | 古野電気株式会社 | 魚群探知機 |
JP2006038819A (ja) * | 2004-06-25 | 2006-02-09 | Macnica Inc | データ通信装置および方法 |
JP5937548B2 (ja) * | 2013-06-26 | 2016-06-22 | 古野電気株式会社 | ディジタルフィルタおよび信号処理装置 |
CN114325671B (zh) * | 2021-12-29 | 2024-04-09 | 南京世海声学科技有限公司 | 改进的基于水下声信标信号多脉冲累积的抗多普勒方法 |
CN115225178B (zh) * | 2022-07-08 | 2023-08-11 | 西南交通大学 | 一种基于频域零点重置的多普勒效应模拟器的实现方法 |
-
1998
- 1998-11-27 JP JP33713998A patent/JP4077092B2/ja not_active Expired - Lifetime
Also Published As
Publication number | Publication date |
---|---|
JP2000162317A (ja) | 2000-06-16 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
US5662115A (en) | Method for determining the velocity-time spectrum of blood flow | |
EP2570782B1 (en) | Doppler measurement instrument and doppler measurement processing method | |
US20100262009A1 (en) | Method and apparatus for enhancing the detecting and tracking of moving objects using ultrasound | |
US6447455B2 (en) | Ultrasound diagnostic apparatus and method for measuring blood flow velocity using doppler effect | |
JP3581676B2 (ja) | 血流速度測定向け超音波診断装置 | |
JP2012247302A (ja) | フーリエ変換でのサンプル数を削減した、短時間信号のピークパワースペクトルを検出する方法及び装置 | |
CN111693136B (zh) | 一种采用回波信号自相关相位谱的声表面波谐振器频率估计算法 | |
CN101342083B (zh) | 多普勒信号频谱计算方法与装置 | |
US6758815B2 (en) | Apparatus and method for indicating mechanical stiffness properties of body tissue | |
CN104169740A (zh) | 用于检测声学信号的设备以及所属的方法 | |
JP4077092B2 (ja) | ドップラ周波数測定方法およびドップラソナー | |
CN113805166A (zh) | 一种雷达物位计的目标跟踪测距方法及系统 | |
JP3773779B2 (ja) | レーダ信号処理装置 | |
JP4386282B2 (ja) | 水中通信システム | |
JP3850950B2 (ja) | ドップラソナー | |
RU183781U1 (ru) | Устройство определения доплеровского сдвига частоты по информационному фазоманипулированному сигналу путем аппроксимации фазового отклонения | |
JP3881078B2 (ja) | 周波数推定方法、周波数推定装置、ドップラソナーおよび潮流計 | |
JP4249332B2 (ja) | 周波数測定方法およびドップラソナー | |
RU186027U1 (ru) | Устройство определения доплеровского сдвига частоты по информационному фазоманипулированному сигналу путем взвешенной аппроксимации фазового отклонения | |
JP3877749B2 (ja) | 信号処理方法、信号処理装置 | |
JP3964095B2 (ja) | 大気温度測定方法および装置 | |
RU2467350C2 (ru) | Способ и устройство обнаружения сигналов при наличии переменного доплеровского эффекта | |
US5216639A (en) | Method for processing a doppler signal | |
JP2563656B2 (ja) | 超音波ドプラ映像装置 | |
JP3391578B2 (ja) | 相関装置および流れ情報表示装置 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20050606 |
|
A977 | Report on retrieval |
Free format text: JAPANESE INTERMEDIATE CODE: A971007 Effective date: 20070522 |
|
A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20070529 |
|
A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20070724 |
|
TRDD | Decision of grant or rejection written | ||
A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20080129 |
|
A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20080131 |
|
R150 | Certificate of patent or registration of utility model |
Free format text: JAPANESE INTERMEDIATE CODE: R150 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20110208 Year of fee payment: 3 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20110208 Year of fee payment: 3 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20120208 Year of fee payment: 4 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20130208 Year of fee payment: 5 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20130208 Year of fee payment: 5 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20140208 Year of fee payment: 6 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
EXPY | Cancellation because of completion of term |