JPH09133750A - 目標運動解析方法 - Google Patents
目標運動解析方法Info
- Publication number
- JPH09133750A JPH09133750A JP29418495A JP29418495A JPH09133750A JP H09133750 A JPH09133750 A JP H09133750A JP 29418495 A JP29418495 A JP 29418495A JP 29418495 A JP29418495 A JP 29418495A JP H09133750 A JPH09133750 A JP H09133750A
- Authority
- JP
- Japan
- Prior art keywords
- target
- initial value
- information
- state quantity
- motion analysis
- 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.)
- Withdrawn
Links
- 230000033001 locomotion Effects 0.000 title claims abstract description 80
- 238000004458 analytical method Methods 0.000 title claims description 90
- 238000004364 calculation method Methods 0.000 claims abstract description 60
- 239000011159 matrix material Substances 0.000 claims abstract description 29
- 238000001514 detection method Methods 0.000 claims abstract description 26
- 238000005259 measurement Methods 0.000 claims description 22
- 238000000034 method Methods 0.000 claims description 13
- 238000012545 processing Methods 0.000 claims description 5
- 238000011156 evaluation Methods 0.000 abstract description 14
- 230000005855 radiation Effects 0.000 abstract 1
- 238000010586 diagram Methods 0.000 description 5
- 238000007796 conventional method Methods 0.000 description 1
- 239000004615 ingredient Substances 0.000 description 1
- 238000011835 investigation Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
- 230000017105 transposition Effects 0.000 description 1
Landscapes
- Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)
Abstract
(57)【要約】
【課題】 非線形方程式を解いて目標の運動解析を行な
う際の解の収束を速める。 【解決手段】 目標の放射音に対応する受信信号が、方
位情報算出部２２と周波数情報算出部２３に入力され、
観測方位角θm （t)と観測周波数成分νkm（t)が算出さ
れ、初期値算出部２４及びマニューバ検出部２５に渡さ
れる。マニューバ検出部２５は、マニューバ発生時点で
の内部状態量の推定値VX0 * の速度に関する成分をゼロ
とした内部状態量の推定値VXdo* と推定情報行列Σ
VXd0* -1を、初期値算出部２４に送る。初期値算出部
は、推定値VXdo* と行列ΣVXd0* -1を先見情報として、
その情報と目標に対する観測量とで構成された擬似線形
方程式を解き、解析初期値VX00 * を求めて内部状態量推
定部２６に送る。内部状態量推定部２６は、初期値VX00
* を適用したと評価関数で目標の運動解析を行う。
う際の解の収束を速める。 【解決手段】 目標の放射音に対応する受信信号が、方
位情報算出部２２と周波数情報算出部２３に入力され、
観測方位角θm （t)と観測周波数成分νkm（t)が算出さ
れ、初期値算出部２４及びマニューバ検出部２５に渡さ
れる。マニューバ検出部２５は、マニューバ発生時点で
の内部状態量の推定値VX0 * の速度に関する成分をゼロ
とした内部状態量の推定値VXdo* と推定情報行列Σ
VXd0* -1を、初期値算出部２４に送る。初期値算出部
は、推定値VXdo* と行列ΣVXd0* -1を先見情報として、
その情報と目標に対する観測量とで構成された擬似線形
方程式を解き、解析初期値VX00 * を求めて内部状態量推
定部２６に送る。内部状態量推定部２６は、初期値VX00
* を適用したと評価関数で目標の運動解析を行う。
Description
【０００１】
【発明の属する技術分野】本発明は、目標から放射され
る音を、運動可能な艦船等の航走体に取付けられた受波
器センサアレイで受信し、雑音に乱された観測量から、
相対移動する目標の位置と速度等に関する内部状態量の
推定を行う目標運動解析方法に関するものである。
る音を、運動可能な艦船等の航走体に取付けられた受波
器センサアレイで受信し、雑音に乱された観測量から、
相対移動する目標の位置と速度等に関する内部状態量の
推定を行う目標運動解析方法に関するものである。
【０００２】
【従来の技術】従来、このような分野の技術としては、
例えば次の文献に示されたものがあった。 文献１；ACASSP,1989 IEEE、J.M.Passrieux,D.Pillion,
P.Blanc-Benon and C.Jauffret著“Target Motion Anal
ysis with Bearing and Frequencies Measurements via
Instrumental Variable Estimator”P.2645-2648 文献２；GEC JOURNAL OF RESEARCH 、３［３］（１９８
５）、M.S.Woolfson著“An Evaluation of Manoeuvre D
etector Algorithms”P.181-190 図２は、従来の目標運動解析方法を示す幾何学的説明図
であり、観測系と運動系が示されている。図２における
（Ｘ，Ｙ）は原点ｏの固定座標系であり、その固定座標
系中に航走体１と目標２とが示されている。同図中のVx
1 （t)は時刻ｔにおける航走体１の位置ベクトル、Vx2
（t)は時刻ｔにおける目標の位置ベクトル、ｒ（t)は時
刻ｔにおける航走体１と目標２間の距離‖Vx2 （t)−Vx
1 （t)‖、θ（t)は時刻ｔにおける航走体１から見た目
標２のＹ軸を基準とする方位角をそれぞれ示している。
ここで“‖‖”はベクトルのノルムを表わす。目標運動
解析は、目標２が等速直線運動を行っていると仮定し、
雑音に乱された音源の目標２に対する方位角と受信音の
周波数の観測時系列から、目標２の位置と速度を同定す
るものである。従来技術を説明する前に、目標運動解析
方法の動作原理について説明する。時刻ｔにおける方位
角θ（t)の観測値及び音源の第ｋ固有周波数成分ｆk の
ドップラシフト周波数νk (t) の観測雑音を含んだ観測
値を、次の（１）式のように表わす。
例えば次の文献に示されたものがあった。 文献１；ACASSP,1989 IEEE、J.M.Passrieux,D.Pillion,
P.Blanc-Benon and C.Jauffret著“Target Motion Anal
ysis with Bearing and Frequencies Measurements via
Instrumental Variable Estimator”P.2645-2648 文献２；GEC JOURNAL OF RESEARCH 、３［３］（１９８
５）、M.S.Woolfson著“An Evaluation of Manoeuvre D
etector Algorithms”P.181-190 図２は、従来の目標運動解析方法を示す幾何学的説明図
であり、観測系と運動系が示されている。図２における
（Ｘ，Ｙ）は原点ｏの固定座標系であり、その固定座標
系中に航走体１と目標２とが示されている。同図中のVx
1 （t)は時刻ｔにおける航走体１の位置ベクトル、Vx2
（t)は時刻ｔにおける目標の位置ベクトル、ｒ（t)は時
刻ｔにおける航走体１と目標２間の距離‖Vx2 （t)−Vx
1 （t)‖、θ（t)は時刻ｔにおける航走体１から見た目
標２のＹ軸を基準とする方位角をそれぞれ示している。
ここで“‖‖”はベクトルのノルムを表わす。目標運動
解析は、目標２が等速直線運動を行っていると仮定し、
雑音に乱された音源の目標２に対する方位角と受信音の
周波数の観測時系列から、目標２の位置と速度を同定す
るものである。従来技術を説明する前に、目標運動解析
方法の動作原理について説明する。時刻ｔにおける方位
角θ（t)の観測値及び音源の第ｋ固有周波数成分ｆk の
ドップラシフト周波数νk (t) の観測雑音を含んだ観測
値を、次の（１）式のように表わす。
【０００３】
【数１】 ここで、目標２に対する航走体１上の観測者から見た相
対位置ベクトルVx（t)を次の（２）式とし、各時刻ｔ＝
ｔ1 ，ｔ2 ，…，ｔn における内部状態量VX（t)を
（３）式のように定義する。但し、Ｔは転置を表わし、
“’”は時間微分を表す。即ち、内部状態量VX（t)は、
目標２の位置及び速度と、目標２の音源固有周波数成分
とを含んでいる。 Vx(t) ＝Vx2 (t) −Vx1 (t) ＝［ｒx (t) ｒy (t) ］T ・・・（２） VX（t)＝［Vx（t)T Vx’(t) T １／ｆ1 … １／ｆp ］T ・・・（３） 時刻ｔ＝ｔ0 における内部状態量VX0 をVX(t0 ）とおく
とき、VX(t）とVX0 の関係は、次の（４）に示す状態遷
移式で与えられる。
対位置ベクトルVx（t)を次の（２）式とし、各時刻ｔ＝
ｔ1 ，ｔ2 ，…，ｔn における内部状態量VX（t)を
（３）式のように定義する。但し、Ｔは転置を表わし、
“’”は時間微分を表す。即ち、内部状態量VX（t)は、
目標２の位置及び速度と、目標２の音源固有周波数成分
とを含んでいる。 Vx(t) ＝Vx2 (t) −Vx1 (t) ＝［ｒx (t) ｒy (t) ］T ・・・（２） VX（t)＝［Vx（t)T Vx’(t) T １／ｆ1 … １／ｆp ］T ・・・（３） 時刻ｔ＝ｔ0 における内部状態量VX0 をVX(t0 ）とおく
とき、VX(t）とVX0 の関係は、次の（４）に示す状態遷
移式で与えられる。
【０００４】
【数２】 内部状態量VX（t)の推定は、時刻ｔ＝ｔ1 ，ｔ2 ，…，
ｔn における観測値の組から構成される次の（５）式の
評価関数Ｌ（VX0 ）を最小とするVX0 を求めることによ
って行う。
ｔn における観測値の組から構成される次の（５）式の
評価関数Ｌ（VX0 ）を最小とするVX0 を求めることによ
って行う。
【０００５】
【数３】 つまり、航走体１及び目標２の運動は、音速ｃに比べて
はるかに小さいと仮定している。また、“Ｅ［］”は、
アンサンブル（集合）平均を示している。
はるかに小さいと仮定している。また、“Ｅ［］”は、
アンサンブル（集合）平均を示している。
【０００６】Ｌ（VX0 ）の最少化において、次の（７）
なる非線形方程式を繰り返し法を用いて解く必要があ
り、そのための初期値VX00を定める必要がある。よっ
て、観測誤差が微小であるとして、（８）の方程式を導
く。 ∂L ( VX0 ) / ∂VX0 ＝ 0 ・・・（７） Vb ＝ MB VX0 ＋ MW ( VX(t) ) VN ・・・（８）
なる非線形方程式を繰り返し法を用いて解く必要があ
り、そのための初期値VX00を定める必要がある。よっ
て、観測誤差が微小であるとして、（８）の方程式を導
く。 ∂L ( VX0 ) / ∂VX0 ＝ 0 ・・・（７） Vb ＝ MB VX0 ＋ MW ( VX(t) ) VN ・・・（８）
【数４】
【数５】 文献１では、この方程式の重み行列MW（VX(t))の未知要
素を適当な設定値で置き換えた近似重み行列MW* を用い
て、VX0 に関する擬似線形方程式を解くことで、目標運
動解析の解析初期値VX00 * VX00 * ＝（Σ-1/2MW* -1MB )+ Σ-1/2MW* -1Vb ・・・（１０） を定めている。ここで、“+ ”は、一般逆行列を表し、
Σ-1/2はΣのコレスキー分解で定まる上三角行列の逆行
列を示している。続く目標運動解析では、（１０）式で
与えられるVX00 * を初期値として（５）式の評価関数Ｌ
（VX0 ）を最小化する最適な内部状態量VX0 * を、
（７）式の非線形方程式を解くことによって求める。ま
た、求められた内部状態量の推定値VX0 * の推定情報行
列ΣVX0* -1は次の（１１）式で定義され、時刻ｔにおけ
るその推定情報行列ΣVX0* -1は、（４）式の関係から
（１２）式となる。
素を適当な設定値で置き換えた近似重み行列MW* を用い
て、VX0 に関する擬似線形方程式を解くことで、目標運
動解析の解析初期値VX00 * VX00 * ＝（Σ-1/2MW* -1MB )+ Σ-1/2MW* -1Vb ・・・（１０） を定めている。ここで、“+ ”は、一般逆行列を表し、
Σ-1/2はΣのコレスキー分解で定まる上三角行列の逆行
列を示している。続く目標運動解析では、（１０）式で
与えられるVX00 * を初期値として（５）式の評価関数Ｌ
（VX0 ）を最小化する最適な内部状態量VX0 * を、
（７）式の非線形方程式を解くことによって求める。ま
た、求められた内部状態量の推定値VX0 * の推定情報行
列ΣVX0* -1は次の（１１）式で定義され、時刻ｔにおけ
るその推定情報行列ΣVX0* -1は、（４）式の関係から
（１２）式となる。
【０００７】
【数６】 文献２によれば、目標が変針或いは変速等のマニューバ
を行った場合には、目標運動解析の仮定条件が成立しな
くなるため、（１３）式の推定残差 Vεの統計分布を調
べてマニューバの検出を行う。 Vε＝VZm −VZ（VX0 * ） ・・・（１３） 目標の状態量VX0 に関して何らかの値やその精度等の先
見情報が得られている場合には、これら先見状態量及び
その状態量の情報行列をそれぞれVXd0 * 及びΣVX0* -1と
して、次の（１４）式の先見情報を含んだ評価関数を最
小化することにより、先見情報を活かした内部状態量VX
（t)の推定を行う。 Ｈ（VX0 ）＝（VZm −VZ（VX0 ））T Σ-1（VZm −VZ（VX0 ）） ＋（VX0 −VXd0 * ）ΣVX0* -1（VX0 −VXd0 * ） ・・・（１４） この（１４）式を最小化する際の初期値はVXd0 * を用い
るが、先見情報が未知の成分については適当な値を仮定
する。例えば時刻ｔ＝ｔd でマニューバが検出された場
合を想定すると、時刻ｔd における目標２の相対位置な
らびに固有周波数に関する推定状態量及びその状態量の
推定情報行列は、時刻ｔd 以降の運動解析における先見
情報として考えられる。よって、推定状態量及び該推定
情報行列の速度に関する成分をゼロと仮定したものをそ
れぞれVXd0 * 及びΣVX0* -1と置き、時刻ｔd 以降の観測
値を用いて運動解析を行うように初期化する。即ち、マ
ニューバの検出時刻ｔd を改めてｔ＝ｔ0 とし、VXd0 *
を初期値として、（１４）式を最少化するVX0 を求め、
そして内部状態量VX（t ）を推定する。
を行った場合には、目標運動解析の仮定条件が成立しな
くなるため、（１３）式の推定残差 Vεの統計分布を調
べてマニューバの検出を行う。 Vε＝VZm −VZ（VX0 * ） ・・・（１３） 目標の状態量VX0 に関して何らかの値やその精度等の先
見情報が得られている場合には、これら先見状態量及び
その状態量の情報行列をそれぞれVXd0 * 及びΣVX0* -1と
して、次の（１４）式の先見情報を含んだ評価関数を最
小化することにより、先見情報を活かした内部状態量VX
（t)の推定を行う。 Ｈ（VX0 ）＝（VZm −VZ（VX0 ））T Σ-1（VZm −VZ（VX0 ）） ＋（VX0 −VXd0 * ）ΣVX0* -1（VX0 −VXd0 * ） ・・・（１４） この（１４）式を最小化する際の初期値はVXd0 * を用い
るが、先見情報が未知の成分については適当な値を仮定
する。例えば時刻ｔ＝ｔd でマニューバが検出された場
合を想定すると、時刻ｔd における目標２の相対位置な
らびに固有周波数に関する推定状態量及びその状態量の
推定情報行列は、時刻ｔd 以降の運動解析における先見
情報として考えられる。よって、推定状態量及び該推定
情報行列の速度に関する成分をゼロと仮定したものをそ
れぞれVXd0 * 及びΣVX0* -1と置き、時刻ｔd 以降の観測
値を用いて運動解析を行うように初期化する。即ち、マ
ニューバの検出時刻ｔd を改めてｔ＝ｔ0 とし、VXd0 *
を初期値として、（１４）式を最少化するVX0 を求め、
そして内部状態量VX（t ）を推定する。
【０００８】図３は、従来の目標運動解析方法に用いら
れる目標運動解析装置の機能ブロック図である。この機
能ブロック図における目標運動解析装置は、入力端子１
１を介して図示しない受波器アレイからの受信信号を入
力する方位情報算出部１２及び周波数情報算出部１３を
備えている。方位情報算出部１２は、目標２の放射音の
到来方位を受信信号から算出するものであり、周波数情
報算出部１３は、目標２の放射音の周波数成分を算出す
るものである。方位情報算出部１２と周波数情報算出部
１３の出力側は、共に初期値算出部１４とマニューバ検
出部１５に接続されている。初期値算出部１４は、方位
情報と周波数情報から初期状態量を算出する機能を有
し、マニューバ検出部１５は、方位情報算出部１２及び
周波数情報算出部１３か与えられた方位情報と周波数情
報から、目標２のマニューバ検出を行う機能を有してい
る。初期値算出部１４とマニューバ検出部１５には、内
部状態量推定部１６が接続されている。内部状態量推定
部１６は目標の運動を解析するものであり、その解析結
果は、出力端子１７を介して出力されるようになってい
る。次に、マニューバ検出時の目標２の位置情報を先見
情報とする場合における図３の装置の動作を説明する 目標２の放射した音が受波器センサアレイで受信され、
その受信信号が入力端子１１から入力される。受信信号
は方位情報算出部１２と周波数情報算出部１３とに入力
される。方位情報算出部１２では、受信信号から目標２
に対する観測方位角θm （t)を算出し、周波数情報算出
部１３ではトップラシフト観測周波数成分νkm（t)（k=
1,2,…,p）を算出し、それら算出結果が、初期値算出部
１４及びマニューバ検出部１５に送られる。各時刻ｔ＝
ｔ1 ，ｔ2 ，…，ｔn におけるｎ組の観測方位角θ
m （t)とドップラシフト観測周波数成分νkm（t)（k=1,
2,…,p）が、初期値算出部１４に入力されると、該初期
値算出部１４は（１０）を用いて、運動解析の初期値VX
00 * を求め、その初期値VX00 * を内部状態量推定部１６
に送る。内部状態量推定部１６は、そのVX00 * を運動解
析の初期値とし、（１４）式の非線形最小二乗間題を解
いて内部状態量に対する推定値VX0 * を求める。そし
て、内部状態量推定部１６は、目標運動の内部情報量Vx
2 （t)の推定値Vx2 * （t)を算出して、結果を出力端
子１７から出力するともに、内部状態量に対する推定値
VX0 * をマニューバ検出部１５に送る。
れる目標運動解析装置の機能ブロック図である。この機
能ブロック図における目標運動解析装置は、入力端子１
１を介して図示しない受波器アレイからの受信信号を入
力する方位情報算出部１２及び周波数情報算出部１３を
備えている。方位情報算出部１２は、目標２の放射音の
到来方位を受信信号から算出するものであり、周波数情
報算出部１３は、目標２の放射音の周波数成分を算出す
るものである。方位情報算出部１２と周波数情報算出部
１３の出力側は、共に初期値算出部１４とマニューバ検
出部１５に接続されている。初期値算出部１４は、方位
情報と周波数情報から初期状態量を算出する機能を有
し、マニューバ検出部１５は、方位情報算出部１２及び
周波数情報算出部１３か与えられた方位情報と周波数情
報から、目標２のマニューバ検出を行う機能を有してい
る。初期値算出部１４とマニューバ検出部１５には、内
部状態量推定部１６が接続されている。内部状態量推定
部１６は目標の運動を解析するものであり、その解析結
果は、出力端子１７を介して出力されるようになってい
る。次に、マニューバ検出時の目標２の位置情報を先見
情報とする場合における図３の装置の動作を説明する 目標２の放射した音が受波器センサアレイで受信され、
その受信信号が入力端子１１から入力される。受信信号
は方位情報算出部１２と周波数情報算出部１３とに入力
される。方位情報算出部１２では、受信信号から目標２
に対する観測方位角θm （t)を算出し、周波数情報算出
部１３ではトップラシフト観測周波数成分νkm（t)（k=
1,2,…,p）を算出し、それら算出結果が、初期値算出部
１４及びマニューバ検出部１５に送られる。各時刻ｔ＝
ｔ1 ，ｔ2 ，…，ｔn におけるｎ組の観測方位角θ
m （t)とドップラシフト観測周波数成分νkm（t)（k=1,
2,…,p）が、初期値算出部１４に入力されると、該初期
値算出部１４は（１０）を用いて、運動解析の初期値VX
00 * を求め、その初期値VX00 * を内部状態量推定部１６
に送る。内部状態量推定部１６は、そのVX00 * を運動解
析の初期値とし、（１４）式の非線形最小二乗間題を解
いて内部状態量に対する推定値VX0 * を求める。そし
て、内部状態量推定部１６は、目標運動の内部情報量Vx
2 （t)の推定値Vx2 * （t)を算出して、結果を出力端
子１７から出力するともに、内部状態量に対する推定値
VX0 * をマニューバ検出部１５に送る。
【０００９】次ステップにおいて、各時刻ｔ＝ｔn+1 ，
ｔn+2 ，…，ｔm （＞ｔn ) の（ｍ−ｎ）組の観測方位
角θm （t)及びドップラシフト観測周波数成分νkm（t)
（k=1,2,…,p) が初期値算出部１４及びマニューバ検出
部１５に入力されると、そのマニューバ検出部１５で
は、（１３）式で定義される推定算差の統計分布を調べ
る。その結果、例えばマニューバの発生が検出される。
マニューバ検出部１５でマニューバが検出された場合に
は、該マニューバ検出部１５は（１２）式を用い、マニ
ューバ検出時点（ｔ＝ｔd ）の推定情報行列を算出し、
そのうちの速度に関する成分をゼロとした内部状態量の
推定値及び該推定情報行列をそれぞれVXdo * 及びΣVX0*
-1とし、内部状態量推定部１６に送る。内部状態量推定
部１６は、それらマニューバ検出部１５から送られた推
定値VXdo * 及び推定情報行列ΣVX0* -1を先見情報とし
て、その推定値VXdo * を運動解析の初期値として用い、
時刻ｔ＝ｔd+1 ，ｔd+2 ，…，ｔm の（ｍ−ｄ）組の観
測方位角θm （t)及びドップラシフト観測周波数成分ν
km（t)（k ＝1,2,…,p）に対して、（１４）式の評価関
数を最小とする目標運動解析を行う。また、マニューバ
検出部１５でマニューバが検出されなかった場合、各時
刻ｔ＝ｔ1 ，ｔ2 ，…，ｔn で求めた内部状態量の推定
値VX0 * を運動解析の初期値として、ｔ＝ｔn+1 ，ｔ
n+2 ，…，ｔm （＞ｔn ）における（ｍ−ｎ）組の観測
方位角θm （t)及びドップラシフト観測周波数成分νkm
（t)（k ＝1 ,2, …,p）に対して、（５）式の評価関数
を最小とする目標運動解析を行う。
ｔn+2 ，…，ｔm （＞ｔn ) の（ｍ−ｎ）組の観測方位
角θm （t)及びドップラシフト観測周波数成分νkm（t)
（k=1,2,…,p) が初期値算出部１４及びマニューバ検出
部１５に入力されると、そのマニューバ検出部１５で
は、（１３）式で定義される推定算差の統計分布を調べ
る。その結果、例えばマニューバの発生が検出される。
マニューバ検出部１５でマニューバが検出された場合に
は、該マニューバ検出部１５は（１２）式を用い、マニ
ューバ検出時点（ｔ＝ｔd ）の推定情報行列を算出し、
そのうちの速度に関する成分をゼロとした内部状態量の
推定値及び該推定情報行列をそれぞれVXdo * 及びΣVX0*
-1とし、内部状態量推定部１６に送る。内部状態量推定
部１６は、それらマニューバ検出部１５から送られた推
定値VXdo * 及び推定情報行列ΣVX0* -1を先見情報とし
て、その推定値VXdo * を運動解析の初期値として用い、
時刻ｔ＝ｔd+1 ，ｔd+2 ，…，ｔm の（ｍ−ｄ）組の観
測方位角θm （t)及びドップラシフト観測周波数成分ν
km（t)（k ＝1,2,…,p）に対して、（１４）式の評価関
数を最小とする目標運動解析を行う。また、マニューバ
検出部１５でマニューバが検出されなかった場合、各時
刻ｔ＝ｔ1 ，ｔ2 ，…，ｔn で求めた内部状態量の推定
値VX0 * を運動解析の初期値として、ｔ＝ｔn+1 ，ｔ
n+2 ，…，ｔm （＞ｔn ）における（ｍ−ｎ）組の観測
方位角θm （t)及びドップラシフト観測周波数成分νkm
（t)（k ＝1 ,2, …,p）に対して、（５）式の評価関数
を最小とする目標運動解析を行う。
【００１０】
【発明が解決しようとする課題】しかしながら、従来の
目標運動解析方法では、次のような課題があった。内部
状態量に関する先見情報が全ての成分について得られて
いる場合を除いて、先見情報が未知の成分については適
当な値を仮定した状態量VXd0 * を、式(14)の評価関数を
最小化する内部状態量VX0 を求める際の初期値として使
わざるを得ず、この仮定した成分の影響により、解の収
束に時間がかかるという問題があった。本発明は、先見
情報が未知の状態量成分に適当な値を仮定せざるを得
ず、そのため解の収束に時間がかかるという問題点を解
決することを目的とする。
目標運動解析方法では、次のような課題があった。内部
状態量に関する先見情報が全ての成分について得られて
いる場合を除いて、先見情報が未知の成分については適
当な値を仮定した状態量VXd0 * を、式(14)の評価関数を
最小化する内部状態量VX0 を求める際の初期値として使
わざるを得ず、この仮定した成分の影響により、解の収
束に時間がかかるという問題があった。本発明は、先見
情報が未知の状態量成分に適当な値を仮定せざるを得
ず、そのため解の収束に時間がかかるという問題点を解
決することを目的とする。
【００１１】
【課題を解決するための手段】第１の発明は、前記課題
を解決するために、目標から放射された音を運動可能な
航走体に取り付けた受波器センサアレイで受信し、該目
標に対する方位角測定結果である方位観測量時系列か
ら、該目標が等速直線運動を行うものと仮定した上で該
目標の運動解析の初期値を定めて、該目標の位置及び速
度に関する状態量を推定する目標運動解析方法におい
て、次のような方法を講じている。即ち、第１の発明の
目標運動解析方法においては、前記目標の状態量に関す
る概略値とその精度からなる先見情報が得られている場
合、その先見情報と、前記方位観測量と該方位観測量の
測角精度と前記航走体の速度成分と該航走体の位置の時
系列とから、観測雑音成分が微小であるとした擬似線形
方程式を構成し、該擬似線形方程式を解くことで前記先
見情報の内容を含む解析初期値を前記運動解析の初期値
として定め、繰り返し演算を行なって非線形最小二乗法
を解くことで前記目標の位置及び速度に関する状態量を
推定するようにしている。
を解決するために、目標から放射された音を運動可能な
航走体に取り付けた受波器センサアレイで受信し、該目
標に対する方位角測定結果である方位観測量時系列か
ら、該目標が等速直線運動を行うものと仮定した上で該
目標の運動解析の初期値を定めて、該目標の位置及び速
度に関する状態量を推定する目標運動解析方法におい
て、次のような方法を講じている。即ち、第１の発明の
目標運動解析方法においては、前記目標の状態量に関す
る概略値とその精度からなる先見情報が得られている場
合、その先見情報と、前記方位観測量と該方位観測量の
測角精度と前記航走体の速度成分と該航走体の位置の時
系列とから、観測雑音成分が微小であるとした擬似線形
方程式を構成し、該擬似線形方程式を解くことで前記先
見情報の内容を含む解析初期値を前記運動解析の初期値
として定め、繰り返し演算を行なって非線形最小二乗法
を解くことで前記目標の位置及び速度に関する状態量を
推定するようにしている。
【００１２】第２の発明は、目標から放射された音を運
動可能な航走体に取り付けた受波器センサアレイで受信
し、該目標に対する方位角測定結果である方位観測量時
系列と該放射された音のドップラシフト周波数測定結果
である周波数観測量時系列とから、該目標の位置及び速
度に関する状態量と目標の音源固有周波数に関する状態
量とを、該目標が等速直線運動を行うものと仮定した上
で該目標の運動解析の初期値を定めて、推定する目標運
動解析方法において、次のような方法を講じている。即
ち、前記目標の状態量に関する概略値とその精度からな
る先見情報が得られている場合、その先見情報と、前記
方位観測量と該方位観測量の測角精度と前記航走体の速
度成分と該航走体の位置の時系列と、前記周波数観測量
と該周波数観測量の測定精度の時系列とから、観測雑音
成分が微小であるとした擬似線形方程式を構成し、該擬
似線形方程式を解くことで前記先見情報の内容を含む解
析初期値を前記運動解析の初期値として定め、繰り返し
演算を行なって非線形最小二乗法を解くことで前記目標
の位置及び速度に関する状態量と目標の音源固有周波数
に関する状態量とを推定するようにしている。
動可能な航走体に取り付けた受波器センサアレイで受信
し、該目標に対する方位角測定結果である方位観測量時
系列と該放射された音のドップラシフト周波数測定結果
である周波数観測量時系列とから、該目標の位置及び速
度に関する状態量と目標の音源固有周波数に関する状態
量とを、該目標が等速直線運動を行うものと仮定した上
で該目標の運動解析の初期値を定めて、推定する目標運
動解析方法において、次のような方法を講じている。即
ち、前記目標の状態量に関する概略値とその精度からな
る先見情報が得られている場合、その先見情報と、前記
方位観測量と該方位観測量の測角精度と前記航走体の速
度成分と該航走体の位置の時系列と、前記周波数観測量
と該周波数観測量の測定精度の時系列とから、観測雑音
成分が微小であるとした擬似線形方程式を構成し、該擬
似線形方程式を解くことで前記先見情報の内容を含む解
析初期値を前記運動解析の初期値として定め、繰り返し
演算を行なって非線形最小二乗法を解くことで前記目標
の位置及び速度に関する状態量と目標の音源固有周波数
に関する状態量とを推定するようにしている。
【００１３】第３の発明は、前記目標の変針及び変速等
のマニューバが検出された場合の第１の発明における前
記先見情報は、前記速度に関する成分をゼロとした該マ
ニューバ検出時点の前記目標の状態量の推定値とその状
態量の推定情報行列とで構成している。そして、前記構
成した先見情報と、前記マニューバ検出以降の前記方位
観測量と該方位観測量の測角精度と前記航走体の速度成
分と該航走体の位置の時系列とから観測雑音成分が微小
であるとした擬似線形方程式を構成し、その擬似線形方
程式を解くことで前記先見情報の内容を含む解析初期値
を前記運動解析の初期値として定め、繰り返し演算を行
なって非線形最小二乗法を解くことで前記目標の位置及
び速度に関する状態量に関する状態量を推定するように
している。
のマニューバが検出された場合の第１の発明における前
記先見情報は、前記速度に関する成分をゼロとした該マ
ニューバ検出時点の前記目標の状態量の推定値とその状
態量の推定情報行列とで構成している。そして、前記構
成した先見情報と、前記マニューバ検出以降の前記方位
観測量と該方位観測量の測角精度と前記航走体の速度成
分と該航走体の位置の時系列とから観測雑音成分が微小
であるとした擬似線形方程式を構成し、その擬似線形方
程式を解くことで前記先見情報の内容を含む解析初期値
を前記運動解析の初期値として定め、繰り返し演算を行
なって非線形最小二乗法を解くことで前記目標の位置及
び速度に関する状態量に関する状態量を推定するように
している。
【００１４】第４の発明は、前記目標の変針及び変速等
のマニューバが検出された場合の第２の発明における前
記先見情報は、前記速度に関する成分をゼロとした該マ
ニューバ検出時点の前記目標の状態量の推定値とその状
態量の推定情報行列とで構成している。そして、前記構
成した先見情報と、前記方位観測量と該方位観測量の測
角精度と前記航走体の速度成分と該航走体の位置の時系
列と、前記周波数観測量と該周波数観測量の測定精度の
時系列とから観測雑音成分が微小であるとした擬似線形
方程式を構成し、その擬似線形方程式を解くことで前記
先見情報の内容を含む解析初期値を前記運動解析の初期
値として定め、繰り返し演算を行なって非線形最小二乗
法を解くことで前記目標の位置及び速度に関する状態量
と目標の音源固有周波数に関する状態量とを推定するよ
うにしている。
のマニューバが検出された場合の第２の発明における前
記先見情報は、前記速度に関する成分をゼロとした該マ
ニューバ検出時点の前記目標の状態量の推定値とその状
態量の推定情報行列とで構成している。そして、前記構
成した先見情報と、前記方位観測量と該方位観測量の測
角精度と前記航走体の速度成分と該航走体の位置の時系
列と、前記周波数観測量と該周波数観測量の測定精度の
時系列とから観測雑音成分が微小であるとした擬似線形
方程式を構成し、その擬似線形方程式を解くことで前記
先見情報の内容を含む解析初期値を前記運動解析の初期
値として定め、繰り返し演算を行なって非線形最小二乗
法を解くことで前記目標の位置及び速度に関する状態量
と目標の音源固有周波数に関する状態量とを推定するよ
うにしている。
【００１５】第１から第４の発明によれば、以上のよう
に目標運動解析方法を構成したので、目標が放射した音
を受波器センサアレイで受信し、受信結果の観測量時系
列から、目標が等速直線運動を行うものと仮定した上で
目標の運動解析の初期値を定めて目標の状態量を推定す
る。ここで、目標の状態量に関する概略値とその精度か
らなる先見情報が得られている場合には、その先見情報
と、方位観測量と方位観測量の測角精度と航走体の速度
成分と航走体の位置の時系列等から、観測雑音成分が微
小であるとした擬似線形方程式が構成される。擬似線形
方程式を解くことで、先見情報の内容を含む解析初期値
が運動解析の初期値として定められる。繰り返し演算で
非線形最小二乗法を解くことで目標の状態量が推定され
る。この時の解析初期値は先見情報の内容を含むので、
非線形最小二乗法を解く上での初期値として適したもの
となる。従って、前記課題を解決できるのである。
に目標運動解析方法を構成したので、目標が放射した音
を受波器センサアレイで受信し、受信結果の観測量時系
列から、目標が等速直線運動を行うものと仮定した上で
目標の運動解析の初期値を定めて目標の状態量を推定す
る。ここで、目標の状態量に関する概略値とその精度か
らなる先見情報が得られている場合には、その先見情報
と、方位観測量と方位観測量の測角精度と航走体の速度
成分と航走体の位置の時系列等から、観測雑音成分が微
小であるとした擬似線形方程式が構成される。擬似線形
方程式を解くことで、先見情報の内容を含む解析初期値
が運動解析の初期値として定められる。繰り返し演算で
非線形最小二乗法を解くことで目標の状態量が推定され
る。この時の解析初期値は先見情報の内容を含むので、
非線形最小二乗法を解く上での初期値として適したもの
となる。従って、前記課題を解決できるのである。
【００１６】
【発明の実施の形態】本実施形態の特徴は、先見情報及
び観測量の時系列から疑似線形方程式を解き、先見情報
を含み観測量の時系列に最適な初期値を設定することで
あり、マニューバが発生した場合にも優れた目標運動解
析を行えるようになっている。まず、本実施形態の目標
運動解析方法の原理を説明する。（８）式の重み状列MW
（VX(t )) の未知要素を適当な推定値で近似した近似重
み行列MW* を用いて、その（８）式を変形すると（１
５）式が得られる。 ( MW* ) -1（Vb−MBVX0 ）＝VN ・・・（１５） （１５）式の両辺に(Σ-1/2 )Tを掛けると（１６）式が
得られる。ここで、Σ-1/2は、Σ-1のコレスキー分解に
より定まる上三角行列であり、次の（１７）式を満た
す。 (Σ-1/2 )T（MW* ）-1（Vb−MBVX0 ）＝(Σ-1/2 )T VN ・・・（１６） Σ＝( Σ1/2 ) T Σ1/2 ・・・（１７） （１６）をVN0 とおくと、 Ｅ［VN0 VN0 T ］＝ (Σ-1/2 )T Ｅ［VN VN T ］( Σ-1/2 ) ＝ (Σ-1/2 )T ΣΣ-1/2 ＝ (Σ-1/2 )T (Σ1 /2 )T Σ1 /2Σ-1/2 ＝Ｉ ・・・（１８） となる。
び観測量の時系列から疑似線形方程式を解き、先見情報
を含み観測量の時系列に最適な初期値を設定することで
あり、マニューバが発生した場合にも優れた目標運動解
析を行えるようになっている。まず、本実施形態の目標
運動解析方法の原理を説明する。（８）式の重み状列MW
（VX(t )) の未知要素を適当な推定値で近似した近似重
み行列MW* を用いて、その（８）式を変形すると（１
５）式が得られる。 ( MW* ) -1（Vb−MBVX0 ）＝VN ・・・（１５） （１５）式の両辺に(Σ-1/2 )Tを掛けると（１６）式が
得られる。ここで、Σ-1/2は、Σ-1のコレスキー分解に
より定まる上三角行列であり、次の（１７）式を満た
す。 (Σ-1/2 )T（MW* ）-1（Vb−MBVX0 ）＝(Σ-1/2 )T VN ・・・（１６） Σ＝( Σ1/2 ) T Σ1/2 ・・・（１７） （１６）をVN0 とおくと、 Ｅ［VN0 VN0 T ］＝ (Σ-1/2 )T Ｅ［VN VN T ］( Σ-1/2 ) ＝ (Σ-1/2 )T ΣΣ-1/2 ＝ (Σ-1/2 )T (Σ1 /2 )T Σ1 /2Σ-1/2 ＝Ｉ ・・・（１８） となる。
【００１７】よって、（８）式を行列MW* を用いてVX0
に関する擬似線形方程式を解くことは、次の（１９）式
で定義された評価関数Ｉ（VX0 ）を最小化する線形最小
二乗解を求めることと等価である。 Ｉ（VX0 ）＝‖(Σ-1/2 )T（MW* ）-1（vb−MBVX0 ）‖2 ・・・（１９） （１９）式を変形すると、次の（２０）式が獲られる。 Ｉ（VX0 ）＝（Vb−MBVX0 ）T （MW* ΣMW* T ）-1（Vb−MBVX0 ） ・・・（２０） 従って、先見情報として、内部状態量の推定状態量VXd0
* 及び該推定状態量の推定情報行列ΣVXd0* -1が知られ
ている場合には、 Ｇ（VX0 ）＝（Vb−MBVX0 ）T （MW*ΣMW*T ）-1 （Vb−MBVX0 ） ＋（VX0 −VXd0 * ）ΣVXd0* -1（VX0 −VXd0 * ）・・・（２１） なる評価関数を最小とする解を求めることで、先見情報
を含み観測量の時系列に適した初期値を得ることができ
る。（２１）式をVX0 で微分しておけば、次の（２２）
の線形方程式が得られ、目標運動解析の解析初期値VX00
は、（２３）式となる。 ｛MBT （MW*ΣMW*T ）-1MB＋ΣVXd0* -1｝VX0 ＝MBT （MW*ΣMW*T ）-1Vb ＋ΣVXd0* -1VXd0 * ・・・（２２） VX00 * ＝｛MBT （MW*ΣMW*T ）-1MB＋ΣVXd0* -1｝-1 ×｛MBT （MWΣMWT ）-1Vb＋ΣVXd0* -1VXd0 * ｝・・・（２３） 今、時刻ｔ＝ｔd においてマニューバが検出された場合
を想定すると、時刻ｔ＝ｔd における推定状態量及び該
状態量の推定情報行列において、速度に関する成分をゼ
ロと置いたものをVXd0 * 及び行列ΣVXd0* -1とし、これ
らを先見情報とする。そして、各時刻ｔ＝ｔd+1 ，ｔ
d+2 ，…，ｔm の（ｍ−ｄ）組の観測方位θm （t)及び
ドップラシフト観測周波数或分νkm（k=1,2,…,p) とか
ら、（２３）式を用いてマニューバ後の目標運動解析の
解析初期値VX00 * を算出する。初期値VX00 * の算出にお
いては、推定情報行列ΣVXd0* -1の速度に関する成分を
ゼロとしているので、推定状態量VXd0 * のうちゼロとし
た仮定した速度成分の影響を除くことができる。また、
該初期値VX00 * はマニューバ検出時の位置に関する先見
情報を含み、かつマニューバ検出以降の観測時系列に対
する疑似線形方程式の解にもなっている。したがって、
事前情報と観測データの両方の情報を活かした初期値と
して、（１４）式の評価関数を最小化する最小二乗解を
得る際に、収束を加速できるのである。
に関する擬似線形方程式を解くことは、次の（１９）式
で定義された評価関数Ｉ（VX0 ）を最小化する線形最小
二乗解を求めることと等価である。 Ｉ（VX0 ）＝‖(Σ-1/2 )T（MW* ）-1（vb−MBVX0 ）‖2 ・・・（１９） （１９）式を変形すると、次の（２０）式が獲られる。 Ｉ（VX0 ）＝（Vb−MBVX0 ）T （MW* ΣMW* T ）-1（Vb−MBVX0 ） ・・・（２０） 従って、先見情報として、内部状態量の推定状態量VXd0
* 及び該推定状態量の推定情報行列ΣVXd0* -1が知られ
ている場合には、 Ｇ（VX0 ）＝（Vb−MBVX0 ）T （MW*ΣMW*T ）-1 （Vb−MBVX0 ） ＋（VX0 −VXd0 * ）ΣVXd0* -1（VX0 −VXd0 * ）・・・（２１） なる評価関数を最小とする解を求めることで、先見情報
を含み観測量の時系列に適した初期値を得ることができ
る。（２１）式をVX0 で微分しておけば、次の（２２）
の線形方程式が得られ、目標運動解析の解析初期値VX00
は、（２３）式となる。 ｛MBT （MW*ΣMW*T ）-1MB＋ΣVXd0* -1｝VX0 ＝MBT （MW*ΣMW*T ）-1Vb ＋ΣVXd0* -1VXd0 * ・・・（２２） VX00 * ＝｛MBT （MW*ΣMW*T ）-1MB＋ΣVXd0* -1｝-1 ×｛MBT （MWΣMWT ）-1Vb＋ΣVXd0* -1VXd0 * ｝・・・（２３） 今、時刻ｔ＝ｔd においてマニューバが検出された場合
を想定すると、時刻ｔ＝ｔd における推定状態量及び該
状態量の推定情報行列において、速度に関する成分をゼ
ロと置いたものをVXd0 * 及び行列ΣVXd0* -1とし、これ
らを先見情報とする。そして、各時刻ｔ＝ｔd+1 ，ｔ
d+2 ，…，ｔm の（ｍ−ｄ）組の観測方位θm （t)及び
ドップラシフト観測周波数或分νkm（k=1,2,…,p) とか
ら、（２３）式を用いてマニューバ後の目標運動解析の
解析初期値VX00 * を算出する。初期値VX00 * の算出にお
いては、推定情報行列ΣVXd0* -1の速度に関する成分を
ゼロとしているので、推定状態量VXd0 * のうちゼロとし
た仮定した速度成分の影響を除くことができる。また、
該初期値VX00 * はマニューバ検出時の位置に関する先見
情報を含み、かつマニューバ検出以降の観測時系列に対
する疑似線形方程式の解にもなっている。したがって、
事前情報と観測データの両方の情報を活かした初期値と
して、（１４）式の評価関数を最小化する最小二乗解を
得る際に、収束を加速できるのである。
【００１８】図１は、本発明の実施形態の目標運動解析
方法に用いられる目標運動解析装置の機能ブロック図で
ある。この目標運動解析装置は、入力端子２１を介し
て、図示しない受波器センサアレイからの受信信号を入
力する方位情報算出部２２及び周波数情報算出部２３を
備えている。方位情報算出部２２は、目標２の放射音の
到来方位を受信信号から算出するものであり、周波数情
報算出部２３は、目標の放射音の周波数成分の算出を行
うものである。方位情報算出部２２と周波数情報算出部
３３の出力側は、共に初期値算出部２４とマニューバ検
出部２５に接続されている。マニューバ検出部２５は、
方位情報算出部２２及び周波数情報算出部２３から与え
られた方位情報と周波数情報から、目標２のマニューバ
の検出を行うと共に、推定状態量VXd0 * 及び推定情報行
列ΣVXd0* -1を出力する機能を有し、それらを初期値算
出部２４に渡すようなっている。初期値算出部２４は、
方位情報と周波数情報と推定状態量VXd0 * 及び推定情報
行列ΣVXd0* -1とから、初期状態量を算出する機能を有
している。初期値算出部２４とマニューバ検出部２５に
は、内部状態量推定部２６が接続されている。内部状態
量推定部２６は目標の運動を解析するものであり、その
解析結果は、出力端子２７を介して出力されるようにな
っている。
方法に用いられる目標運動解析装置の機能ブロック図で
ある。この目標運動解析装置は、入力端子２１を介し
て、図示しない受波器センサアレイからの受信信号を入
力する方位情報算出部２２及び周波数情報算出部２３を
備えている。方位情報算出部２２は、目標２の放射音の
到来方位を受信信号から算出するものであり、周波数情
報算出部２３は、目標の放射音の周波数成分の算出を行
うものである。方位情報算出部２２と周波数情報算出部
３３の出力側は、共に初期値算出部２４とマニューバ検
出部２５に接続されている。マニューバ検出部２５は、
方位情報算出部２２及び周波数情報算出部２３から与え
られた方位情報と周波数情報から、目標２のマニューバ
の検出を行うと共に、推定状態量VXd0 * 及び推定情報行
列ΣVXd0* -1を出力する機能を有し、それらを初期値算
出部２４に渡すようなっている。初期値算出部２４は、
方位情報と周波数情報と推定状態量VXd0 * 及び推定情報
行列ΣVXd0* -1とから、初期状態量を算出する機能を有
している。初期値算出部２４とマニューバ検出部２５に
は、内部状態量推定部２６が接続されている。内部状態
量推定部２６は目標の運動を解析するものであり、その
解析結果は、出力端子２７を介して出力されるようにな
っている。
【００１９】次に、図１の動作を説明するが、ここで
は、マニューバ検出時の目標２の位置情報を先見情報と
した場合を説明する。目標２の放射した音が、図示しな
い受波器センサアレイで受信され、その受信信号が入力
端子２１を介して方位情報算出部２２と周波数情報算出
部２３に送られる。方位情報算出部２２は受信信号から
観測方位角θm (t) を算出し、周波数情報算出部２３は
ドップラシフト観測周波数成分νkm(t)(k=1,2,…,p) を
算出し、それらの結果が初期値算出部２４に送られる。
各時刻ｔ＝ｔ1 ，ｔ2 ，…，ｔn におけるｎ組の観測方
位角θm （t)及びドップラシフト観測周波数成分ν
km(t)(k=1,2,…,p) が、初期値算出部２４に入力される
と、該初期値算出部２４は（１０）式を用いて状態量の
解析初期値VX00 * を求め、その解析初期値VX00 * を内部
状態量推定部２６に送る。内部状態量推定部２６は、解
析初期値VX00 * を運動解析の初期値とし、（５）の評価
関数を最小化する目標運動解析を行って内部状態量の推
定値VX0 * を求める。内部状態量推定部２６は、その推
定値VX0 * をマニューバ検出部２５に送るともに、
（４）式及び（２）式を用いて、目標２の内部状態量Vx
2 （t)の推定値VX2 * （t)を出力端子２７を介して出力
する。
は、マニューバ検出時の目標２の位置情報を先見情報と
した場合を説明する。目標２の放射した音が、図示しな
い受波器センサアレイで受信され、その受信信号が入力
端子２１を介して方位情報算出部２２と周波数情報算出
部２３に送られる。方位情報算出部２２は受信信号から
観測方位角θm (t) を算出し、周波数情報算出部２３は
ドップラシフト観測周波数成分νkm(t)(k=1,2,…,p) を
算出し、それらの結果が初期値算出部２４に送られる。
各時刻ｔ＝ｔ1 ，ｔ2 ，…，ｔn におけるｎ組の観測方
位角θm （t)及びドップラシフト観測周波数成分ν
km(t)(k=1,2,…,p) が、初期値算出部２４に入力される
と、該初期値算出部２４は（１０）式を用いて状態量の
解析初期値VX00 * を求め、その解析初期値VX00 * を内部
状態量推定部２６に送る。内部状態量推定部２６は、解
析初期値VX00 * を運動解析の初期値とし、（５）の評価
関数を最小化する目標運動解析を行って内部状態量の推
定値VX0 * を求める。内部状態量推定部２６は、その推
定値VX0 * をマニューバ検出部２５に送るともに、
（４）式及び（２）式を用いて、目標２の内部状態量Vx
2 （t)の推定値VX2 * （t)を出力端子２７を介して出力
する。
【００２０】次ステップにおいて、各時刻ｔ＝ｔn+1 ，
ｔn+2 ，…，ｔm (>ｔn ) の（ｍ−ｎ）組の観測方位角
θm （t)及びドップラシフト観測周波数成分νkm(t)(k=
1,2,…,p) が、初期値算出部２４及びマニューバ検出部
２５に入力されると、マニューバ検出部２５では、（１
３）式で定義される推定算差の統計分布を調査し、その
調査結果でマニューバの検出を行う。マニューバ検出部
２５でマニューバが検出された場合、そのマニューバ検
出部２５は、（１２）式を用いてマニューバ検出時点ｔ
＝ｔd の推定情報行列を算出し、速度に関する成分をゼ
ロとおいた内部状態量の推定値及び該推定情報行列をそ
れぞれVXd0 * 及びΣVXd0* -1とし、それらを初期値算出
部２４に送る。初期値算出部２４は、内部状態量の推定
値VXd0 * 及び該推定情報行列ΣVXd0* -1を先見情報と
し、各時刻ｔ＝ｔd+1 ，ｔd+2 ，…，ｔm の（ｍ−ｄ）
組の観測方位角θm （t)及びドップラシフト観測周波数
成分νkm（t)（ｋ＝１,2, …,p）に対して（２３）式を
用い運動解析の解析初期値VX00 * を求め、その解析初期
値VX00 * を内部状態量推定部２６に送る。内部状態量推
定部２６は解析初期値VX00 * を運動解析の初期値とし、
（１４）式の評価関数を最少化する目標運動解析を行
う。そして、内部状態量推定部２６は（４）式及び
（２）式を用い、目標２に対する位置ベクトルVx2 (t)
に対応する内部状態量の推定値VX2 * (t) を算出してそ
れを出力端子２７を介して出力する。
ｔn+2 ，…，ｔm (>ｔn ) の（ｍ−ｎ）組の観測方位角
θm （t)及びドップラシフト観測周波数成分νkm(t)(k=
1,2,…,p) が、初期値算出部２４及びマニューバ検出部
２５に入力されると、マニューバ検出部２５では、（１
３）式で定義される推定算差の統計分布を調査し、その
調査結果でマニューバの検出を行う。マニューバ検出部
２５でマニューバが検出された場合、そのマニューバ検
出部２５は、（１２）式を用いてマニューバ検出時点ｔ
＝ｔd の推定情報行列を算出し、速度に関する成分をゼ
ロとおいた内部状態量の推定値及び該推定情報行列をそ
れぞれVXd0 * 及びΣVXd0* -1とし、それらを初期値算出
部２４に送る。初期値算出部２４は、内部状態量の推定
値VXd0 * 及び該推定情報行列ΣVXd0* -1を先見情報と
し、各時刻ｔ＝ｔd+1 ，ｔd+2 ，…，ｔm の（ｍ−ｄ）
組の観測方位角θm （t)及びドップラシフト観測周波数
成分νkm（t)（ｋ＝１,2, …,p）に対して（２３）式を
用い運動解析の解析初期値VX00 * を求め、その解析初期
値VX00 * を内部状態量推定部２６に送る。内部状態量推
定部２６は解析初期値VX00 * を運動解析の初期値とし、
（１４）式の評価関数を最少化する目標運動解析を行
う。そして、内部状態量推定部２６は（４）式及び
（２）式を用い、目標２に対する位置ベクトルVx2 (t)
に対応する内部状態量の推定値VX2 * (t) を算出してそ
れを出力端子２７を介して出力する。
【００２１】また、マニューバ検出部２５においてマニ
ューバが検出されなかった場合、内部状態量の推定値VX
0 * を運動解析の初期値として、各時刻ｔ＝ｔ1 ，
ｔ2 ，…，ｔm のｍ組の観測方位角θm （t)及びドップ
ラシフト観測周波数成分νkm（t)（k ＝１,2, …,p）に
対して（５）式の評価関数を最小とする目標運動解析を
行ない、目標２に対する位置ベクトルVx2 (t) に対応す
る内部状態量の推定値VX2 * (t) を算出して出力端子２
７を介して出力する。以上のように、本実施形態によれ
ば、目標２の状態量に関する先見情報と、観測量及びそ
の観測量の測定精度から、観測雑音成分が微小であると
して擬似線形方程式を構成し、その擬似線形方程式を解
くことで、目標２の状態量に関する先見情報を含んだ運
動解析の解析初期値VX00 * を定めている。そのため、事
前情報と観測データの両方の情報を生かすことが可能に
なり、（１４）式の非線形方程式を解く際に、解の収束
が早くなり、従来の課題を解決できる。
ューバが検出されなかった場合、内部状態量の推定値VX
0 * を運動解析の初期値として、各時刻ｔ＝ｔ1 ，
ｔ2 ，…，ｔm のｍ組の観測方位角θm （t)及びドップ
ラシフト観測周波数成分νkm（t)（k ＝１,2, …,p）に
対して（５）式の評価関数を最小とする目標運動解析を
行ない、目標２に対する位置ベクトルVx2 (t) に対応す
る内部状態量の推定値VX2 * (t) を算出して出力端子２
７を介して出力する。以上のように、本実施形態によれ
ば、目標２の状態量に関する先見情報と、観測量及びそ
の観測量の測定精度から、観測雑音成分が微小であると
して擬似線形方程式を構成し、その擬似線形方程式を解
くことで、目標２の状態量に関する先見情報を含んだ運
動解析の解析初期値VX00 * を定めている。そのため、事
前情報と観測データの両方の情報を生かすことが可能に
なり、（１４）式の非線形方程式を解く際に、解の収束
が早くなり、従来の課題を解決できる。
【００２２】なお、本発明は、上記実施形態に限定され
ず種々の変形が可能である。例えば、上記実施形態で
は、求める内部状態量を、目標２の位置と速度と目標２
の固有音源周波数としているが、目標２の位置と速度だ
けを求める場合にも適用が可能である。この場合には、
各（３），（４），（５，）（９）式から音源固有周波
数及びドップラシフト周波数を取り除けばよい。また、
図１の目標運動解析装置では、方位情報算出部２２、周
波数情報算出部２３、初期値算出部２４、マニューバ検
出部２５、及び内部状態量推定部２６を個別なものとし
て記載しているが、受波器センサアレイから受信信号に
対してオンライン処理で上記実施形態と同様の処理をす
る構成としてもよい。
ず種々の変形が可能である。例えば、上記実施形態で
は、求める内部状態量を、目標２の位置と速度と目標２
の固有音源周波数としているが、目標２の位置と速度だ
けを求める場合にも適用が可能である。この場合には、
各（３），（４），（５，）（９）式から音源固有周波
数及びドップラシフト周波数を取り除けばよい。また、
図１の目標運動解析装置では、方位情報算出部２２、周
波数情報算出部２３、初期値算出部２４、マニューバ検
出部２５、及び内部状態量推定部２６を個別なものとし
て記載しているが、受波器センサアレイから受信信号に
対してオンライン処理で上記実施形態と同様の処理をす
る構成としてもよい。
【００２３】
【発明の効果】以上詳細に説明したように、第１から第
４の発明によれば、目標の状態量に関する先見情報と、
観測量及びその観測量の測定精度から、観測雑音成分が
微小であるとして擬似線形方程式を構成し、その擬似線
形方程式を解くことで、目標の状態量に関する先見情報
を含んだ運動解析の解析初期値を定めている。そのた
め、事前情報と観測データの両方の情報を生かすことが
可能になり、その解析初期値を用いて非線形方程式を解
けば解の収束が速くなり、目標の状態量の推定値が速く
得られる。
４の発明によれば、目標の状態量に関する先見情報と、
観測量及びその観測量の測定精度から、観測雑音成分が
微小であるとして擬似線形方程式を構成し、その擬似線
形方程式を解くことで、目標の状態量に関する先見情報
を含んだ運動解析の解析初期値を定めている。そのた
め、事前情報と観測データの両方の情報を生かすことが
可能になり、その解析初期値を用いて非線形方程式を解
けば解の収束が速くなり、目標の状態量の推定値が速く
得られる。
【図１】本発明の実施形態の目標運動解析方法に用いら
れる目標運動解析装置の機能ブロック図である。
れる目標運動解析装置の機能ブロック図である。
【図２】従来の目標運動解析方法を示す幾何学的説明図
である。
である。
【図３】従来の目標運動解析方法に用いられる目標運動
解析装置の機能ブロック図である。
解析装置の機能ブロック図である。
２２ 方位情報算出部 ２３ 周波数情報算出部 ２４ 初期位置算出部 ２５ マニューバ検出部 ２６ 内部状態量推定部
Claims (4)
- 【請求項１】 目標から放射された音を運動可能な航走
体に取り付けた受波器センサアレイで受信し、該目標に
対する方位角測定結果である方位観測量時系列から、該
目標が等速直線運動を行うものと仮定した上で該目標の
運動解析の初期値を定めて、該目標の位置及び速度に関
する状態量を推定する目標運動解析方法において、 前記目標の状態量に関する概略値とその精度からなる先
見情報が得られている場合、その先見情報と、前記方位
観測量と該方位観測量の測角精度と前記航走体の速度成
分と該航走体の位置の時系列とから、観測雑音成分が微
小であるとした擬似線形方程式を構成し、該擬似線形方
程式を解くことで前記先見情報の内容を含む解析初期値
を前記運動解析の初期値として定め、繰り返し演算を行
なって非線形最小二乗法を解くことまたは該初期値をオ
ンライン処理に用いて前記目標の位置及び速度に関する
状態量を推定することを特徴とする目標運動解析方法。 - 【請求項２】 目標から放射された音を運動可能な航走
体に取り付けた受波器センサアレイで受信し、該目標に
対する方位角測定結果である方位観測量時系列と該放射
された音のドップラシフト周波数測定結果である周波数
観測量時系列とから、該目標の位置及び速度に関する状
態量と目標の音源固有周波数に関する状態量とを、該目
標が等速直線運動を行うものと仮定した上で該目標の運
動解析の初期値を定めて推定する目標運動解析方法にお
いて、 前記目標の状態量に関する概略値とその精度からなる先
見情報が得られている場合、その先見情報と、前記方位
観測量と該方位観測量の測角精度と前記航走体の速度成
分と該航走体の位置の時系列と、前記周波数観測量と該
周波数観測量の測定精度の時系列とから観測雑音成分が
微小であるとした擬似線形方程式を構成し、該擬似線形
方程式を解くことで前記先見情報の内容を含む解析初期
値を前記運動解析の初期値として定め、繰り返し演算を
行なって非線形最小二乗法を解くことまたは該初期値を
オンライン処理に用いて、前記目標の位置及び速度に関
する状態量と目標の音源固有周波数に関する状態量とを
推定することを特徴とする目標運動解析方法。 - 【請求項３】 前記目標の変針及び変速等のマニューバ
が検出された場合の前記先見情報は、前記速度に関する
成分をゼロとした該マニューバ検出時点の前記目標の状
態量の推定値とその状態量の推定情報行列とで構成し、 前記構成した先見情報と、前記マニューバ検出以降の前
記方位観測量と該方位観測量の測角精度と前記航走体の
速度成分と該航走体の位置の時系列とから、観測雑音成
分が微小であるとした擬似線形方程式を構成し、その擬
似線形方程式を解くことで前記先見情報の内容を含む解
析初期値を前記運動解析の初期値として定め、繰り返し
演算を行なって非線形最小二乗法を解くことまたは該初
期値をオンライン処理に用いて前記目標の位置及び速度
に関する状態量に関する状態量を推定することを特徴と
する請求項１記載の目標運動解析方法。 - 【請求項４】 前記目標の変針及び変速等のマニューバ
が検出された場合の前記先見情報は、前記速度に関する
成分をゼロとした該マニューバ検出時点の前記目標の状
態量の推定値とその状態量の推定情報行列とで構成し、 前記構成した先見情報と、前記方位観測量と該方位観測
量の測角精度と前記航走体の速度成分と該航走体の位置
の時系列と、前記周波数観測量と該周波数観測量の測定
精度の時系列とから観測雑音成分が微小であるとした擬
似線形方程式を構成し、その擬似線形方程式を解くこと
で前記先見情報の内容を含む解析初期値を前記運動解析
の初期値として定め、繰り返し演算を行なって非線形最
小二乗法を解くことまたは該初期値をオンライン処理に
用いて前記目標の位置及び速度に関する状態量と目標の
音源固有周波数に関する状態量とを推定することを特徴
とする請求項２記載の目標運動解析方法。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP29418495A JPH09133750A (ja) | 1995-11-13 | 1995-11-13 | 目標運動解析方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP29418495A JPH09133750A (ja) | 1995-11-13 | 1995-11-13 | 目標運動解析方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| JPH09133750A true JPH09133750A (ja) | 1997-05-20 |
Family
ID=17804407
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP29418495A Withdrawn JPH09133750A (ja) | 1995-11-13 | 1995-11-13 | 目標運動解析方法 |
Country Status (1)
| Country | Link |
|---|---|
| JP (1) | JPH09133750A (ja) |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2008268044A (ja) * | 2007-04-23 | 2008-11-06 | Oki Electric Ind Co Ltd | 目標状態量推定方法 |
| JP2009058367A (ja) * | 2007-08-31 | 2009-03-19 | Hitachi Ltd | 目標運動解析方法及び装置 |
-
1995
- 1995-11-13 JP JP29418495A patent/JPH09133750A/ja not_active Withdrawn
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2008268044A (ja) * | 2007-04-23 | 2008-11-06 | Oki Electric Ind Co Ltd | 目標状態量推定方法 |
| JP2009058367A (ja) * | 2007-08-31 | 2009-03-19 | Hitachi Ltd | 目標運動解析方法及び装置 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110749891B (zh) | 一种可估计未知有效声速的自适应水下单信标定位方法 | |
| CN104596514A (zh) | 加速度计和陀螺仪的实时降噪系统和方法 | |
| US5877998A (en) | Recursive method for target motion analysis | |
| Agate et al. | Road-constrained target tracking and identification using a particle filter | |
| JP3095189B2 (ja) | ナビゲーション装置 | |
| Konatowski et al. | A comparison of estimation accuracy by the use of KF, EKF & UKF filters | |
| CN119644321B (zh) | 一种基于椭圆轨迹约束目标建模与跟踪方法 | |
| CN114459484A (zh) | 航迹初始化方法及装置 | |
| Tamas et al. | State estimation based on Kalman filtering techniques in navigation | |
| JP6757227B2 (ja) | 運動パラメータ推定装置、運動パラメータ推定方法及びプログラム | |
| JP3097471B2 (ja) | 追尾装置 | |
| JP5012168B2 (ja) | 目標状態量推定方法 | |
| CN110441748A (zh) | 一种基于幅度信息的α-β滤波方法 | |
| Levy et al. | Unscented Kalman Filter with a Nonlinear Propagation Model for Navigation Applications | |
| JP2739054B2 (ja) | 目標運動解析方法 | |
| JPH0797136B2 (ja) | 多目標追尾方法及びその装置 | |
| JPH1062508A (ja) | 目標運動解析方法 | |
| JP3142489B2 (ja) | 目標運動解析方法及び目標運動解析装置 | |
| Zhang | An adaptive unscented kalman filter for dead reckoning systems | |
| JP3154666B2 (ja) | 目標追尾装置 | |
| JPH10319108A (ja) | センサ姿勢及び位置のバイアス誤差推定装置 | |
| JP3129388B2 (ja) | 目標運動解析方法 | |
| JP3302870B2 (ja) | センサのバイアス誤差推定装置 | |
| JP3218181B2 (ja) | センサのバイアス誤差推定装置 | |
| CN120275950B (zh) | 一种天基雷达光学复合空间目标距离速度测量方法 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A300 | Withdrawal of application because of no request for examination |
Free format text: JAPANESE INTERMEDIATE CODE: A300 Effective date: 20030204 |