JP2713761B2 - 時間空間信号のスペクトル推定方法 - Google Patents

時間空間信号のスペクトル推定方法

Info

Publication number
JP2713761B2
JP2713761B2 JP1101272A JP10127289A JP2713761B2 JP 2713761 B2 JP2713761 B2 JP 2713761B2 JP 1101272 A JP1101272 A JP 1101272A JP 10127289 A JP10127289 A JP 10127289A JP 2713761 B2 JP2713761 B2 JP 2713761B2
Authority
JP
Japan
Prior art keywords
spectrum
signal
time
frequency
estimated
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 - Fee Related
Application number
JP1101272A
Other languages
English (en)
Other versions
JPH02280074A (ja
Inventor
清仁 徳田
敦司 深沢
聡 清水
拓朗 佐藤
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Oki Electric Industry Co Ltd
Original Assignee
Oki Electric Industry Co Ltd
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Oki Electric Industry Co Ltd filed Critical Oki Electric Industry Co Ltd
Priority to JP1101272A priority Critical patent/JP2713761B2/ja
Publication of JPH02280074A publication Critical patent/JPH02280074A/ja
Application granted granted Critical
Publication of JP2713761B2 publication Critical patent/JP2713761B2/ja
Anticipated expiration legal-status Critical
Expired - Fee Related legal-status Critical Current

Links

Landscapes

  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)

Description

【発明の詳細な説明】 (産業上の利用分野) 本発明は、水中における位置計測(音響測位)等を行
うソーナの信号処理等において、音源標定における周波
数スペクトル推定(時間的処理)及び波数スペクトル推
定(空間的処理)を行う時間空間信号のスペクトル推定
方法に関するものである。
(従来の技術) 従来、このような分野の技術としては例えば、沖電気
研究開発、53[4](昭61−10)、似鳥・五十嵐著「水
中音響におけるディジタル信号処理技術」p.53−58(文
献1)、及び尾崎統著「時系列論」(昭63−3−20)放
送大学教材、p.49−55、101−105(文献2)に記載され
るものがあった。
前記文献1に記載されているように、ソーナは、水中
音響を用いて3次元空間内(水中)を移動する目標の探
索、位置計測、類別を行うものである。ソーナは動作様
式により、パッシブソーナとアクティブソーナに分けら
れる。パッシブソーナは目標が放射する音波を用い、ま
たアクティブソーナは目標に向かって音波を当て、その
反射波(エコー)を用いて目標の探索、位置計測、類別
を行う。ソーナで用いられる信号処理、つまり時間空間
信号のスペクトル推定方法は、次のように、信号の時間
的特徴(波形、スペクトル等)を抽出するために用いる
時間的処理と、信号の空間的特徴(位置、形状、移動速
度等)を抽出するために用いる空間的処理に分けられ
る。
時間的処理のうち、整合フィルタおよびウィーナフィ
ルタは、それぞれ定められた波形およびスペクトルを持
つ信号が、既知のスペクトルを持つ雑音に埋れている時
に最大のS/N比を得るフィルタである。スペクトル推定
は、信号の周波数に対する強度を推定するものであり、
雑音に埋れた周期的信号(線スペクトル)を検出し、そ
の信号を放射する目標を類別するために用いられる。
空間的処理のうち、ビームフォーミング(BF)は、受
信アレイを構成する多数の受波器の出力信号を用い、空
間を伝搬する波動の方向性を利用して信号のS/Nの改
善、信号の入射方向、および強度(空間スペクトル)の
推定等を行う。遅延、および位相推定は、少数の受波器
で受信される信号の間に生ずる遅延または位相の差を推
定する問題であり、主に目標の位置計測のために用いら
れる。これは、ビームフォーミングの簡約化であると考
えられる。
(発明が解決しようとする課題) しかしながら、上記いずれのスペクトル推定方法であ
っても、受信アレー数が少ない時には、入力データ数が
少ないため、ビームフォーミングの空間的分解能が劣下
し、正しく波数スペクトル(方位)が推定できない。そ
の上、FFT(高速フーリエ変換、Fast Fourier Transfor
m)を用いて周波数スペクトルを推定するので、短時間
では、サンプリングデータ数が少ないため、FFTの特性
から、周波数分解能が劣下するという問題があり、それ
らを解決することが困難であった。
本発明は前記従来技術が持っていた課題として、周波
数分解能と空間的分解能の劣下という点について解決し
た時間空間信号のスペクトル推定方法を提供するもので
ある。
(課題を解決するための手段) 前記課題を解決するために、第1の発明は、任意の周
波数及び空間位置の音響信号を複数の受信アレーで受信
し、そのアレー受信信号から、周波数スペクトル及び波
数スペクトルを推定するパッシブ受波方式における時間
空間信号のスペクトル推定方法において、前記各アレー
受信信号をそれぞれ独立な多変量データとみなし、その
各多変量データ間の相関度を考慮した線形予測を行う多
変量線形予測分析方法を用いて前記周波数スペクトルを
推定する。さらに、その周波数スペクトル推定結果から
導出されるコヒーレンシー(位相差量を評価する関連度
関数)を用いて位相差分析することにより、前記波数ス
ペクトルを推定するようにしたものである。
第2の発明は、空間位置が既知である音響信号を用い
て第1の発明を適用することにより、アクティブ受波方
式の信号処理を行うようにしている。
(作 用) 第1、第2の発明によれば、以上のように時間空間信
号のスペクトル推定方法を構成したので、多変量線形予
測分析では、各アレー受信信号を独立な信号、即ち多変
量データとみなし、各々互いの相関度を考慮して線形予
測することにより、少ない受信アレー数で、しかも短時
間かつ高分解能の周波数スペクトル推定を行なう。ま
た、コヒーレンシーを用いた位相差分析では、短い受信
アレー間隔で、高分解能の波数スペクトル推定を可能に
させる。従って、前記課題を解決できるのである。
(実施例) 第1図は、本発明の実施例を示すもので、多変量線形
予測分析を用いた方位(波数スペクトル)推定の機能ブ
ロック図である。
この機能ブロックでは、パッシブ受波あるいはアクテ
ィブ受波を行う受信アレー部10を備え、その出力側には
ディジタル・シグナル・プロセッサ(DSP)等で構成さ
れる時間的処理帯域制限フィルタ部20、時間的処理スペ
クトル分析部30及び空間的処理方位決定部40が接続され
ている。
受信アレー部10は、任意の周波数及び空間位置の音響
信号、あるいは空間位置が既知である音響信号を受信し
てそれを電気信号に変換する複数の受信アレー10−1〜
10−4と、各受信アレー出力をディジタル信号に変換す
る複数のアナログ/ディジタル変換器(以下、A/D変換
器という)11−1〜11−4とを備えている。時間的処理
帯域制限フィルタ部20は、A/D変換器11−1〜11−4の
出力側に接続された複数の帯域制限フィルタ(以下、BP
Fという)21−1〜21−4で構成されている。BPF21−1
〜21−4は、ディジタルフィルタ等で構成される。時間
的処理スペクトル分析部30は、BPF21−1〜21−4の出
力を入力して多変量線形予測分析を行う機能を有し、さ
らに、空間的処理方位決定部40は、多変量線形予測分析
結果を入力して位相差分析を行い、推定方位kを出力す
る機能を有している。
次に、動作を説明する。
音源からの音響信号、例えば時間空間信号s(t,)
は、受信アレー10−1〜10−4で受信されて電気信号に
変換され、その電気信号がA/D変換器11−1〜11−4に
より、サンプリング周波数fsでサンプリングされ、時間
的離散値s(n,)となる。ここで、nはサンプリン
グしたデータの数、は受信アレー10−1〜10−4の3
次元位置ベクトル、mは受信アレー10−1〜10−4の番
号であり、n=1,2,…,N、m=1,2,…,M(但し、第1図
では受信アレー数M=4の場合を示している)となる。
受信アレー部10の出力信号は、BPF21−1〜21−4に
おいて各々帯域制限された後、時間的処理スペクトル分
析部30へ入力される。
時間的処理スペクトル分析部30では、以下に述べる多
変量線形予測分析を行う。なお、多変量時系列に関する
理論は、前記文献2に記載されている。
多変量線形予測分析を行う場合、まずM個の受信アレ
ー出力を離散的多変量時系列とみなし、受信アレー出力
行列(n)を、 (n)=[s(n,),s(n,), …,s(n,)] ……(1) 但し、n=1,2,…,N(N;データ数)t;転置行列 とする。また受信アレー出力行列(n)は、過去P個
の時系列行列(n−1),…,(n−P)と、下記
(2)式に示す白色雑音時系列行列(予測誤差行列とも
いう)(n)とを用いて、下記(3)式に示す多変量
自己回帰モデル(以下、M次元AR(P)モデルという)
を満足するものとする。
(n)=[e(n,),e(n,), …,e(n,)] ……(2) 但し、(p);予測係数行列 p=0,1,…,P (3)式において、予測係数行列(p)は次式のよ
うになる。
(0)=− 但し、;単位行列 即ち、 また、予測誤差行列(n)は次式のようになる。
E[(n)]= 但し、E[ ];アンサンブル (集合)平均 今、(3)式を要素表示すると、一般に、 但し、q=1,2,…,M となる。時間平均が特定の時刻に無関係で時間差τi,j
のみに依存する定常場では、予測誤差の相互相関関数Re
(τi,j:ri,rj)を(5)式のように定義できる。
但し、−T〜T;時間フレーム τi,j=τ−τ ここで、相互相関関数を Rs(τi,j)= E[s(n+τi,j) ・s(n,ξ)] ……(6) と定義すると、(5)式は(5−1)式のようになる。
さらに、(5−1)式の相互相関関数Rs(τi,j
i,)と周波数クロススペクトル(入出力間のスペ
クトルの関係)の間に、ウィナー・キンチン(Wiener−
Khintchine)の定理を用いると、クロススペクトルPs
q,ξ(ω)は し、周波数ω=2πf Rs(τi,ji,) …自己相関となり、(7)式は(8)式のようになる。
ここで、係数Ai,q(ω)を、フーリエ変換式を用い
とおくと、定常場のクロススペクトル間には、 但し、*;複素共役 の関係が成り立つ。この(10)式を行列表示するため
に、クロススペクトル行列e,s,Aを次式のようにお
くと、 (10)式は、 e=・s・*t となる。従って、クロススペクトル行列(入力信号のパ
ワースペクトル) s=()-1・e・(*t-1 ……(11) により、スペクトルが推定できる。
今、空間上の任意の点に置かれた1受信アレー
(例えば、10−1)に平面波が入射した時、その周波数
スペクトルSm(ω)は、理論的には S(ω,rm)=EXP[−j] ・δ(ω−ω) ……(12) 但し、;任意波数ベクトル(方位) と表わされ、周波数ω=ωの時には、 S(ω0,)=EXP[−j) ……(13) となる。
以上のような多変量線形予測分析により、ωを任意
周波数として積分δ(ω−ω)の近似スペクトルS
(ω)が求まり、(13)式は S(ω0,)= S(ω)EXP[−j) ……(14) となる。各受信アレー10−1〜10−4で近似スペクトル
S(ω)は共通であるから、時間的処理スペクトル分
析部30の出力として、受信アレー10−1〜10−4ごと
に、(13)式の周波数スペクトルS(ω0,)が出力
され、それが空間的処理方位決定部40に入力される。
空間的処理方位決定部40では、各周波数スペクトルS
(ω0,)について、以下のように、コヒーレンシー
値を用いて位相差を算出し、推定方位を求める。
即ち、(11)式のクロススペクトル行列sの要素か
らコヒーレンシーを、 Cij(ω)=Psij(ω)/Psii(ω) ……(15) 但し、i,j=1,2,…,M と定義し、一般にクロススペクトルPs(ω)を Ps(ω)=As(ω)e−jθ(ω) ……(16) 但し、 As(ω)=sqrt{real[Ps(ω)]+ imag[Ps(ω)]} θ(ω)=atan{imag[Ps(ω)]/ real[Ps(ω)]} とすると、(15)式のコヒーレンシーは(17)式のよう
になる。
Cij(ω)=[Asij(ω)/Asii(ω)] ・EXP[−j(θij(ω)−θii(ω))] ……(17) この(17)式中の振幅成分 Asij(ω)/Asij(ω) を無視して位相成分のみを考えると、位相差Δij(ω)
は Δij(ω)=θij(ω)−θii(ω) ……(18) で求まる。ここで、周波数ω、互いに異なる空間位置
i,におけるコヒーレンシーの理論式は、 Cij(ω)=EXP[−jΘij] ……(19) 但し、ijij ij となる。今、空間方位(波数スペクトル)を推定する
ために、i=1、j=2,3,4の場合を考え、実際の位相
差は、多変量線形予測分析により、(18)式の位相差Δ
ij(ω)で求まったものとすると、 ijij=Δij(ω) ……(20) 但し、i=1、j=2,3,4 が成り立つ。波数ベクトルのxyz成分を(k0x,k0y,k
0z)とすると、入射平面波の音源方位を意味する波形ベ
クトルのxyz成分は、 として求まり、音源の方位が推定できる。
以上のように、本実施例では、時間的処理スペクトル
分析部30により、多変量線形予測分析で周波数スペクト
ルを推定するので、少ない受信アレー数(例えば、最低
4個)で、しかも短時間推定が可能となる。さらに、空
間的処理方位決定部40により、位相差分析で波数スペク
トル(方位)kを推定するので、各受信アレー間隔をλ
/2(但し、λ;入射波の波長、k;波数、λ=2π/)
以下にすることが可能となり、装置の小型化が図れる。
なお、本発明は上記の実施例に限定されず、第1図の
各部20,30,40をマイクロコンピュータを用いたプログラ
ム制御等で実行したり、あるいは本発明を空中の音源標
定等に適用する等、種々の変形と応用が可能である。
(発明の効果) 以上詳細に説明したように、第1、第2の発明によれ
ば、多変量線形予測分析により周波数スペクトルを推定
するので、少ない受信アレー数で、しかも短時間かつ高
分解能の推定が可能となる。その上、前記周波数スペク
トル推定結果から導出されるコヒーレンシーを用いて位
相差分析することにより、波数スペクトルの推定を行う
ので、高分解能を維持しつつ受信アレー間隔を、入射波
の波長λの1/2以下にすることが可能となり、装置の小
型化が図れるばかりか、高精度な時間空間信号スペクト
ル分析方法を提供することができる。
【図面の簡単な説明】
第1図は本発明の実施例を示す多変量線形予測分析を用
いた方位推定の機能ブロック図である。 10……受信アレー部、10−1〜10−4……受信アレー、
11−1〜11−4……A/D変換器、20……時間的処理帯域
制限フィルタ部、21−1〜21−4……BPF、30……時間
的処理スペクトル分析部、40……空間的処理方位決定
部。
フロントページの続き (72)発明者 佐藤 拓朗 東京都港区虎ノ門1丁目7番12号 沖電 気工業株式会社内 (56)参考文献 特開 平2−268282(JP,A) 特開 平2−17481(JP,A) 特開 平2−201186(JP,A) 特開 平2−74886(JP,A) 特開 昭61−2087(JP,A) 特公 平1−38270(JP,B2)

Claims (2)

    (57)【特許請求の範囲】
  1. 【請求項1】任意の周波数及び空間位置の音響信号を複
    数の受信アレーで受信し、そのアレー受信信号から、周
    波数スペクトル及び波数スペクトルを推定する時間空間
    信号のスペクトル推定方法において、 前記各アレー受信信号をそれぞれ独立な多変量データと
    みなし、その各多変量データ間の相関度を考慮した線形
    予測を行う多変量線形予測分析方法を用いて前記周波数
    スペクトルを推定し、 その周波数スペクトル推定結果から導出されるコヒーレ
    ンシーを用いて位相差分析することにより、前記波数ス
    ペクトルを推定することを特徴とする時間空間信号のス
    ペクトル推定方法。
  2. 【請求項2】音響信号は空間位置が既知である請求項1
    記載の時間空間信号のスペクトル推定方法。
JP1101272A 1989-04-20 1989-04-20 時間空間信号のスペクトル推定方法 Expired - Fee Related JP2713761B2 (ja)

Priority Applications (1)

Application Number Priority Date Filing Date Title
JP1101272A JP2713761B2 (ja) 1989-04-20 1989-04-20 時間空間信号のスペクトル推定方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP1101272A JP2713761B2 (ja) 1989-04-20 1989-04-20 時間空間信号のスペクトル推定方法

Publications (2)

Publication Number Publication Date
JPH02280074A JPH02280074A (ja) 1990-11-16
JP2713761B2 true JP2713761B2 (ja) 1998-02-16

Family

ID=14296249

Family Applications (1)

Application Number Title Priority Date Filing Date
JP1101272A Expired - Fee Related JP2713761B2 (ja) 1989-04-20 1989-04-20 時間空間信号のスペクトル推定方法

Country Status (1)

Country Link
JP (1) JP2713761B2 (ja)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102866383A (zh) * 2012-10-22 2013-01-09 重庆大学 一种基于空域自适应滤波的波达方向估计方法

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH05196717A (ja) * 1992-01-20 1993-08-06 Tech Res & Dev Inst Of Japan Def Agency ソーナー信号検出処理方式
JP2779435B2 (ja) * 1995-10-31 1998-07-23 防衛庁技術研究本部長 ソーナー信号検出処理装置
CN110441731A (zh) * 2019-07-27 2019-11-12 西南电子技术研究所(中国电子科技集团公司第十研究所) 窄带电磁辐射源目标定位方法
CN112034440B (zh) * 2020-09-01 2023-12-22 中国科学院声学研究所 一种基于波数谱能量累积分布的目标深度辨识方法及系统

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102866383A (zh) * 2012-10-22 2013-01-09 重庆大学 一种基于空域自适应滤波的波达方向估计方法

Also Published As

Publication number Publication date
JPH02280074A (ja) 1990-11-16

Similar Documents

Publication Publication Date Title
CN103713276B (zh) 基于最小互熵谱分析的波达方向估计方法
US7773458B2 (en) Systems and methods for detection and analysis of amplitude modulation of underwater sound
US5251186A (en) Preprocessor and adaptive beamformer for linear-frequency modulation active signals
US6446025B1 (en) Multiple propagation wave parameter measuring method and apparatus and machine-readable recording medium recording multiple propagation wave parameter measuring program
JP2713761B2 (ja) 時間空間信号のスペクトル推定方法
US6654315B1 (en) Sonar display system and method
JP7056212B2 (ja) 方位推定方法および装置
GB2104218A (en) Detecting harmonically-rich acoustic sources
Touzé et al. Double-Capon and double-MUSICAL for arrival separation and observable estimation in an acoustic waveguide
JPH10332818A (ja) 目標位置局限方法
Grondin et al. A study of the complexity and accuracy of direction of arrival estimation methods based on GCC-PHAT for a pair of close microphones
CN113189570B (zh) 一种基于复域压缩感知的阵列信号处理方法及系统
Wu et al. Switched-element direction finding
CN110361696B (zh) 基于时间反转技术的封闭空间声源定位方法
CN113281755A (zh) 一种基于深度和声速扫描的高精度单波束测深方法
JP2690606B2 (ja) 音源数決定方法
CN112666558A (zh) 一种适用于汽车fmcw雷达的低复杂度music测向方法及装置
JP2779435B2 (ja) ソーナー信号検出処理装置
JP3030332B2 (ja) ノーマルモード固有値推定方法
JP2685632B2 (ja) 音源特徴抽出方法
Wang et al. On performance improvement of tone frequency estimation
JPH08339200A (ja) 到来受信信号の周波数分析方法
Hioka et al. Multiple-speech-source localization using advanced histogram mapping method
Akins et al. Modal-MUSIC: A passive mode estimation algorithm for partially spanning arrays
RU2470315C1 (ru) Способ компьютерно-интерферометрического обнаружения-пеленгования радиосигналов с расширенным спектром

Legal Events

Date Code Title Description
LAPS Cancellation because of no payment of annual fees