(57)【特許請求の範囲】
【請求項１】
センサによって検知された目標移動体の観測データを用いて前記目標移動体の距離、方位、速力および針路を含む運動パラメータを解析する目標運動解析方法であって、距離、速力、針路のうちの少なくとも１つの運動パラメータについてオペレータ入力を受け付け、前記オペレータ入力された運動パラメータを除く次元縮退した状態ベクトルを設定し、設定された前記状態ベクトルを用いるとともに前記オペレータ入力された運動パラメータを真として解析処理を行い、得られた解析結果を表示装置上に表示することを特徴とするオペレータ介在型目標運動解析方法。
【請求項２】
センサによって検知された目標移動体の観測データを用いて前記目標移動体の距離、方位、速力および針路を含む運動パラメータを解析する目標運動解析装置であって、距離、速力、針路のうちの少なくとも１つの運動パラメータについてオペレータ入力を受け付ける手段と、前記オペレータ入力された運動パラメータを除く次元縮退した状態ベクトルを設定する手段と、設定された前記状態ベクトルを用いるとともに前記オペレータ入力された運動パラメータを真として解析処理をする手段と、得られた解析結果を表示装置上に表示する手段とを有することを特徴とするオペレータ介在型目標運動解析装置。
【請求項３】
前記表示する手段は、さらに前記目標移動体の逐次解析時刻における解析位置履歴データと、最新解析時刻における前記目標移動体の解析結果に基づき等速直進性を仮定して過去に遡った位置計算履歴データとを並列して表示する手段を含むことを特徴とする請求項２記載のオペレータ介在型目標運動解析装置。
【請求項４】
前記表示する手段は、さらに前記目標移動体の逐次解析時刻における解析方位履歴データと観測方位履歴データとを並列して表示する手段を含むことを特徴とする請求項２記載のオペレータ介在型目標運動解析装置。
【請求項５】
前記表示装置において、前記目標移動体の逐次時刻における解析方位履歴データ又は解析位置履歴データのうち、異常データと判断される解析方位履歴データ又は解析位置履歴データを指定することにより、該当する解析時刻の観測データを排除して、再度解析処理を実施する手段を備えることを特徴とする請求項２記載のオペレータ介在型目標運動解析装置。
【請求項６】
前記表示装置において、前記目標移動体の運動パラメータの解析結果が、オペレータ入力を反映して解析されたものである場合に、それをオペレータに告知するための告知手段を備えることを特徴とする請求項２記載のオペレータ介在型目標運動解析装置。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
本発明は、船舶、航空機等、目標移動体の運動諸元を解析する技術に関する。
【０００２】
【従来の技術】
目標移動体から放射される音波又は電波をパッシブセンサにより時系列的に測定した観測データを用いて、目標移動体の運動諸元を解析する技術が知られている。
【０００３】
図１は、目標運動解析の原理を示す図である。観測船１０１は、パッシブセンサを使って目標１０２が放射する音波又は電波を受信することにより、到来音波又は電波の方位~By及び周波数~fを観測する。パッシブセンサから得られる観測情報は次のとおりである：
~By = By + σBy
~f  = f  + σf
ここで、~By：観測方位、By：方位の真値、σBy：方位誤差、~f：観測周波数、f：観測周波数の真値、σf：周波数誤差である。
【０００４】
観測船１０１の目標運動解析装置は、観測された方位及び周波数、そして観測時刻及び観測位置の時系列データに基づいて、目標１０２の状態ベクトル^X（距離r、方位By、速力Mt、針路Ct、放射音波又は電波の周波数f0）を解析する。その解析手法は、状態ベクトル^Xから推定される各観測時刻の解析方位By及び解析周波数fと、観測方位~By及び周波数~fとの残差２乗を観測誤差分散σ~By2及びσ~f2で正規化し、これらを時系列的に合計した評価関数Ф(^Xn)が最小となるような状態ベクトル^Xを求める技術である。
【０００５】
このような非線形最適化問題を解いて目標運動諸元を解析する目標運動解析技術に関連する従来技術として、特開平９−１３３７５０号公報「目標運動解析方法」（特許文献１）又は特許第３１２９３８８号「目標運動解析方法」（特許文献２）に記載された技術が知られている。
【０００６】
特許文献１は、目標運動諸元の概略値とその精度からなるオペレータ先見情報が入力されている場合に、推定状態量とこの推定状態量の推定情報行列を加えて擬似線形方程式を構成し、この擬似線形方程式を解いて目標運動諸元の初期値を解析する技術について記載する。
【０００７】
特許文献２は、目標速力が指定値以下である拘束条件の下で擬似線形方程式を解いて目標運動諸元の初期値を解析する技術について記載する。
【０００８】
【特許文献１】
特開平９−１３３７５０号公報
【特許文献２】
特許第３１２９３８８号
【０００９】
【発明が解決しようとする課題】
特許文献１は、目標運動諸元の概略値とその精度に関するオペレータ先見情報を必要とする。しかし通常、この概略値の精度に関するオペレータ先見情報が得られることはないため、適当な精度を仮定せざるをえない。従って観測データの測定誤差の状況により概略値の精度に関する仮定値が適切でない場合、目標運動諸元の解析値がオペレータ先見情報として入力した概略値に収束しない、あるいは収束時間が遅延し、オペレータ先見情報が目標運動解析に反映されない、あるいは反映されるのに時間がかかるという問題がある。また概略値が誤った情報である場合に、オペレータ先見情報が誤りであることを認知判断する方法が提供されていないという問題がある。また特許文献２は、目標速力下限値および目標距離上下限値に関するオペレータ先見情報を活用した解析について考慮されていない。
【００１０】
本発明の目的は、前述した従来技術の問題を解決し、オペレータ先見情報を反映して目標運動解析の収束時間を短縮することにある。
【００１１】
【課題を解決するための手段】
本発明によれば、前記目的は、センサによって検知された目標移動体の観測データを用いて目標移動体の運動パラメータを解析する目標運動解析方法において、距離、速力、針路の運動パラメータのうち少なくとも１つの運動パラメータについてオペレータ入力を受け付け、オペレータ入力された運動パラメータを除く次元縮退した状態ベクトルを設定し、設定された状態ベクトルを用いて解析する手段をもつことにより達成される。すなわち、オペレータ入力の精度を仮定しなくとも、オペレータ入力された運動パラメータ以外の残りの未知の運動パラメータに関する解析を実施することにより収束時間の短縮を可能とする。
【００１２】
また本発明によれば、目標移動体に関する距離上下限値及び速力上下限値のオペレータ入力を受け付ける手段を設け、ラグランジュ乗数法及びフィルタゲインの調節を行う処理を設けることにより、距離上下限値及び速力上下限値の制約条件の下で解析を実施し収束時間の短縮を可能とする。
【００１３】
また本発明によれば、解析結果を表示する表示装置において、目標移動体の逐次解析時刻における解析位置履歴データと、現在時刻における目標移動体の解析結果に基づき等速直進性を仮定して過去に遡った位置計算履歴データとの比較表示、および目標移動体の逐次解析時刻における解析方位履歴データと、観測データベース記憶手段における観測方位履歴データとの比較表示を実施することにより、オペレータ入力が誤りであること、または異常データが混在していることを認知判断することができ、オペレータ入力を再度実施すること、またはオペレータが異常データを排除することにより収束時間の短縮を可能とする。
【００１４】
さらに本発明によれば、解析結果を表示する表示装置において、オペレータ入力無しの場合の解析位置履歴データまたは解析方位履歴データに加えてオペレータ入力有りの場合の解析位置履歴データ又は解析方位履歴データを表示することにより、オペレータ入力が誤りであることを認知判断し、オペレータ入力を再度実施することにより収束時間の短縮を可能とする。
【００１５】
【発明の実施の形態】
以下、本発明の実施形態について図面を用いて詳細に説明する。なお以下に説明する本発明の実施形態は、対象目標を船舶として音響センサにより目標音波の到来方位を時系列的に観測し、オペレータ先見情報を反映して目標運動諸元を解析する場合について説明する。しかし本発明は、船舶のみならず航空機等の音波又は電波を放射しながら移動する目標の音波又は電波の到来方位と受信周波数を時系列的に観測する場合にも適用可能である。
【００１６】
図２は、本実施形態の目標運動解析装置２１３の構成を示す図である。目標運動解析装置２１３には、パッシブセンサ２０１、航海センサ２０２及び位置演算器２０３が接続される。パッシブセンサ２０１は、目標音源からの音波の到来方位及び受信周波数を受信する。航海センサ２０２は、ログ及びジャイロを有し、目標の位置を計測する。位置演算器２０３は、観測船の位置を計測する。
【００１７】
目標運動解析装置２１３は、ＣＰＵ、メモリ、記憶装置、入力装置２０８及び表示装置２０６を有する専用計算機である。観測データベース２０４は、記憶装置上に格納され、パッシブセンサ２０１、航海センサ２０２及び位置演算器２０３から入力される観測データを時系列に蓄積する。解析結果データベース２１０は、記憶装置上に格納され、解析結果のデータを格納する。入力装置２０８は、目標運動諸元に関するオペレータ先見情報を入力し、また観測データベース２０４上の観測データの削除のための入力などを行う。表示装置２０６は、解析結果のデータを種々の形式で表示する。
【００１８】
解析処理部２０５、状態ベクトル設定部２０７、時系列残差演算部２０９、位置計算履歴演算部２１１及び位置履歴残差演算部２１２は、メモリに格納され、ＣＰＵによって実行されるプログラムである。状態ベクトル設定部２０７は、入力装置２０８からオペレータ先見情報を入力し、解析に用いる状態ベクトルや制約条件を設定する。解析処理部２０５は、観測データベース２０４上の観測データと状態ベクトル設定部２０７からの状態ベクトルや制約条件を入力し、これらのデータに最小２乗法や適応型フィルタの手法を適用して目標運動諸元に関する非線形解析を実施し、その解析結果を表示装置２０６に表示するとともに解析結果データベース２１０に格納する。
【００１９】
時系列残差演算部２０９は、解析処理部２０５が解析した状態ベクトルに基づいて目標の等速直進性を仮定して過去に遡った各観測時刻ごとの解析方位と観測データベース２０４に蓄積されている観測方位との方位残差を計算し、表示装置２０６に表示する。位置計算履歴演算部２１１は、解析処理部２０５が解析した状態ベクトルに基づいて目標の等速直進性を仮定して過去に遡った各観測時刻ごとの位置計算履歴を計算し、表示装置２０６に表示するとともに、位置履歴残差演算部２１２へ送出する。位置履歴残差演算部２１２は、位置計算履歴演算部２１１で計算された各観測時刻ごとの位置計算履歴と解析目標位置履歴との残差を計算し、表示装置２０６に表示する。
【００２０】
なおオペレータは、表示装置２０６の表示情報に基づいて入力された先見情報の妥当性を評価し、誤値であると判断した場合は、再度入力装置２０８からオペレータ先見情報を入力することができる。また表示装置２０６の表示情報に基づいて異常データの混在が認知された場合、オペレータが入力装置２０８から異常データを指定することにより、該当時刻の観測データを観測データベース２０４から検索排除し、解析処理部２０５により再度解析を実施することができる。
【００２１】
図３は、観測船３０１から見た目標移動体３０２の運動諸元を説明する図である。図で各記号の意味は次のとおりである。
【００２２】
座標原点：最新観測時刻tnにおける観測船の位置
^Rn     ：最新観測時刻tnにおける観測船からの解析目標距離
^Byn    ：最新観測時刻tnにおける観測船からの解析目標方位
^XTn    ：最新観測時刻tnにおける観測船からの解析目標X軸方向位置
^YTn    ：最新観測時刻tnにおける観測船からの解析目標Y軸方向位置
^Ctn    ：最新観測時刻tnにおける解析目標針路
^Mtn    ：最新観測時刻tnにおける解析目標速度
^VXn    ：最新観測時刻tnにおける解析目標X軸方向速度
^VYn    ：最新観測時刻tnにおける解析目標Y軸方向速度
解析対象となる目標移動体（以下目標という）の運動パラメータは、観測船から見た目標の距離Rと方位By、目標速力Mt及び目標針路Ctであるが、これらの運動パラメータは、^XT、^YT、^VX、^VYの４つのパラメータで表現できる。そこでこれらのパラメータを要素とする状態ベクトルとして（数１）を定義する。また目標が等速直進運動をして時刻tnに至ったと仮定したときの運動モデルは、（数３）に示す線形化された式で表現される。n-1の付いたパラメータは時刻tn-1の数値、(XOn,YOn)は時刻tnの観測船の位置座標である。以下^のついた記号は解析値、~のついた記号は誤差が含まれる可能性のある観測データを示す。
【００２３】
状態ベクトル設定部２０７は、入力装置２０８を介して入力されたオペレータ先見情報のタイプとその数値データを受け取り、入力されたタイプに応じた状態ベクトルをメモリに設定する。オペレータ先見情報のタイプとして、本実施例では、距離上下限（または上限のみ）、速力上下限（または上限のみ）、距離、速力、針路、距離と速力、距離と針路、および速力と針路の８種が可能である。（表１）に各オペレータ先見情報について適用される状態ベクトルを表現する数式を示す。解析手法として最小２乗法を用いる場合と、適応型フィルタを用いる場合についての状態ベクトルを示す。
【００２４】
【表１】
0003881967.tif 000002
図４は、最小２乗法を用いる解析処理部２０５の解析処理手順を示すフローチャートである。この解析手法は、マーカート法及びラグランジュ乗数法を利用する。以下に示す処理手順に適用する数式は、オペレータ先見情報が入力されない場合、距離上下限および速力上下限が入力された場合についての数式である。他のタイプのオペレータ先見情報が入力された場合に適用する数式については、後述する。解析処理部２０５は、観測データベース２０４に蓄積された全観測時刻の観測データをメモリに読み込む（ステップ４０６）。このときメモリ中の観測データのうち、オペレータが指定した異常データ及び音響センサの感度が劣化する方向の観測データを排除する。次に解析処理部２０５は、有効な観測データを用いて擬似線形解析処理を行い、初期解析値を算出する（ステップ４０７）。
【００２５】
擬似線形解析処理の手順は次のとおりである。解析処理部２０５は、前回の時刻について実行された非線形解析処理による解析値があるか判定する（ステップ４０２）。前回の解析値がある場合には、その解析値を初期解析値とし、ステップ４０３，４０４をスキップして非線形解析処理に戻る。前回の解析値がない場合には、まず擬似線形解析の準備処理を行う（ステップ４０３）。使用する状態ベクトルは、オペレータ先見情報が入力されない場合と、距離上下限および速力上下限が入力された場合には（数１）である。（数１）の数値から（数１２）に示す観測方位を計算するが、観測誤差が微小であることを仮定し、（数１２）の左辺は（数１３）のように展開される。（数１４）は、（数１２）および（数１３）から求めたものである。（数１５）は（数１４）を行列表記したものである。すなわち非線形問題は、（数１６）および（数１７）に示される線形モデルの最小化問題に置き換えられる。解析処理部２０５は、（数１４）（数１５）から係数行列H及び定数行列Gを計算する。次に解析処理部２０５は、（数１８）により初期解析値^X0nを算出する（ステップ４０４）。
【００２６】
次に非線形解析処理に戻り、解析処理部２０５は、（数５）に示す非線形解析に必要なパラメータを初期化する（ステップ４０８）。次にループ計算に入り、状態ベクトル^Xnを用いて（数４）に示す最小２乗法による評価関数値Ф（^Xn）を算出する（ステップ４０９）。次にФ（^Xn）の繰り返し計算のために状態ベクトルの修正量ΔXを算出する（ステップ４１０）。（数２）はベクトルΔXの要素を示す。まず（数４）を最小化するための（数６）（数７）（数８）を算出する。修正量ΔXは、通常（数９）により導出されるが、ここでは反復修正を安定して行うための縮小因子λを加えて（数１０）により計算する。
【００２７】
オペレータ先見情報として目標速力上限値および目標距離上限値が入力された場合には、（数１９）で示す制約条件が存在する。MAX(Mt)は速力の上限値、MAX(R)は距離の上限値である。この場合には（数４）の代わりにこの式を拡張した（数２０）に示す評価関数Lを適用する。ただし制約条件個数j=2とする。また修正量ΔXは、（数２１）（数２２）（数２３）を算出し、（数２４）により計算する。オペレータ先見情報として目標速力の上下限値および目標距離の上下限値が入力された場合には、制約条件として（数１９）を用い、（数２０）の評価関数Lの計算には制約条件個数j=4とし、（数２４）により修正量ΔXを計算すればよい。
【００２８】
次に解析処理部２０５は、状態ベクトル^Xnに修正量ΔXを加えて（数１１）に示す新たな状態ベクトル^Xnを算出する（ステップ４１１）。ここで^Xknは今回の^Xn、^Xk-1nは前回の^Xn又は初期解析値である。次に状態ベクトル^Xnの修正量ΔXが所定値以下となり状態ベクトル^Xnがほとんど修正されない状態となるか、あるいは縮小因子λの調整に伴う繰り返しループ（ループ１又はループ２）の回数が上限に達しているかを判定する（ステップ４１２）。条件を満足している場合には、ステップ４１３へ進む。
【００２９】
条件を満足していない場合には、前回の評価関数Фの値と今回の評価関数Фの値とを比較判定する（ステップ４１８）。評価関数の値が減少している場合には縮小因子λを倍にし、修正量ΔXを増加する方向へ調整し（ステップ４１９）、ステップ４０９に戻る。評価関数の値が変わらないか増加している場合には、縮小因子λを半減し、修正量ΔXを減少する方向へ調整し（ステップ４２０）、ステップ４０９に戻る。
【００３０】
ステップ４１２の条件を満足している場合、解析処理部２０５は、現在の状態ベクトル^Xnを用いて（数３１）により制約条件を満足していない条件式番号Igと、制約条件逸脱最大値gmaxを算出する（ステップ４１３）。次に制約条件逸脱最大値gmaxが前回の制約条件逸脱最大値c以下か、あるいはループ３の回数が上限に達しているかを判定する（ステップ４１４）。条件を満足している場合には、ステップ４１５へ進む。条件を満足していない場合には、（数３２）に示すようにペナルティパラメータriを更新し（ステップ４２２）、ステップ４０９に戻る。
【００３１】
ステップ４１４の条件を満足している場合、解析処理部２０５は、（数３３）に示すように制約条件逸脱最大値cを今回の制約条件逸脱最大値gmaxに置き換える（ステップ４１５）。次に制約条件逸脱最大値cが制約条件逸脱許容値εより小さくなっているか、あるいはループ３の回数が上限に達しているかを判定する（ステップ４１６）。条件を満足している場合には、現在の状態ベクトル^Xnから解析結果を算出する（ステップ４１７）。すなわち（数２５）により目標運動諸元（距離R、方位By、速力Mt、針路Ct）を算出する。また（数２６）（数２７）（数２８）（数２９）（数３０）により解析標準偏差を算出する。解析処理部２０５は、算出された解析結果を解析結果データベース２１０に格納し、また図形化して表示装置２０６上に表示する。ステップ４１６の条件を満足していない場合、（数３４）によりラグランジュ乗数uiを更新し（ステップ４２１）、ステップ４２２の処理を行った後にステップ４０９に戻る。
【００３２】
なお誤差をともなって観測される情報から内部の状態モデルを推定する最小２乗法及びマーカート法の詳細については、以下の文献に記載されている。
「最小二乗法による実験データ解析」中川徹  他  東京大学出版会
また内部の状態モデルに制約条件がある場合の最小２乗解の計算方法の詳細については以下の文献に記載されている。
「数値計算ハンドブック」大野他  オーム社。
【００３３】
オペレータ先見情報として距離上下限値及び速力上下限値が入力された場合には、上記処理手順に従って解析することにより、これらの入力値を反映した解析が可能となる。
【００３４】
オペレータ先見情報として目標距離Ropが入力された場合、入力された時刻の状態ベクトルを（数３８）に従って設定し、上記の擬似線形解析を実施する。（数３８）のr/nは、観測データ1,2,…nのうちオペレータ入力のあった時刻のデータを示す。次に（数３５）に従って置き換えられた状態ベクトルを用いて上記の非線形解析を実施する。すなわち（数１）の^XTr/n，^YTr/nを^Byr/nに置き換え、３次元ベクトルに縮退した状態ベクトルを用いる。非線形解析は、運動モデル（数３７）に基づいて、拡張評価関数L（数２０）のfm(^Xn)に（数３９）、g1(^Xn)、g2(^Xn)に（数４０）を適用し、（数２４）に（数４１）及び（数４２）を適用して状態ベクトル修正量ΔX（数３６）を算出することにより実行される。解析結果は（数４３）によって算出される。これにより、目標距離に関するオペレータ先見情報Ropを基準とした目標方位、速力、針路に関する解析が可能となる。
【００３５】
オペレータ先見情報として目標速力Mtopが入力された場合、状態ベクトルを（数４７）に従って設定し、上記の擬似線形解析を実施する。次に（数４４）に従って置き換えられた状態ベクトルを用いて上記の非線形解析を実施する。すなわち（数１）の^VXn，^VYnを^Ctnに置き換え、３次元ベクトルに縮退した状態ベクトルを用いる。非線形解析は、運動モデル（数４６）に基づいて、拡張評価関数L（数２０）のfm(^Xn)に（数４８）、g1(^Xn)、g2(^Xn)に（数４９）を適用し、（数２４）に（数５０）及び（数５１）を適用して状態ベクトル修正量ΔX（数４５）を算出することにより実行される。解析結果は（数５２）によって算出される。これにより、目標速力に関するオペレータ先見情報Mtopを基準とした目標距離、方位、針路に関する解析が可能となる。
【００３６】
オペレータ先見情報として目標針路Ctopが入力された場合、状態ベクトルを（数５３）に従って設定し、運動モデル（数５５）に基づいて（数５６）（数５７）から（数５８）の擬似線形方程式をたてて擬似線形解析を実施する。次に（数５３）の状態ベクトルを用いて上記の非線形解析を実施する。すなわち（数１）の^VXn，^VYnを^Mtnに置き換え、３次元ベクトルに縮退した状態ベクトルを用いる。非線形解析は、運動モデル（数５５）に基づいて、拡張評価関数L（数２０）のfm(^Xn)に（数５９）、g1(^Xn)、g2(^Xn)に（数６０）を適用し、（数２４）に（数６１）及び（数６２）を適用して状態ベクトル修正量ΔX（数５４）を算出することにより実行される。解析結果は（数６３）によって算出される。これにより、目標針路に関するオペレータ先見情報Ctopを基準とした目標距離、方位、速力に関する解析が可能となる。（数５３）の状態ベクトルは、By又はCtを含まないので、擬似線形解析においてもそのまま使用できる。
【００３７】
オペレータ先見情報として目標距離Ropと速力Mtopが入力された場合、入力された時刻の状態ベクトルを（数６７）に従って設定し、上記の擬似線形解析を実施する。次に（数６４）に従って置き換えられた状態ベクトルを用いて上記の非線形解析を実施する。すなわち（数１）の^XTr/n，^YTr/nを^Byr/nに置き換え、^VXr/n，^VYr/nを^Ctr/nに置き換え、２次元ベクトルに縮退した状態ベクトルを用いる。非線形解析は、運動モデル（数６６）に基づいて、拡張評価関数L（数２０）のfm(^Xn)に（数６８）、g1(^Xn)、g2(^Xn)に（数６９）を適用し、（数２４）に（数７０）及び（数７１）を適用して状態ベクトル修正量ΔX（数６５）を算出することにより実行される。解析結果は（数７２）によって算出される。これにより、目標距離と速力に関するオペレータ先見情報Rop、Mtopを基準とした目標方位、針路に関する解析が可能となる。
【００３８】
オペレータ先見情報として目標距離Ropと針路Ctopが入力された場合、入力された時刻の状態ベクトルを（数７６）に従って設定し、上記の擬似線形解析を実施する。次に（数７３）に従って置き換えられた状態ベクトルを用いて上記の非線形解析を実施する。すなわち（数１）の^XTr/n，^YTr/nを^Byr/nに置き換え、^VXr/n，^VYr/nを^Mtr/nに置き換え、２次元ベクトルに縮退した状態ベクトルを用いる。非線形解析は、運動モデル（数７５）に基づいて、拡張評価関数L（数２０）のfm(^Xn)に（数７７）、g1(^Xn)、g2(^Xn)に（数７８）を適用し、（数２４）に（数７９）及び（数８０）を適用して状態ベクトル修正量ΔX（数７４）を算出することにより実行される。解析結果は（数８１）によって算出される。これにより、目標距離と針路に関するオペレータ先見情報Rop、Ctopを基準とした目標方位、速力に関する解析が可能となる。
【００３９】
オペレータ先見情報として目標速力Mtopと針路Ctopが入力された場合、状態ベクトルを（数８２）に従って設定し、運動モデル（数８４）に基づいて（数８５）（数８６）から（数８７）の擬似線形方程式をたてて擬似線形解析を実施する。次に（数８２）の状態ベクトルを用いて上記の非線形解析を実施する。ここで（数８２）は、（数１）のベクトル要素のうち^XTn，^YTnのみを含む２次元ベクトルに縮退している。非線形解析は、運動モデル（数８４）に基づいて、拡張評価関数L（数２０）のfm(^Xn)に（数８８）、g1(^Xn)、g2(^Xn)に（数８９）を適用し、（数２４）に（数９０）及び（数９１）を適用して状態ベクトル修正量ΔX（数８３）を算出することにより実行される。解析結果は（数９２）によって算出される。これにより、目標速力、針路に関するオペレータ先見情報Mtop、Ctopを基準とした目標距離、方位に関する解析が可能となる。（数８２）の状態ベクトルは、By又はCtを含まないので、擬似線形解析においてもそのまま使用できる。
【００４０】
図５は、適応型フィルタを用いる解析処理部２０５の解析処理手順を示すフローチャートである。この解析手法は、カルマンフィルタを利用する。以下に示す解析処理に適用する状態ベクトルは、（数９８）に従うものであり、オペレータ先見情報が入力されない場合、距離上下限および速力上下限が入力された場合の状態ベクトルである。他のタイプのオペレータ先見情報が入力された場合に適用される状態ベクトルは、（表１）に挙げるものに置き換わる。以下の処理手順は、観測データベース２０４に蓄積されている各観測時刻ごとのデータを１つ１つ使用して状態ベクトル及び共分散行列の初期値を逐次更新することにより解析値を収束させる。
【００４１】
前回の時刻について実行された解析値がない初期状態の場合、解析処理部２０５は、（数９９）により状態ベクトルの初期値を設定し（ステップ５０２）、（数１００）により共分散行列の初期値を設定する（ステップ５０３）。ただし次元縮退している場合の（数９９）（数１００）は、各々の次元の形式に縮退する。
【００４２】
次に解析処理部２０５は、観測データベース２０４に蓄積された全観測時刻の観測データをメモリに読み込む（ステップ５０６）。このときメモリ中の観測データのうち、オペレータが指定した異常データ及び音響センサの感度が劣化する方向の観測データを排除する。次に解析処理部２０５は、（数１０１）に従い状態ベクトルの初期値を観測データの観測時刻まで時間遷移させる（ステップ５０７）。次に（数１０２）に従い共分散行列の初期値を観測データの観測時刻まで時間遷移させる（ステップ５０８）。さらにゲイン調節係数μを１．０に初期化する（ステップ５０９）。
【００４３】
次に解析処理部２０５は、（数１０３）により観測データに基づいて状態ベクトルを修正更新する（ステップ５１０）。次に（数１０４）により共分散行列を修正更新する（ステップ５１１）。次に（数１０５）により状態ベクトルを目標運動諸元に変換し、（数１０６）により共分散行列を目標運動諸元解析値の解析誤差標準偏差に変換する（ステップ５１２）。次に解析処理部２０５は、算出された目標運動諸元がオペレータ先見情報として入力された距離上下限および速力上下限の制約条件を満足しているか、あるいはμが−０．１以下かを判定する（ステップ５１３）。オペレータ先見情報として距離上限および速力上限のみが入力された場合には、各々上限のみの制約条件となる。満足していない場合はゲイン調節係数μを０．１だけ減算し（ステップ５１５）、ステップ５１０から５１３までの処理を繰り返す。満足している場合はすべての観測データを処理済みか判定する（ステップ５１４）。未処理データがある場合にはステップ５０６に戻り、未処理データがなくなったとき解析結果を解析結果データベース２１０に格納し、図形化して表示装置２０６上に表示し、処理を終了する。
【００４４】
なお誤差をともなって観測される情報から内部の状態モデルを推定するカルマンフィルタの詳細については以下の文献に記載されている。
「応用カルマンフィルタ」片山徹  朝倉書店。
【００４５】
上記のようにオペレータ先見情報が入力された場合に、（数１）又は（数９８）で示される目標の状態ベクトルの位置に関する要素又は速度に関する要素が別の運動パラメータに置き換えられ、３次元又は２次元の状態ベクトルに次元縮退する。この別の運動パラメータは、オペレータ入力された運動パラメータと連立すると、位置に関する要素又は速度に関する要素が決定されるような運動パラメータである。ただし目標速力Mtopと針路Ctopが入力された場合には、MtopとCtopによって速度に関する要素が決まるから、速度に関する要素が除去された状態ベクトルが設定される。
【００４６】
上記のようにオペレータ先見情報が入力された場合に、（数１）又は（数９８）で示される４次元の状態ベクトルが３次元又は２次元の状態ベクトルに次元縮退するか、または解析処理に制約条件が課せられるために、解析処理中の探索範囲が狭まり、結果として解析収束までの処理時間を短縮することができる。
【００４７】
時系列残差演算部２０９は、解析結果データベース２１０上の解析結果を参照し、（数９３）により解析方位を算出し、（数９４）により観測方位と解析方位との方位残差を求め、その結果を表示装置２０６上に表示する。図６は、縦軸６０２を時間とし横軸６０１を方位残差として、観測方位ごとの方位残差の表示例を示す。方位残差６０３はオペレータ先見情報入力なしの場合の解析結果に基づき計算された方位残差、方位残差６０４はオペレータ先見情報入力ありの場合の解析結果に基づき計算された方位残差を示している。この例ではMtに関するオペレータ先見情報が入力された場合の例を示している。方位残差データが残差０の中心線に近接しているほど解析結果が正しいと判断することができる。この例によれば、オペレータは、先見情報入力前後の比較、解析収束判断、異常データ６０５の混在判断を容易に行える。
【００４８】
位置計算履歴演算部２１１は、解析結果データベース２１０上の解析結果を参照し、（数９７）により目標の等速直進性を仮定して過去に遡った目標の位置履歴を計算し、（数９６）により各解析時刻の解析目標位置を計算する。ベクトルDn(m)は、現在時刻nから遡った時刻mの位置ベクトルを示す。位置履歴残差演算部２１２は、これらの計算結果から解析目標位置履歴の的針線残差を計算してその結果を表示装置２０６上に表示する。図７は、縦軸７０２を目標位置のトレースラインとし、横軸７０１を的針線残差として、観測データごとの残差の表示例を示す。７０３はオペレータ先見情報入力なしの場合の解析結果に基づき計算された解析目標位置残差、７０４はオペレータ先見情報入力ありの場合の解析結果に基づき計算された解析目標位置残差である。解析目標位置残差が残差０の中心線（トレースライン）に近接しているほど解析結果が正しいと判断することができる。この例によれば、オペレータは、先見情報入力前後の比較、解析収束判断、異常データ７０５の混在判断を容易に行える。
【００４９】
位置計算履歴演算部２１１は、解析結果データベース２１０上の解析結果を参照し、（数９３）により解析方位を計算し、観測方位ξ(m)とともに、表示装置２０６上に表示する。図８は、縦軸を時間とし横軸を方位として、時系列の方位履歴の表示例を示す。８０２は観測方位、８０３は解析方位を示す。また８０７はオペレータ先見情報入力なしの場合の解析結果に基づき計算された解析方位のトレースライン、８０４はオペレータ先見情報入力ありの場合の解析結果に基づき計算された解析方位８０３のトレースラインである。観測方位と解析方位が近接しているほど解析結果が正しいと判断することができる。この例によれば、オペレータは、先見情報入力前後の比較、解析収束判断を容易に行える。
【００５０】
位置計算履歴演算部２１１は、解析結果データベース２１０上の解析結果を参照し、（数９５）により目標の等速直進性を仮定して過去に遡った目標の位置計算距離^Rn(m)を計算し、解析距離^Rmとともに表示装置２０６上に表示する。図９は、縦軸を距離とし横軸を時間として、時系列の解析距離及び位置計算距離の表示例を示す。９０２，９０７は解析距離、９０３は位置計算距離、９０４は位置計算距離９０３のトレースラインである。９０７はオペレータ先見情報入力なしの場合の解析距離、９０２はオペレータ先見情報入力ありの場合の解析距離である。解析距離９０２，９０７と位置計算距離９０３が近接しているほど解析結果が正しいと判断することができる。この例によれば、オペレータは、先見情報入力前後の比較、解析収束判断を容易に行える。
【００５１】
位置計算履歴演算部２１１は、解析結果データベース２１０上の解析結果を参照し、（数９７）により最新解析時刻における解析結果に基づいて目標の等速直進性を仮定して過去に遡った目標位置を計算し、（数９６）により各解析時刻における解析目標位置を計算し、表示装置２０６上に表示する。図１０は、横軸をX方向、縦軸をY方向として、時系列の位置計算による位置および解析目標位置の表示例を示す。１００９は観測船位置、１００８は最新観測時刻における観測船の位置を示す。１００２，１００７は解析目標位置、１００３は目標位置計算による位置、１００４は位置１００３のトレースラインである。１００７はオペレータ先見情報入力なしの場合の解析目標位置、１００２はオペレータ先見情報入力ありの場合の解析目標位置である。１００５，１００６は最新の位置データである。解析目標位置１００２と位置１００３とが近接しているほど解析結果が正しいと判断することができる。この例によれば、オペレータは、先見情報入力前後の比較、解析収束判断を容易に行える。
【００５２】
図１１は、シミュレーションによって従来技術による解析結果と本発明による解析結果を比較する例を示す。両解析結果は、目標速力に関するオペレータ先見情報を反映している。１１００は観測船の位置、１１０２は従来技術による解析目標位置、１１０３は本発明による解析目標、１１０１は目標運動ベクトル（真値）を示す。解析目標１１０３は、解析目標１１０２に比べて早い段階から目標運動ベクトル１１０１に収束していることがわかる。また１１０４は本発明による目標距離解析値の時間変化、１１０５は従来技術による目標距離解析値の時間変化を示す。１１０６は本発明による目標速力解析値の時間変化、１１０７は本発明による目標速力解析値の時間変化を示す。１１０８は本発明による目標針路解析値の時間変化、１１０９は従来技術による目標針路解析値の時間変化を示す。この例では、目標速力解析値１１０６は目標速力解析値１１０７と差がないが、目標距離解析値１１０４、目標針路解析値１１０８は、目標距離解析値１１０５、目標針路解析値１１０９に比べて早い段階から収束することがわかる。これは、従来技術が状態ベクトルを（数１）または（数９８）と定義してオペレータ先見情報を反映しているため、目標距離解析値１１０５と目標針路解析値１１０９との間に、評価関数Ф（^Xn）（数４）が最小となる相関関係が生成されることが原因である。本発明によれば、状態ベクトルを（数４４）のように設定し、オペレータ先見情報を解析対象パラメータから縮退することによってこの問題を解決する。
【００５３】
図１２は、目標速力に関するオペレータ先見情報として誤値を入力した場合の解析結果をシミュレーション例で示したものである。解析目標位置（先見情報：誤値）１２０２と位置計算位置（最新解析結果に基づくトレースライン）１２０１との乖離状況、目標距離解析値（先見情報：誤値）１２０４と目標距離（最新解析結果に基づくトレースライン）１２０３との乖離状況、および方位残差（先見情報：誤値）１２０５の乖離状況からオペレータ先見情報が誤値であることが判断可能である。
【００５４】
図１３は、目標速力に関するオペレータ先見情報として真値を入力した場合の解析結果をシミュレーション例で示したものである。解析目標位置（先見情報：真値）１３０２と位置計算位置（最新解析結果に基づくトレースライン）１３０１との合致状況、目標距離解析値（先見情報：真値）１３０４と目標距離（最新解析結果に基づくトレースライン）１３０３との合致状況、および方位残差（先見情報：真値）１３０５の合致状況からオペレータ先見情報が真値であることが判断可能である。
【００５５】
以下に本実施形態で参照した（数１）〜（数１０６）を挙げる。
【００５６】
【数１】
0003881967.tif 000003
【数２】
0003881967.tif 000004
【数３】
0003881967.tif 000005
【数４】
0003881967.tif 000006
【数５】
0003881967.tif 000007
【数６】
0003881967.tif 000008
【数７】
0003881967.tif 000009
【数８】
0003881967.tif 000010
【数９】
0003881967.tif 000011
【数１０】
0003881967.tif 000012
【数１１】
0003881967.tif 000013
【数１２】
0003881967.tif 000014
【数１３】
0003881967.tif 000015
【数１４】
0003881967.tif 000016
【数１５】
0003881967.tif 000017
【数１６】
0003881967.tif 000018
【数１７】
0003881967.tif 000019
【数１８】
0003881967.tif 000020
【数１９】
0003881967.tif 000021
【数２０】
0003881967.tif 000022
【数２１】
0003881967.tif 000023
【数２２】
0003881967.tif 000024
【数２３】
0003881967.tif 000025
【数２４】
0003881967.tif 000026
【数２５】
0003881967.tif 000027
【数２６】
0003881967.tif 000028
【数２７】
0003881967.tif 000029
【数２８】
0003881967.tif 000030
【数２９】
0003881967.tif 000031
【数３０】
0003881967.tif 000032
【数３１】
0003881967.tif 000033
【数３２】
0003881967.tif 000034
【数３３】
0003881967.tif 000035
【数３４】
0003881967.tif 000036
【数３５】
0003881967.tif 000037
【数３６】
0003881967.tif 000038
【数３７】
0003881967.tif 000039
【数３８】
0003881967.tif 000040
【数３９】
0003881967.tif 000041
【数４０】
0003881967.tif 000042
【数４１】
0003881967.tif 000043
【数４２】
0003881967.tif 000044
【数４３】
0003881967.tif 000045
【数４４】
0003881967.tif 000046
【数４５】
0003881967.tif 000047
【数４６】
0003881967.tif 000048
【数４７】
0003881967.tif 000049
【数４８】
0003881967.tif 000050
【数４９】
0003881967.tif 000051
【数５０】
0003881967.tif 000052
【数５１】
0003881967.tif 000053
【数５２】
0003881967.tif 000054
【数５３】
0003881967.tif 000055
【数５４】
0003881967.tif 000056
【数５５】
0003881967.tif 000057
【数５６】
0003881967.tif 000058
【数５７】
0003881967.tif 000059
【数５８】
0003881967.tif 000060
【数５９】
0003881967.tif 000061
【数６０】
0003881967.tif 000062
【数６１】
0003881967.tif 000063
【数６２】
0003881967.tif 000064
【数６３】
0003881967.tif 000065
【数６４】
0003881967.tif 000066
【数６５】
0003881967.tif 000067
【数６６】
0003881967.tif 000068
【数６７】
0003881967.tif 000069
【数６８】
0003881967.tif 000070
【数６９】
0003881967.tif 000071
【数７０】
0003881967.tif 000072
【数７１】
0003881967.tif 000073
【数７２】
0003881967.tif 000074
【数７３】
0003881967.tif 000075
【数７４】
0003881967.tif 000076
【数７５】
0003881967.tif 000077
【数７６】
0003881967.tif 000078
【数７７】
0003881967.tif 000079
【数７８】
0003881967.tif 000080
【数７９】
0003881967.tif 000081
【数８０】
0003881967.tif 000082
【数８１】
0003881967.tif 000083
【数８２】
0003881967.tif 000084
【数８３】
0003881967.tif 000085
【数８４】
0003881967.tif 000086
【数８５】
0003881967.tif 000087
【数８６】
0003881967.tif 000088
【数８７】
0003881967.tif 000089
【数８８】
0003881967.tif 000090
【数８９】
0003881967.tif 000091
【数９０】
0003881967.tif 000092
【数９１】
0003881967.tif 000093
【数９２】
0003881967.tif 000094
【数９３】
0003881967.tif 000095
【数９４】
0003881967.tif 000096
【数９５】
0003881967.tif 000097
【数９６】
0003881967.tif 000098
【数９７】
0003881967.tif 000099
【数９８】
0003881967.tif 000100
【数９９】
0003881967.tif 000101
【数１００】
0003881967.tif 000102
【数１０１】
0003881967.tif 000103
【数１０２】
0003881967.tif 000104
【数１０３】
0003881967.tif 000105
【数１０４】
0003881967.tif 000106
【数１０５】
0003881967.tif 000107
【数１０６】
0003881967.tif 000108
【発明の効果】
以上説明したように、本発明によれば、入力されたオペレータ先見情報を利用して目標運動解析の収束時間を短縮することができる。
【図面の簡単な説明】
【図１】目標運動解析の原理を示す図である。
【図２】実施形態の目標運動解析装置の構成を示す図である。
【図３】目標移動体の運動諸元を説明する図である。
【図４】実施形態の最小２乗法を用いる解析処理手順を示すフローチャートである。
【図５】実施形態の適応型フィルタを用いる解析処理手順を示すフローチャートである。
【図６】時系列方位残差の表示例を示す図である。
【図７】解析目標位置履歴残差の表示例を示す図である。
【図８】観測方位と解析方位の時系列データの表示例を示す図である。
【図９】解析距離の時系列履歴データと最新の解析結果に基づく位置計算履歴データの表示例を示す図である。
【図１０】解析目標位置履歴データと最新の解析結果に基づく位置計算履歴データの対応関係を観測船の位置履歴データとともに示す表示例の図である。
【図１１】目標速力に関するオペレータ先見情報を入力した場合の従来技術による解析結果と本発明技術による解析結果の比較をシミュレーション例で示す図である。
【図１２】目標速力に関するオペレータ先見情報として誤値を入力した場合の解析結果をシミュレーション例で示す図である。
【図１３】目標速力に関するオペレータ先見情報として真値を入力した場合の解析結果をシミュレーション例で示す図である。
【符号の説明】
２０４：観測データベース、２０５：解析処理部、２０６：表示装置、２０７：状態ベクトル設定部、２０８：入力装置、２０９：時系列残差演算部、２１０：解析結果データベース、２１１：位置計算履歴演算部、２１２：位置履歴残差演算部、２１３：目標運動解析装置。
【図１】
0003881967.tif 000109
【図２】
0003881967.tif 000110
【図３】
0003881967.tif 000111
【図４】
0003881967.tif 000112
【図５】
0003881967.tif 000113
【図６】
0003881967.tif 000114
【図７】
0003881967.tif 000115
【図８】
0003881967.tif 000116
【図９】
0003881967.tif 000117
【図１０】
0003881967.tif 000118
【図１１】
0003881967.tif 000119
【図１２】
0003881967.tif 000120
【図１３】
0003881967.tif 000121