(57)【要約】【課題】目標移動体の目標運動パラメータを高精度に求めることができる目標運動解析方法及び目標運動解析装置を提供する。【解決手段】目標運動解析装置１００は、水上を運動する複数の目標移動体から放射される音波を、パッシブセンサ１によって取得し、取得された音波の観測方位と周波数を観測データとして時系列に記憶部１０２に記憶する。解析処理部２０３は、記憶部１０２に格納された時系列の観測データを用いて、入力装置１０３を介してオペレータが指定した複数の目標の中の少なくとも二つの特定の目標の距離、方位、速度および針路を含む目標運動パラメータを、特定の目標相互の位置関係を規定する制約条件を用いて、それぞれの目標運動パラメータを算出する。【選択図】図１
2014025786.tif 000002
【特許請求の範囲】
【請求項１】
水上を運動する複数の目標移動体から放射される音波もしくは電磁波の信号を、パッシブセンサによって取得し、該パッシブセンサにより取得された音波もしくは電磁波の信号を観測データとして時系列に格納する記憶部と、
該記憶部に格納された時系列の前記観測データを用いて前記複数の目標移動体の中のオペレータによって指定された少なくとも二つの特定の目標移動体の距離、方位、速度および針路を含む目標運動パラメータを、前記特定の目標移動体相互の位置関係を規定する情報を用いて、それぞれの目標運動パラメータを算出する解析処理部と、
を備えることを特徴とする目標運動解析装置。
【請求項２】
前記特定の目標移動体相互の位置関係を規定する情報とは、前記算出される目標運動パラメータのうちの前記距離に対する最大距離制限値及び最小距離制限値の制約条件であることを特徴とする請求項１に記載の目標運動解析装置。
【請求項３】
前記制約条件は、自艦が潜水艦である場合に、該潜水艦の現在位置する海域において、対潜艦艇がアクティブソナーを用いて潜水艦の探索をする際に、音波の屈折により潜水艦の探索が可能な対潜艦艇自身からの有効水平方向半径距離のデータに基づいて算出される前記最大距離制限値及び最小距離制限値であり、
該有効水平方向半径距離のデータに基づいて前記制約条件を算出設定する制約条件算出部を備えることを特徴とする請求項２に記載の目標運動解析装置。
【請求項４】
センサによって取得された目標移動体の観測データを用いて、前記目標移動体の距離、方位、速度および針路を含む目標運動パラメータを解析する目標運動解析方法であって、
前記目標移動体が複数であり、
オペレータによって前記複数の目標移動体のうちの少なくとも二つの目標移動体を特定目標移動体としての選択をうける目標移動体選択工程と、
該２つの特定目標移動体同士の間の相対位置関係に係る目標相対位置関係データと、前記センサによって検知された前記特定目標移動体の時系列の前記観測データと、自身の位置及び針路のデータとを用いて、少なくとも前記特定目標移動体の位置、速度、針路を含む目標状態ベクトルをそれぞれ推定演算する目標状態ベクトル推定演算工程と、
前記特定目標移動体の新たな観測データと、観測データを取得したタイミングに対応する推定演算された前記目標状態ベクトルとの差分に対する評価値を算出する評価値算出工程と、
目標状態ベクトル推定演算工程で演算された前記目標状態ベクトルの各要素で偏微分した微係数行列により算出された前記目標状態ベクトルの更新量と、繰り返し計算における前回の前記評価値と今回の前記評価値の差分と、前記目標相対位置関係データに対応して設定された最大距離制限値及び最小距離制限値の制約条件と、に基づいて、前記目標状態ベクトルが収束しておらず前記目標状態ベクトル推定演算工程、前記評価値算出工程を繰り返すか、前記目標状態ベクトルが収束したかを判定する目標状態ベクトル収束判定工程と、
該目標状態ベクトル収束判定工程において、前記目標状態ベクトルが収束したと判定されたとき、前記２つの特定目標移動体の目標運動パラメータを出力する目標運動パラメータ出力工程と、を含むことを特徴とする目標運動解析方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、船舶、航空機等、目標移動体の目標運動パラメータを解析する技術に関する。
【背景技術】
【０００２】
センサにより移動する目標（目標移動体）の方位、もしくは目標の方位と目標の放射する音響、電磁波等の周波数とが観測されている状況において、これらの誤差を含んだ情報を用いて、観測によっては直接的に得られない目標までの距離、目標の移動速度、目標の移動方向（針路）等の目標移動体の目標運動パラメータを計算する方法として、最小二乗法を用いた方法が特許文献１に示されている。
【先行技術文献】
【特許文献】
【０００３】
【特許文献１】特開２００５−００９９５４号公報
【発明の概要】
【発明が解決しようとする課題】
【０００４】
しかしながら、特許文献１に記載された技術は、単一の目標移動体を対象としたものであり、二つ以上の目標移動体が、例えば、船団や編隊を組んだグループとして行動している場合であっても、個別の目標移動体毎に解析を行うものである。
そして、センサによって検知された目標移動体の観測データを用いて目標移動体の距離、方位、移動速度及び針路を含む目標運動パラメータを解析するに当たって、目標移動体の距離、移動速度、針路のうちの少なくとも１つの目標運動パラメータについて先験情報に基づいてオペレータ入力を受け付け、オペレータ入力された目標運動パラメータを除く次元縮退した目標状態ベクトルを設定し、設定された次元縮退された目標状態ベクトルを用いて解析処理を行う技術である。
【０００５】
従って、目標移動体の距離、移動速度、針路のうちの少なくとも１つの目標運動パラメータを先験情報に基づいてオペレータ入力するとした場合、又はオペレータ入力された目標運動パラメータに大きな誤差がある場合、目標移動体の距離、移動速度及び針路の目標運動パラメータの内の残りの二つの目標運動パラメータの推定精度が低下したり、オペレータ入力された目標運動パラメータの補正を行う工程を伴うために、精度の高い目標移動体の距離、方位、移動速度及び針路の目標運動パラメータを得るのに時間を要したりする問題があった。
【０００６】
本発明は、前記した従来の課題を解決するものであり、目標移動体の目標運動パラメータを高精度に求めることができる目標運動解析方法及び目標運動解析装置を提供することを目的とする。
【課題を解決するための手段】
【０００７】
前記課題を解決するために、本発明に係る目標運動解析装置は、複数の目標移動体に対する観測データを入力情報とし、複数の目標間の距離に関する制約条件を考慮したものを用いて、目標運動パラメータの最適解を算出する。
【０００８】
本発明によれば、グループを形成する複数の目標移動体は、相互に一定の位置関係を保って行動するのが普通である。その複数の目標移動体間の相対位置の関係、特に相対距離は、先験情報として推定可能なものである。前記制約条件下で複数の目標移動体からの観測データに関連を持たせた最適化問題を解くことにより、当該複数目標の目標運動パラメータを単独で推定するよりも高精度に求めることができる。
【発明の効果】
【０００９】
本発明によれば、目標移動体の目標運動パラメータを高精度に求めることができる目標運動解析方法及び目標運動解析装置を提供することができる。
【図面の簡単な説明】
【００１０】
【図１】本発明に係る目標運動解析装置の機能構造ブロック図である。
【図２】目標運動解析装置における目標移動体の目標運動パラメータの説明図である。
【図３】一つのグループとして編隊を組んで行動する複数の目標移動体、例えば、艦船の編隊における２隻の艦船の相対位置関係の説明図である。
【図４】目標運動解析装置における、例えば、２隻の目標艦船に対する目標運動パラメータ解析計算の全体の流れを示すフローチャートである。
【図５】目標運動解析装置の制約条件算出部が行う処理の説明図である。
【図６】解析処理部における目標運動パラメータ解析計算の流れを示すフローチャートである。
【図７】解析処理部の出力する目標運動パラメータの説明図である。
【発明を実施するための形態】
【００１１】
以下、本発明の実施形態に係る目標運動解析方法及び目標運動解析装置１００について図面を参照して説明する。
本実施形態では、例えば、編隊を組んで航行する水上艦艇を潜水艦が目標移動体として捉える場合を例に、潜水艦に搭載される目標運動解析装置１００について説明する。
【００１２】
図１は、本発明に係る目標運動解析装置の機能構造ブロック図である。
図１に示すように目標運動解析装置１００には、パッシブセンサ１が検出した信号を処理するセンサ信号処理装置２からの処理結果のデータ（観測データ）、航海センサ５からの信号に基づいて自艦の位置を演算する位置演算部６からの自艦位置情報、水温センサ７が検出した海水温を示す信号等が入力される。
【００１３】
ここで、パッシブセンサ１は、例えば、潜水艦の艦首部に設けられた多数のハイドロフォン（受波器）の球形状に配置されたアレイから構成された艦首部パッシブソナー、潜水艦の艦体の左右両側面に配置された受波器のアレイから構成された艦体側面パッシブソナー、潜水艦の艦首部から艦体の左右両側面に配置された受波器のアレイから構成されたコンフォーマルパッシブソナー、潜水艦の艦尾から繰り出されて曳航される曳航式パッシブソナー（ＴＡＳＳ(Towed Array Sonar System)）のうちのいずれか一つ又はその組み合わせである。
【００１４】
ちなみに、ＴＡＳＳは、曳航ケーブルに沿って複数個所にハイドロフォンを配置したものであり、ハイドロフォンが曳航ケーブルの周方向に複数個配置されて、音源の方位を検出可能な多局、例えば、４局のものである。
【００１５】
パッシブセンサ１からの信号は、センサ信号処理装置２において処理され、個々の水上艦艇からのものと考えられる音源からの音波を、水平方位角で走査し、同一方位角で周波数スペクトル分布（音紋）が同一のものを１つの目標からの音波信号（音紋）として同定する。各目標の音紋は、例えば、水上艦のエンジン音、スクリュー音が組み合わさったものであり、その目標の水上艦艇の主機関がデーゼルエンジンかガスタービンエンジンかによっても異なり、又、スクリューの回転速度やスクリューの翼の大きさや形状、推進軸が何本かによっても異なり、その特徴から編隊航行をする艦艇の種類の差として判別できる。このようなセンサ信号処理装置２における音紋の分析により、鯨等の発する音声と水上艦艇からの音波を区別して、センサ信号処理装置２の表示装置（図示せず）に、音紋の周波数スペクトルを表示したり、音源の方位角を表示したり、抽出した音紋をヘッドフォンに出力してソナー要員に目標音源を確認させたりすることができるようになっている。
【００１６】
センサ信号処理装置２は、音源の自艦からの方位角（観測方位）を演算するとともに、算出した方位角の誤差（標準偏差）をも算出し、目標運動解析装置１００に出力する。また、センサ信号処理装置２は、音源の音紋から、例えば、強度の一番強い周波数を選択してその音源の代表周波数（以下「観測周波数」と称する）として特定し、特定された観測周波数とその観測周波数の算出誤差（標準偏差）も目標運動解析装置１００に出力する。
【００１７】
この際、センサ信号処理装置２は、前記したパッシブセンサ１を構成している艦首部パッシブソナー、コンフォーマルパッシブソナー、艦体側面パッシブソナー、曳航式パッシブソナー等のそれぞれからの信号を総合的に組み合わせ評価して、音源の自艦からの方位角（観測方位）を演算するとともに、算出した方位角の誤差（標準偏差）をも算出し、目標運動解析装置１００に出力する。また、センサ信号処理装置２は、音源の音紋から、強度の一番強い周波数を選択してその音源の観測周波数として特定し、特定された観測周波数とその観測周波数の算出誤差（標準偏差）も目標運動解析装置１００に出力するという精度の高い観測データを出すこともできる。
【００１８】
その他に、センサ信号処理装置２は、艦首部パッシブソナー、コンフォーマルパッシブソナー、艦体側面パッシブソナー、曳航式パッシブソナー等のそれぞれが取得した独自の信号により、音源の自艦からの方位角（観測方位）を演算する。その際に算出した方位角の誤差（標準偏差）をも算出し、目標運動解析装置１００に出力する。また、音源の音紋から、強度の一番強い周波数を選択してその音源の観測周波数として特定し、特定された観測周波数とその観測周波数の算出誤差（標準偏差）も目標運動解析装置１００に出力する。
【００１９】
艦首部パッシブソナー、コンフォーマルパッシブソナー、艦体側面パッシブソナー、曳航式パッシブソナーのそれぞれごとに、前記した方位角（観測方位）、方位角の誤差（標準偏差）、観測周波数、その観測周波数の算出誤差（標準偏差）を目標運動解析装置１００に出力する場合は、センサ信号処理装置２は、艦首部パッシブソナー、コンフォーマルパッシブソナー、艦体側面パッシブソナー、曳航式パッシブソナーの配置されている相対配置の距離から、特定の音紋を有する音源に対し、三角測量の原理で音源の自艦からの方位角と距離とを、誤差は大きいが概算的にセンサ信号処理装置２において算出でき、その算出結果をも目標運動解析装置１００に出力するようにしても良い。
【００２０】
ここで、センサ信号処理装置２は、センサ信号処理装置２を操作するオペレータ用の図示しない入力装置、表示装置を有している。目標運動解析装置１００を操作するオペレータ席とセンサ信号処理装置２を操作するオペレータ席とは、並列して配置され、センサ信号処理装置２の表示装置と目標運動解析装置１００の後記する表示装置１０４も並列して配置されている。このような配置とすることにより、この２つのセンサ信号処理装置２と目標運動解析装置１００のそれぞれのオペレータが並列して座り、オペレータは、互いにコミュニケーションをして共同作業が行える配置であることが好ましい。
【００２１】
航海センサ５は自艦の運動を計測し、位置演算部６に入力する。位置演算部６は、自艦の位置、針路、速度等を算出して自艦運動情報として目標運動解析装置１００に出力する。針路とは、自艦の進行方向を、真北（真の北の方位）方向を基準に右回りに角度表示の方向で示したものである。ちなみに、磁気コンパスによる北方向は、地球磁場の各地域での偏差があるので、その偏差データを基に、磁気コンパスの示す真北方位を補正して針路を算出する。
航海センサ５は、例えば、慣性航法装置(ＩＮＳ：Inertial Navigation System)用のジャイロスコープ（ジャイロセンサとも呼ばれる）や、磁気コンパスや、水中での対水速力を検出するログ（速力計）等を含んでおり、それらの信号に基づいて位置演算部６において前記した自艦運動情報が算出される。
【００２２】
なお、位置演算部６が、慣性航法の機能やログによる補正機能を有していても、ログは水に対する相対的速度を示すだけなので、位置演算部６で算出した自艦の位置には、慣性航法の演算誤差、海流による偏流誤差等が蓄積される。したがって、位置演算部６は、図示しないが、自艦（潜水艦）がマストを露頂させた際に、ＧＰＳ（Global Positioning System）により、自艦の位置、針路等を補正する。
【００２３】
水温センサ７は、海水温を測定するセンサであり、水中音速を推定するために用いられるものである。水温センサ７としては、例えば、温度計、海水の電気抵抗測定プローブ等がある。
ちなみに、水中音速は、海水温度により主に支配されるが、海水の塩分濃度の影響も受ける。
【００２４】
目標運動解析装置１００は、演算部であるＣＰＵ１０１、ハードディスク装置や半導体記憶装置等で構成された記憶部１０２、入力装置１０３、表示装置１０４、ＲＡＭ等を有する計算機で構成されている。
記憶部１０２には、目標運動解析装置１００用の解析プログラム等を格納したプログラム１０２ａ、観測データを格納する領域である観測データベース１０２ｂ、自艦運動情報を格納する領域である自艦運動データベース１０２ｃ、解析の途中結果をも含めた解析結果を格納する領域である解析結果データベース１０２ｄ、海域音響特性データベース１０２ｅを含んでいる。
【００２５】
ここで、海域音響特性データベース１０２ｅとは、音響測定艦が各海域の季節ごとに海水温度データ分布、塩分濃度データ等を測定し、対潜戦闘艦艇（以下、単に「対潜艦艇」と称する）がアクティブソナーを用いて潜水艦の探索をする際に、音波の屈折により潜水艦の探索が可能な対潜艦艇自身からの有効水平方向半径距離のデータ等である。
【００２６】
ＣＰＵ１０１は、記憶部１０２からプログラム１０２ａを読み込んで実行する機能部としての、データ取得部２０１、制約条件算出部２０２、解析処理部２０３、入力制御部２０５、表示制御部２０７を含んでいる。
データ取得部２０１は、センサ信号処理装置２から入力される目標移動体である音源の自艦からの観測方位、観測方位の標準偏差、音源の観測周波数、観測周波数の標準偏差を観測データとして時系列的に取得し、観測データを取得したタイミングにおける時刻情報を付加して記憶部１０２の観測データベース１０２ｂに時系列的に記憶させる。
【００２７】
また、データ取得部２０１は、位置演算部６から、観測データを取得したタイミングにおける自艦位置、自艦の速度、自艦の針路等のデータを含む自艦運動情報を取得して、時刻情報を付加して記憶部１０２の自艦運動データベース１０２ｃに時系列的に記憶させる。
更に、データ取得部２０１は、適宜、水温センサ７からの信号を取得して、観測データベース１０２ｂに時刻情報を付加して時系列的に記憶させる。
【００２８】
入力制御部２０５は、入力装置１０３からの指示に従って、制約条件算出部２０２における処理や解析処理部２０３における処理を制御したり、表示制御部２０７を介して表示装置１０４に表示させる内容を制御したりする。ちなみに、入力装置１０３は、キーボードやマウスに限定されるものではなく、表示装置１０４の表示画面に設けられたタッチパネル等の画面入力装置をも含む。
【００２９】
（制約条件算出部２０２）
次に、制約条件算出部２０２におけるオペレータの入力による制約条件算出機能について説明する。
入力制御部２０５は、検出された複数の目標移動体（船団、艦艇の編隊等）の中からオペレータの先験情報に基づく入力装置１０３の操作による解析の対象とする自艦にとって現時点で最大脅威の目標である基準目標の選定、基準目標（目標Ａ（目標移動体））に対して所定の距離の間隔を取って編隊行動している僚艦目標（目標Ｂ（目標移動体）の選定を受け付ける。この際の先験情報とは、パッシブセンサ１の受信信号をセンサ信号処理装置２で処理した結果に基づいて目標移動体が、対潜艦艇ならば、アクティブソナーを定期的に用いて、潜水艦の存在の有無を確認していることから容易に対潜艦艇であると判断でき、アクティブソナーを用いていない場合は、その音紋から艦隊防空用の艦艇、護衛艦艇により護衛されている船団又は編隊の艦船や船舶等（揚陸艦、航空母艦、輸送船等）であるとの判定ができる。
【００３０】
目標運動解析装置１００のオペレータは、目標Ａと目標Ｂとが、先験情報からともに対潜艦艇と判定した場合、入力装置１０３を操作して制約条件算出部２０２に目標Ｂの目標運動パラメータの内の自艦からの最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎ等の制約条件を算出させる。
【００３１】
制約条件算出部２０２は、目標Ａ，Ｂが共に対潜艦艇であるとの先験情報を基に、前記した現在自艦が所在している海域の海域音響特性データベース１０２ｅに基づいて、潜水艦の探索が可能な対潜艦艇自身からの有効水平方向半径距離のデータを検索し、目標Ａ，Ｂは、目標Ｂの目標運動パラメータに対する最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎ等の制約条件（図５参照）を算出する。これは、対潜艦艇は、その対潜艦艇が護衛艦艇として護衛している船団や編隊の外周に位置し、潜水艦の探索が可能な対潜艦艇自身からの有効水平方向半径距離が互いに一部オーバーラップさせて、護衛艦艇の警戒網をかいくぐって潜水艦の侵入を許さないように、互いの距離を所定の距離に保つように警戒配置で航行することが通常であるからである。
【００３２】
そして、目標運動解析装置１００のオペレータは、例えば、自艦により近い側の対潜艦艇を目標Ａとし、その近傍には、自艦からの距離が目標Ａよりも遠いが、対潜艦艇の目標Ｂが存在するとして、目標Ｂを選択する。自艦により近い側の対潜艦艇である目標Ａは、音源の強度が高いことや、概算による三角測量により算出した目標までの距離が近いことから容易に選択できる。
この最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎの定義については、後記する。
【００３３】
なお、目標Ａが対潜艦艇で無く、護衛されている艦艇、船舶であり、目標Ｂが対潜艦艇であるとオペレータが判定した場合や、目標Ａ、Ｂ共に対潜艦艇ではないと判定した場合は、潜水艦の探索が可能な対潜艦艇自身からの有効水平方向半径距離のデータを用いず次のように最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎを算出設定する。センサ信号処理装置２の複数の目標のそれぞれの音紋と三角測量の原理による概算の距離、方位による船団又は編隊の陣形から、目標Ａを中心とした所定の半径から算出する目標Ｂの距離に対する最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎ等の制約条件を設定する。
例えば、複数の目標移動体の編隊が単縦陣型の場合は、比較的接近して先行艦の艦尾が目視できる距離、先行艦からの発光信号を確認可能な距離である。
大規模な艦艇部隊や船団（艦隊）が対空警戒配置で航行している場合は、先頭艦や編隊後備の艦艇は、艦隊防空艦として艦隊の中央に対して相当の距離を置いて配置される。そのような船団や艦隊に対しては、目標運動解析装置１００のオペレータは、入力装置１０３を介して、目標Ａに対する警戒配置距離（前衛又は後衛）を設定する。これから、制約条件算出部２０２は、目標Ｂに対する目標運動パラメータの最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎ等の制約条件を算出する。
以下では、目標Ａ、目標Ｂが共に対潜作戦を行う対潜艦艇である場合を例に説明を行う。
【００３４】
（解析処理部２０３）
次に解析処理部２０３について説明する。
解析処理部２０３は、先ず、目標Ａに対して観測データベース１０２ｂに蓄積された観測データに、制約条件無しで公知の最小二乗法等の最適化手法を用いて、目標運動パラメータを算出する。解析処理部２０３で算出された目標運動パラメータは、表示装置１０４に表示される。
図２は、目標運動解析装置における目標移動体の目標運動パラメータの説明図である。図２において、座標原点には、自艦位置を半円形の符号で示す。縦軸は、方位角で真北をＹ軸の＋符号で示し、横軸は方位角で真東をＸ軸の＋符号で示す。当然、Ｙ軸の−符号は真南を示し、Ｘ軸の−符号は真西を示す。
【００３５】
図２中の各記号の意味は次のとおりである。
＾Ｒ  ：最新観測時刻ｔnにおける自艦から目標までの推定距離（図２中、符号“＾”
を符号“Ｒ”の上に記載して表示）
＾Ｂｙ：最新観測時刻ｔnにおける自艦から見た目標の推定方位角（推定方位）を真北
Ｙからの時計回り角度で表示（図２中、符号“＾”を符号“Ｂ”の上に記載して
表示）
＾ｘ  ：最新観測時刻ｔnにおける自艦を基準とした目標の解析位置のＸ軸成分（図２
中、符号“＾”を符号“ｘ”の上に記載して表示）
＾ｙ  ：最新観測時刻ｔnにおける自艦を基準とした目標の解析位置のＹ軸成分（図２
中、符号“＾”を符号“ｙ”の上に記載して表示）
＾Ｃｔ：最新観測時刻ｔnにおける目標の推定針路（図２中、符号“＾”を符号“Ｃ”
の上に記載して表示）。真北方向を基準に時計回りの角度で表示。
＾Ｍｔ：最新観測時刻ｔnにおける目標の推定速度（図２中、符号“＾”を符
号“Ｍ”の上に記載して表示し、ベクトルであることを下線で表示）
＾νx：目標の推定速度ベクトルのＸ軸成分（図２中、符号“＾”を符号“ν”の上に
記載して表示）
＾νy：目標の推定速度ベクトルのＹ軸成分（図２中、符号“＾”を符号“ν”の上に
記載して表示）
【００３６】
ちなみに、自艦からの艦首方向を中心にして、時計の文字盤表示で１２時の方向を自艦の艦首方向として、目標の方位角や目標の速度ベクトルの方向を示す方法もあるが、潜水艦の場合、例えば、真北をＹ軸とした海底地形図上に自艦の位置や進路、速度を表示し、その同じ座標系上で目標の水上艦艇等の位置や針路等を表示しないと、目標を追尾している自艦（潜水艦）が座礁する危険があり、ここでは、目標移動体の目標運動パラメータ、自艦からの方位、距離、目標の速度、針路も同じ座標系上で表すとしている。
また、符号“＾（ハット）”は、解析結果、又は演算結果を示す記号として示し、後記するように符号“￣（オーバーライン）”は、観測データを示す符号として用いる。
【００３７】
そして、解析処理部２０３は、オペレータが入力装置１０３を介して選択した目標Ｂに対して、制約条件算出部２０２において設定された最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎを用いて、目標Ｂの目標運動パラメータ（自艦からの距離、方位、速度、針路）を算出する。この目標Ｂの目標運動パラメータの算出の方法については、図６のフローチャートの説明の中で詳細に説明する。
【００３８】
（全体フローチャート）
次に、図３から図５を参照しながら複数の目標移動体の目標運動パラメータを算出する全体フローチャートについて説明する。
図３は、一つのグループとして編隊を組んで行動する複数の目標移動体、例えば、艦船の編隊における２隻の艦船の相対位置関係の説明図である。符号４０１を付した目標Ａに対し、符号４０２Ａを付した目標Ｂは、目標Ａと一定の位置関係を保持して編隊運動をする。図３では、目標Ｂは常に目標Ａから距離ｒ以内に位置するように運動していることを示している。ちなみに、目標Ｂは常に目標Ａから距離ｒ以内に位置するということは、符号４０２Ａの位置に限定されなく、符号４０２Ｂを付した仮想線の位置に目標Ｂが位置していることもある。
【００３９】
図４は、目標運動解析装置における、例えば、２隻の目標艦船に対する目標運動パラメータ解析計算の全体の流れを示すフローチャートである。
先ず、ステップＳ５０１では、解析処理部２０３が、目標運動解析装置１００のオペレータの入力装置１０３を介した入力操作により、目標Ａが選定され、記憶部１０２の観測データベース１０２ｂの目標Ａに対応する観測データを読み出して、目標Ａに対する目標運動パラメータ（自艦からの目標Ａの距離及び方位、並びに目標Ａの速度及び針路）の解析値を公知の方法で、解析処理部２０３において算出する（「目標Ａの目標運動パラメータの算出」）。
【００４０】
次に、ステップＳ５０２では、目標運動解析装置１００のオペレータの入力装置１０３を介した入力操作により、目標Ａの目標運動パラメータの解析値に基づく自艦からの目標Ａの位置を基準として、目標Ｂに対する最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎを算出設定する（「目標Ｂに対する制約条件算出」）。
ここで、ステップＳ５０１，Ｓ５０２が、特許請求の範囲に記載の「目標移動体選択工程」に対応する。最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎは、特許請求の範囲に記載の「該２つの特定目標移動体同士の間の相対位置関係に係る目標相対位置関係データ」に対応する。
【００４１】
図５を参照しながらステップＳ５０２における最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎの算出設定方法について詳細に説明する。図５は、目標運動解析装置の制約条件算出部が行う処理の説明図である。符号４０１を付した解析位置（＾ｘA，＾ｙA）の目標Ａに対し、符号４０２を付した目標Ｂは、目標Ａから一定の距離ｒ内に位置する。
また、自艦を符号７０３で示してある。符号７０５を付した解析位置（＾ｘB，＾ｙB）の目標Ｂは、目標Ｂに関し制約条件無しで解析処理部２０３において従来と同様の公知の方法で算出したと想定したものである。符号４０２を付した目標Ｂは、目標Ｂに関し後記するように制約条件を付けて解析処理部２０３において解析された結果の位置を示している。
【００４２】
制約条件算出部２０２は、入力装置１０３から、目標Ａに対する基準目標の指定、目標Ｂの選択指定、目標Ａから目標Ｂまでの相対位置の半径ｒの入力を受け、解析処理部２０３から目標Ａの解析位置（＾ｘA，＾ｙA）と目標Ｂの最新の観測方位のデータを取得する。そして、制約条件算出部２０２は、目標Ａの解析位置（＾ｘA，＾ｙA）を中心とする半径ｒの円７０６を描き、これと観測方位線７０７との交点を求める。一般に、交点は２点求まるが、このうち自艦７０３からの距離が遠い側の交点７０８までの距離を最大距離制限値Ｒｍａｘ、自艦７０３からの距離が近い側の交点７０９までの距離を最小距離制限値Ｒｍｉｎとし、解析処理部２０３へ出力する。
ここで、半径ｒの円７０６は、制約条件算出部２０２が、海域音響特性データベース１０２ｅから、自艦（潜水艦）が現在存在する海域における前記した対潜艦艇からの有効水平方向半径距離のデータから容易に得られる。
【００４３】
具体的に図５に従ってより具体的に説明する。自艦７０３を基準とした目標Ａの位置ベクトルをＲA、目標Ｂの最新の観測方位を￣Ｂｙnとする。目標Ａの推定位置ベクトル＾ＲAは、次式（１）のように目標Ａの推定位置ベクトル＾ＲAのＸ軸成分＾ｘA、Ｙ軸成分＾ｙAにより表される。
【数１】
2014025786.tif 000003
符号４０１で示した目標Ａから観測方位線７０７に下ろした垂線の長さｈは、次式（２）から求まる。
【００４４】
【数２】
2014025786.tif 000004
そして、観測方位線７０７に下ろした垂線の足７１０から、自艦７０３までの距離ｒhを次式（３）から求める。
【数３】
2014025786.tif 000005
【００４５】
半径ｒの円７０６と観測方位線７０７との交点７０８，７０９それぞれから垂線の足７１０までの距離ｋは、次式（４Ａ）により求める。
また、ｒhとｋの関係を用いて、次式（４Ｂ），（４Ｃ）により最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎを求める。
【数４】
2014025786.tif 000006
【００４６】
図４に戻って、ステップＳ５０３では、ステップＳ５０２において設定した制約条件に従って目標Ｂに対する目標運動パラメータ（自艦からの目標Ｂの距離及び方位、並びに目標Ａの速度及び針路）を、解析処理部２０３において算出する（「目標Ｂの目標運動パラメータの算出」）を算出する。
【００４７】
次に、図６、図７を参照し、適宜、図１を参照しながらステップＳ５０３における詳細な処理について説明する。図６は、解析処理部における目標運動パラメータ解析計算の流れを示すフローチャートである。図７は、解析処理部の出力する目標運動パラメータの説明図であり、パッシブセンサ（センサ）１の観測データを用いた目標運動解析の原理を示す図である。
【００４８】
図７に示すようにセンサ信号処理装置２から得られる観測データは、時刻ｔi(ｉ＝１〜ｎ)における音波の到来方位である観測方位￣Ｂｙi（観測方位情報）と観測周波数￣ｆi（観測周波数情報）、及び観測方位の誤差である標準偏差σByi（図示せず）、観測周波数の誤差である標準偏差σfi（図示せず）である。この間に自艦も符号７０３Ａで示した位置（ownｘ1，ownｙ1）から軌跡７０４に示すように符号７０３Ｂで示した位置（ownｘn，ownｙn）にまで移動している。この自艦の各時刻ｔi(ｉ＝０〜ｎ)の運動軌跡と針路は、前記した記憶部１０２の自艦運動データベース１０２ｃに格納されており、自艦からの相対的な真北Ｙ軸を基準線とした時計回りの角度で示す目標Ｂの観測方位、￣Ｂｙiに基づき、目標Ｂが一定の速度で同一方向に航行を続けていると仮定すれば、解析処理部２０３において解析による目標Ｂの推定方位＾Ｂｙi、推定距離＾Ｒは、容易に算出できる。
【００４９】
また、自艦の目標Ｂに対する相対速度が変化することで、ドップラー効果により観測周波数￣ｆiが変化する。このことから、自艦の針路と速度を加味することで、目標Ｂの推定速度＾Ｍｔ、推定針路＾Ｃｔの算出も可能である。
【００５０】
要約すると、目標運動解析装置１００は、これら時系列の観測方位情報、周波数情報に基づき、最新観測時刻ｔnにおける自艦から目標Ｂまでの推定距離＾Ｒ、目標Ｂの推定方位＾Ｂｙ、推定速度＾Ｍｔ、推定針路＾Ｃｔから構成される目標Ｂの状態ベクトル（目標状態ベクトル）＾ＳBを求めるものである（後記する式（８）参照）。
ここで、ベクトルや行列を示す符号には、下線（ ）を付して、明示する。
【００５１】
そして、目標運動解析装置１００では、解析処理部２０３において、前記した目標Ｂの状態ベクトル＾ＳBの値から時刻ｔiの推定方位＾Ｂｙi、推定周波数（「推定音源周波数」とも称する）＾ｆiを求め、同時刻の目標Ｂの観測方位￣Ｂｙi、観測周波数￣ｆiとの差を求める。これを標準偏差σByi、標準偏差σfiで規格化した上で二乗し、観測データ数ｎ個分を合計したものを算出し、更に最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎの制約条件からの目標Ｂの自艦からの推定距離＾Ｒに基づいて算出される逸脱量ｇk（＾ＳB(m)），（ｋ＝１，２）に応じたペナルティを加算した評価関数Фを計算する。
【００５２】
ここで、推定方位＾Ｂｙiと観測方位￣Ｂｙiとの差を標準偏差σByiで規格化した上で二乗し、推定周波数＾ｆiと観測周波数￣ｆiとの差を標準偏差σfiで規格化した上で二乗し、観測データ数ｎ個分を合計したものを算出するとは、次式（５）による。
ここで、解析値を示す符号“＾”は、次式（５）では、符号“ｆi”、符号“Ｂｙi”の上に付してある。以下、同様である。
【００５３】
【数５】
2014025786.tif 000007
解析処理部２０３は、繰返し処理によって、評価関数Фがより小さな値となるように目標Ｂの状態ベクトル＾ＳB値を更新していき、最終的にФを最小とする目標Ｂの状態ベクトル＾ＳBを最良推定値として出力する。
以下、図６のフローチャートにしたがって、解析処理部２０３が、目標Ｂの状態ベクトル＾ＳBを最良推定値として出力する処理について詳細に説明する。
【００５４】
ステップＳ６０１では、観測データベース１０２ｂに格納された目標Ｂの観測データを読み込む。この目標Ｂの観測データは、目標Ａに対する観測データと同じ時間帯のものであっても良いし、目標Ａの針路、速力は、状態ベクトルとして既に算出済みであるため、任意の時刻の目標Ａの位置を外挿計算で算出しても良い。つまり、目標Ｂの観測データの取得時刻が目標Ａの観測データの取得時刻と食い違っていても、目標Ａの状態ベクトルを用いて外挿計算することによって、目標Ｂの観測データの取得時刻における目標Ａの位置を求め、それを基準に前記した制約条件を算出することができる。
目標Ｂの観測データとしては、パッシブセンサ１を構成する艦首部パッシブソナー、コンフォーマルパッシブソナー、艦体側面パッシブソナー、曳航式パッシブソナー等それぞれからの信号を総合的に組み合わせ評価して、目標Ｂに対応する音源の自艦からの方位角（観測方位）￣Ｂｙi（式（６）参照）、観測方位の誤差（標準偏差σByi）、観測周波数￣ｆi（式（７）参照）、観測周波数の算出誤差（標準偏差σfi）を含んでいる。
更に、ステップＳ６０１では、自艦運動データベース１０２ｃから自艦運動情報を読み込む。
【００５５】
【数６】
2014025786.tif 000008
ここで、￣Ｂｙiは、時刻ｔiにおける観測方位を示し、時刻ｔ1から最新時刻ｔn（ｉ＝１〜ｎ）までの観測方位を示す。Ｂｙは、ｎ個の観測方位￣Ｂｙiを並べたベクトル（観測方位ベクトル）である。
【００５６】
【数７】
2014025786.tif 000009
ここで、￣ｆiは、時刻ｔiにおける観測周波数を示し、時刻ｔ1から最新時刻ｔn（ｉ＝１〜ｎ）までの観測周波数を示す。ｆは、ｎ個の観測周波数￣ｆiを並べたベクトル（観測周波数ベクトル）である。
【００５７】
ステップＳ６０２では、目標Ｂの状態ベクトルの繰り返し計算の推定値の初期値＾ＳB(0)（式（８）参照）を設定する（「状態ベクトル初期値計算」）。ここで初期値＾ＳB(0)、及びその各パラメータ＾ｘ(0)，＾ｙ(0)，＾νx(0)，＾νy(0)，＾ｆ(0)には、最新時刻ｔnに対応する添字「ｎ」を付すべきであるが、数式表現が複雑になるのでここでは省略する。
【数８】
2014025786.tif 000010
式（８）では、初期状態として繰り返しの引数ｍ＝０としている。
【００５８】
ステップＳ６０３では、解析処理部２０３は、制約条件算出部２０２で設定された最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎを読み込む（「設定された制約条件の読み込み」）。
ステップＳ６０４では、ステップＳ６０３で読み込んだ制約条件に従って、解析処理の制御を行うパラメータであるラグランジュ乗数ｕ、ペナルティ乗数ｐの初期値を設定する（「ペナルティ乗数初期値設定」）。
【００５９】
ラグランジュ乗数ｕ及びペナルティ乗数ｐは、繰り返し計算の中で必要に応じて更新設定される。具体的には、最大距離制限値Ｒｍａｘに対応して、ｕ1(s)，ｐ1(s)が設定されるとともに、最小距離制限値Ｒｍｉｎに対応してｕ2(s)，ｐ2(s)が設定される。
ここで、引数ｓ＝０のときのｕ1(s),ｕ2(s),ｐ1(s),ｐ2(s)の初期値としては、次式（９）に示すように経験値として設定される。
ここで、制約条件をつけて最適解を得る場合にラグランジュ乗数とペナルティ乗数とを組み合わせて導入する手法は、安定的に制約条件付の最適解を得る方法として公知の技術である。このラグランジュ乗数及びペナルティ乗数は、繰り返し計算によって得られる目標Ｂの状態ベクトル＾ＳB(m)に対する制約条件である最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎに基づいて、後記する式（１８），（１９）によって算出される逸脱量ｇk（＾ＳB(m)），（ｋ＝１，２）に対する後記する判定式（３２Ａ），（３２Ｂ），（３３）の結果に応じて、後記する式（３４）〜（３６）により初期値から更新設定される。
【００６０】
【数９】
2014025786.tif 000011
【００６１】
ステップＳ６０５は、ステップＳ６０２で設定した目標Ｂの状態ベクトルの初期値＾ＳB(0)に対する評価値Φ(0)を計算する。
評価値Φ(m)は、後記する式（２３）で示すように、方位情報評価値φBy(m)、周波数情報評価値φf(m)、制約条件である最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎからの逸脱量ｇk（＾ＳB(m)），（ｋ＝１，２）に対するペナルティφg1(m)，φg2(m)の和として算出される。
以下に、評価値Φ(m)の算出の詳細について説明する。
ここで、次式（１０）〜（２３）では、繰り返し計算回数の引数ｍをｍ＝０とせず、ｍのまま一般化して表示してある。
観測時刻ｔiにおける自艦から見た目標Ｂの推定方位＾Ｂｙi(0)を、目標Ｂの状態ベクトル＾ＳB(0)の値（＾ｘ(0)，＾ｙ(0)，＾νx(0)，＾νy(0)）と自艦の時刻ｔi及び最新時刻ｔnにおける位置情報である座標（ownｘi，ownｙi），（ownｘn，ownｙn）から、次式（１０）によって計算する。
【００６２】
【数１０】
2014025786.tif 000012
ここで、目標Ｂの推定位置ベクトル＾ＲB(m)（＝（＾ｘ(m)，＾ｙ(m)））は、時刻ｔnにおける自艦位置（ownｘn，ownｙn）を基準にした時刻ｔi（ｉ＝１〜ｎ）における目標Ｂの相対位置である。式（１０）を用いて行なっている計算は、時刻ｔiにおける目標Ｂの推定方位＾Ｂｙi(m)を求めるため、時刻ｔiにおける自艦位置（ownｘi，ownｙi）を基準にした目標Ｂの相対位置を求めている。そのため、時刻ｔnにおける目標Ｂの推定位置ベクトル＾ＲB(m)（＝（＾ｘ(m)，＾ｙ(m)））に基づいて、時刻ｔiから時刻ｔnまでの目標Ｂ及び自艦の移動量を割り戻す計算を行なっている。そして、式（１０）で用いられている＾νx(m)，＾νy(m)は、時刻ｔnにおける目標Ｂの推定速度ベクトルのＸ軸成分とＹ軸成分である。
【００６３】
観測方位ベクトル￣Ｂｙ(m)と推定方位ベクトル＾Ｂｙ(m)の差分に対する評価値φBy(m)を、次式（１１）によって計算する。
【数１１】
2014025786.tif 000013
ここで、推定方位＾Ｂｙi(m)をｎ個並べた推定方位ベクトル＾Ｂｙ(m)は、次式（１２）で定義され、観測方位誤差の標準偏差σByiのｉ＝１〜ｎの二乗の逆数の対角行列ＷByは、次式（１３）で定義される。
【００６４】
【数１２】
2014025786.tif 000014
【数１３】
2014025786.tif 000015
【００６５】
時刻ｔiにおける目標Ｂの推定周波数＾ｆi(m)を、目標Ｂの状態ベクトル＾ＳB(m)の値（＾νx(m)，＾νy(m)，＾ｆ(m)）、推定方位＾Ｂｙi(m)、自艦の真北を示すＹ軸から時計回りで示した角度表示の針路  ownＣi、自艦の速度ownＶi及び水中音速Ｖｓから、次式（１４）〜（１６）によって計算する。
ちなみに、水中音速Ｖｓは、記憶部１０２の観測データベース１０２ｂに記憶された当該海域の海水温度から容易に算出される。
【数１４】
2014025786.tif 000016
【００６６】
観測周波数￣ｆiと推定周波数＾ｆi(m)の差分に対する評価値φf(m)を、次式（１７）によって計算する。
【数１５】
2014025786.tif 000017
【００６７】
最大距離制限値Ｒｍａｘと目標Ｂの状態ベクトル＾ＳB(m)の距離成分＾Ｒ(m)との差分量（逸脱量）ｇ1（＾ＳB(m)）、最小距離制限値Ｒｍｉｎと目標Ｂの状態ベクトル＾ＳB(m)の距離成分＾Ｒ(m)との差分量（逸脱量）ｇ2（＾ＳB(m)）を式（１８），（１９），（２０）に基づいて計算する。
【００６８】
【数１６】
2014025786.tif 000018
【００６９】
目標Ｂの状態ベクトル＾ＳB(m)の距離成分＾Ｒ(m)と距離制限値（最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎ）の差分量（逸脱量）ｇk（＾ＳB(m)）に対するペナルティφｇk(m)を次式（２１），（２２）によって計算する。ここで、引数ｋは、１又は２である。
【００７０】
【数１７】
2014025786.tif 000019
【００７１】
目標Ｂの状態ベクトル＾ＳB(m)（繰り返し計算回数ｍ回目の推定値）に対する評価値Φ(m)を、方位情報評価値φBy(m)、周波数情報評価値φf(m)、制約条件である最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎからの逸脱量に対するペナルティφg1(m)，φg2(m)の和として算出する（式（２３）参照）。
【数１８】
2014025786.tif 000020
【００７２】
ステップＳ６０６では、次のステップＳ６０７の処理において目標Ｂの状態ベクトル＾ＳB(m)の更新量算出に用いる、推定方位＾Ｂｙi(m)の微係数、推定周波数＾ｆi(m)の微係数、及び目標Ｂの状態ベクトル＾ＳB(m)の距離成分＾Ｒ(m)の最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎからの逸脱量ｇ1（＾ＳB(m)）,ｇ2（＾ＳB(m)）に対するペナルティの微係数を計算する（「微係数計算」）。
その詳細を次式（２４）〜（２７）で詳細に説明する。
【００７３】
自艦からみた目標Ｂの推定方位＾Ｂｙi(m)を計算するのに使用した式（１０）を、目標Ｂの状態ベクトル＾ＳB(m)の各要素で偏微分した推定方位＾Ｂｙi(m)の微係数である微係数行列Ｈ(m)を、次式（２４）に基づき求める。
【数１９】
2014025786.tif 000021
【００７４】
次に、自艦からみた推定周波数＾ｆi(m)を計算するのに使用した式（１４）を、目標Ｂの状態ベクトル＾ＳB(m)の各要素で偏微分した推定周波数＾ｆi(m)の微係数である微係数行列Ｋ(m)を、次式（２５）に基づき求める。
【００７５】
【数２０】
2014025786.tif 000022
【００７６】
更に、目標Ｂの状態ベクトル＾ＳB(m)の距離成分＾Ｒ(m)と、最大距離制限値Ｒｍａｘ及び最小距離制限値Ｒｍｉｎとの差分量（逸脱量）ｇ1（＾ＳB(m)），ｇ2（＾ＳB(m)）を、目標Ｂの状態ベクトル＾ＳB(m)の各要素で偏微分した微係数の行列Ｇ1(m)，Ｇ2(m)を、次式（２６），（２７）に基づき求める。
【００７７】
【数２１】
2014025786.tif 000023
【００７８】
ステップＳ６０７では、ステップＳ６０５または繰り返し計算の前回のステップＳ６０８において目標Ｂの状態ベクトル＾ＳB(m)に対する評価値Φ(m)（ｍ＝０〜繰り返し計算が終了する最後引数）を計算する過程で求めた推定方位＾Ｂｙi(m)、推定周波数＾ｆi(m)、及びステップＳ６０４において設定又は後記するステップＳ６１８で更新されたペナルティ乗数と、ステップＳ６０６で求めた推定方位＾Ｂｙi(m)の微係数、推定周波数＾ｆi(m)の微係数、及び目標Ｂの状態ベクトル＾ＳB(m)の距離成分＾Ｒ(m)の最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎからの逸脱量に基づいて算出される逸脱量Ｄ(m)対するペナルティの微係数とから、次式（２８）によって目標Ｂの状態ベクトル＾ＳB(m)の更新量ΔＳBを求める。
【００７９】
【数２２】
2014025786.tif 000024
【００８０】
ここで、Ｉは、単位行列であり、は、目標Ｂの状態ベクトル＾ＳBの探索状態に応じ更新量ΔＳBのスケールを制御するための縮小因子である。
【００８１】
また、式（２８）における第１項の逆行列の部分は、２階微分の項に相当し、それに乗じられる第２項の部分は、１階微分の項に相当しており、最適解を得るためのニュートン法を改良したレーベンバーグ・マーカート法（修正ニュートン法）として知られている。
ニュートン法では、繰り返し計算の中で最適解に接近する途中で、突然最適解から遠い解にジャンプする場合があるので、ここでは安定的に最適解を得るため、レーベンバーグ・マーカート法を用いている。
そして、ステップＳ６０７では、更に、目標Ｂの状態ベクトル＾ＳB(m)の更新量ΔＳBを用いて仮更新の目標Ｂの状態ベクトル＾ＳB new(m)を次式（２９）に示すように算出する（「状態ベクトル更新量計算」）。
【００８２】
【数２３】
2014025786.tif 000025
【００８３】
ステップＳ６０８では、ステップＳ６０７において求めた仮更新の目標Ｂの状態ベクトル＾ＳB new(m)に対し、ステップＳ６０５と同様の処理を行い、評価値Φ(m)newを計算する。
【００８４】
ステップＳ６０９では、ステップＳ６０８において算出された仮更新の目標Ｂの状態ベクトル＾ＳB new(m)に対する評価値Φ(m)newを、仮更新前の目標Ｂの状態ベクトル＾ＳB(m)に対する評価値Φ(m)と比較する。そして、（評価値Φ(m)＞（評価値Φ(m)new）か否かを判定する（「（前回Φ）＞（今回Φ）？」）。評価値Φ(m)の値が評価値Φ(m)newより大の場合（Ｙｅｓ）は、ステップＳ６１０へ進み、そうでない場合（Ｎｏ）は、ステップＳ６１１へ進む。
【００８５】
ステップＳ６０９においてＹｅｓでステップＳ６１０に進むと、縮小因子λをそれまでの０．５倍とし、ステップＳ６１２へ進む。
ステップＳ６０９においてＮｏでステップＳ６１１に進むと、縮小因子λをそれまでの２倍にして、ステップＳ６０７に戻り、ステップＳ６０７の状態ベクトル更新量計算を繰り返す。ここで、ステップＳ６１０，Ｓ６１１において縮小因子λの値を変更するのは、確実に目標Ｂの状態ベクトル＾ＳB(m)を最適解に収束させるためであり、縮小因子λの値を最初は大きくし、繰り返し計算の中で、ステップＳ６０９でＹｅｓの場合は、次々に縮小因子λの値を小さくし、安定的に目標Ｂの状態ベクトル＾ＳB(m)の最適解を得るようにする。
【００８６】
続いて、ステップＳ６１２では、ステップＳ６０７において算出した目標Ｂの状態ベクトル＾ＳB(m)の更新量ΔＳBによって目標Ｂの状態ベクトル＾ＳB(m+1)に更新する（式（３０）参照）。
その際に、これに対応する評価値Φ(m)newもΦ(m+1)として解析結果データベース１０２ｄ（図１参照）に記憶させる（式（３１）参照）。
【数２４】
2014025786.tif 000026
【００８７】
ステップＳ６１３では、ステップＳ６０７において算出した目標Ｂの状態ベクトル＾ＳB(m)の更新量ΔＳBが予め設定された規定値以下か否かを判定する（「状態ベクトル更新量が規定値以下？」）。更新量ΔＳBが規定値以下の場合（Ｙｅｓ）は、ステップＳ６１４へ進み、そうでない場合（Ｎｏ）は、ステップＳ６０６へ戻り、再度計算を行う。
ここで、更新量ΔＳBが規定値以下か否かの判定は、更新量ΔＳBのノルムが予め設定された規定値以下か否かで判定する。
【００８８】
ステップＳ６１４では、ステップＳ６１２において更新された目標Ｂの状態ベクトル＾ＳB(m+1)の距離成分＾Ｒ(m+1)が最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎの制約条件を満足するか否かを判定する（「状態ベクトルが制約条件を満足？」）。目標Ｂの状態ベクトル＾ＳB(m+1)の距離成分＾Ｒ(m+1)が、最大距離制限値Ｒｍａｘ以下、最小距離制限値Ｒｍｉｎ以上の条件を満足している場合（Ｙｅｓ）は、この繰り返し計算回数の段階の目標Ｂの状態ベクトル＾ＳB(m+1)を解析処理部２０３における最適解析値として採用して繰り返し計算を終了する。そして、ステップＳ５０１で算出された目標Ａの目標運動パラメータ並びに目標Ｂの状態ベクトル＾ＳB(m+1)のデータ（目標Ｂの目標運動パラメータ）を、表示制御部２０７（図１参照）を介して、表示装置１０４に出力表示する（特許請求の範囲に記載の「目標運動パラメータ出力工程」に対応する）。ここで、繰り返し計算回数ｍ＋１は、ステップＳ６１３，Ｓ６１４の判定を共に満足したときの繰り返し計算回数ｍの最後の引数値である。
目標Ｂの状態ベクトル＾ＳB(m+1)の距離成分＾Ｒ(m+1)が、最大距離制限値Ｒｍａｘ以下、最小距離制限値Ｒｍｉｎ以上の条件を満足していない場合（Ｎｏ）は、ステップＳ６１５へ進む。
【００８９】
ステップＳ６１５では、ステップＳ６０７において更新された目標Ｂの状態ベクトル＾ＳB(m+1)の距離成分＾Ｒ(m+1)の制約条件（最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎ）との差分量（逸脱量）ｇ1（＾ＳB(m+1)），ｇ2（＾ＳB(m+1)）を計算する（「制約条件逸脱量計算」）。
この制約条件（最大距離制限値Ｒｍａｘ、最小距離制限値Ｒｍｉｎ）との差分量（逸脱量）ｇ1（＾ＳB(m+1)），ｇ2（＾ＳB(m+1)）の計算は、式（１８），（１９），（２０）に基づいて、繰り返し計算回数の引数ｍをｍ＋１に読み替えて行う。次いで、ステップＳ６１５では、最大逸脱量gmax、最大逸脱量gmax,k（ｋ＝１，２）を、次式（３２Ａ），（３２Ｂ）により計算する。
【００９０】
式（３２Ａ）に従って、最大距離制限値Ｒｍａｘからの差分量（逸脱量）ｇ1（＾ＳB(m+1)）と、ラグランジュ乗数ｕ1(s)とペナルティ乗数ｐ1(s)の比（正値）との内の第１の最大値の絶対値、小距離制限値Ｒｍｉｎからの差分量（逸脱量）ｇ2（＾ＳB(m+1)）と、ラグランジュ乗数ｕ2(s)とペナルティ乗数ｐ2(s)の比（正値）との内の第２の最大値の絶対値とを比較し、第１の最大値の絶対値及び第２の最大値の絶対値の内の大きい方の値を逸脱量ｇmaxとして算出する。最大逸脱量gmax,k  （ｋ＝１，２）については、式（３２Ｂ）に従って算出されるが、ほぼ同様であり重複する説明を省略する。
ここで、ラグランジュ乗数ｕ1(s)，ｕ2(s)、ペナルティ乗数ｐ1(s)，ｐ2(s)の添字ｓは、この繰り返し計算回数ｍにおけるラグランジュ乗数ｕ1(s)，ｕ2(s)、ペナルティ乗数ｐ1(s)，ｐ2(s)それぞれの更新回数を示す引数である。
ステップＳ６１５で算出された最大逸脱量ｇmaxの値は、最大逸脱量記録値Ｄ(s+1)の候補値して保持される。
【数２５】
2014025786.tif 000027
【００９１】
そして、ステップＳ６１５では、今回の最大逸脱量ｇmax,k（ｋ＝１，２）が、前回の最大逸脱量ｇmaxの値である最大逸脱量記録値Ｄ(s)に所定値βを乗じた結果の値を超えているか否かを次式（３３）によって判定する。制約条件を逸脱している場合は、その式番号ｋ（ｋ＝１or２）を解析結果データベース１０２ｄに記憶させるとともに更新回数ｓにおいて逸脱していることを示すフラグを解析結果データベース１０２ｄに記憶させる。
【数２６】
2014025786.tif 000028
ここで、Ｄ(s)は、前回の更新回数ｓにおける記憶された最大逸脱量記録値Ｄ(s)を意味し、βは逸脱判定許容係数である。逸脱判定許容係数βの値としては、経験値として例えば、０．２５の値を用いる。
なお、最大逸脱量記録値Ｄ(s)のｓ＝０のときの初期値としては、大きな値である所定の経験値が設定される。
【００９２】
ステップＳ６１６では、ステップＳ６１５において算出された最大逸脱量記録値の候補値Ｄ(s+1)が、前回の最大逸脱量記録値Ｄ(s)から減少しているか否かを判定する（「今回逸脱量＜前回逸脱量？」）。最大逸脱量記録値の候補値Ｄ(s+1)が、前回の最大逸脱量記録値Ｄ(s)から減少している場合（Ｙｅｓ）は、最大逸脱量記録値の候補値Ｄ(s+1)を正式の最大逸脱量記録値Ｄ(s+1)として解析結果データベース１０２ｄに記憶させて、ステップＳ６１７へ進む。そうでない場合（Ｎｏ）は、前回の最大逸脱量記録値Ｄ(s)を今回の最大逸脱量記録値Ｄ(s+1)として解析結果データベース１０２ｄに記憶させて、ステップＳ６１８へ進む。
【００９３】
ステップＳ６１７では、ラグランジュ乗数ｕ1(s+1)，を次式（３４），（３５）に従って更新する（「ラグランジュ乗数更新」）。
【数２７】
2014025786.tif 000029
【００９４】
ステップＳ６１６においてＮｏでステップＳ６１８へ進む場合は、今回の繰り返し計算回数ｍにおいて使用したラグランジュ乗数ｕk(s)  （ｋ＝１，２）の値を、次回の繰り返し計算回数ｍ＋１におけるラグランジュ乗数の更新回数の引数（ｓ+１）のラグランジュ乗数ｕk(s+1)  （ｋ＝１，２）の値として用いる。
【００９５】
ステップＳ６１８では、ステップＳ６１５において式（３３）の判定により、制約条件を逸脱している式番号ｋ（ｋ＝１or２）のペナルティ乗数ｐk(s)について、次式（３６）に従って、ペナルティ乗数ｐk(s)をｐk(s+1)に更新する（「ペナルティ乗数更新」）。ステップＳ６１５において式（３３）の判定により、制約条件を逸脱していない式番号ｋ（ｋ＝１or２）のペナルティ乗数ｐk(s)については、前回のペナルティ乗数ｐk(s)の値をそのまま用いてペナルティ乗数ｐk(s+1)として更新する（「ペナルティ乗数更新」）。
【００９６】
【数２８】
2014025786.tif 000030
ここで、ペナルティ乗数ｐk(s)  （ｋ＝１ｏｒ２）に乗じられる定数αは、経験的に設定されたペナルティ乗数更新係数である。
ステップＳ６１８の後、ステップＳ６０６へ戻る。
以上で、一連の目標Ｂの目標運動パラメータの算出の処理が完了する。
ここで、フローチャートのステップＳ６０２，Ｓ６１２が特許請求の範囲に記載の「目標状態ベクトル推定演算工程」に対応し、ステップＳ６０５，Ｓ６０８が、特許請求の範囲に記載の「評価値算出工程」に対応し、ステップＳ６０９，Ｓ６１３，Ｓ６１４が、特許請求の範囲に記載の「目標状態ベクトル収束判定工程」に対応する。
【００９７】
以上のステップＳ６０５〜Ｓ６１８の繰り返し計算を、ステップＳ６１３，Ｓ６１４の判定を共に満足するまで繰り返すことにより、目標Ｂの状態ベクトル＾ＳB(m+1)の距離成分＾Ｒ(m+1)が最大距離制限値Ｒｍａｘと最小距離制限値Ｒｍｉｎの間に収まる目標Ｂの状態ベクトル＾ＳB(m)を精度良く求めることができる。ここで、繰り返し計算回数ｍ＋１は、ステップＳ６１３，Ｓ６１４の判定を共に満足したときの繰り返し計算回数ｍの最後の引数値である。
特に、潜水艦による船団や艦隊に対する魚雷襲撃行動においては、船団や艦隊の針路に対して斜め前方位置に占位して、目標の針路と速度と自艦の発射する魚雷の速度（雷速）を考慮して目標の未来位置を予測し、その予測した未来位置に向けて魚雷を発射したとき、その魚雷にとって目標の面積が大きくなる目標の側面（舷側）に突入するように魚雷襲撃行動を行うのが通常である。
そのため、図５に示すように自艦から見ると目標Ａの位置に対して、目標Ｂは目標Ａの進行方向の斜め前後に一定の距離をおいて位置することが多いのが普通であり、前記したように目標Ａを中心に半径ｒの円を先験情報に基づいて設定し、最大距離制限値Ｒｍａｘと最小距離制限値Ｒｍｉｎを設定して、目標Ｂの状態ベクトル＾ＳB(m)の収束を図ることは、妥当である。
【００９８】
本実施形態によれば、最大脅威の目標である目標Ａの目標運動パラメータを目標運動解析装置１００により取得した際に、目標Ａの近くに潜在的に存在する可能性のある脅威度が目標Ａに次いで高い目標Ｂに対して、目標Ａとの相対位置関係情報が正確に把握された状態で目標Ｂの目標運動パラメータを取得することができる。
その結果、自艦である潜水艦から複数の水上艦艇に対して、例えば、魚雷による襲撃行動を行う場合に、目標Ａだけにとらわれることなく目標Ｂに対しても、魚雷攻撃が行える範囲か否かの判断が容易にでき、魚雷の発射時点において目標Ａ，Ｂに魚雷が達する時点における未来位置を、複数の発射魚雷に個別に精度良く入力することができる。
【００９９】
現在の魚雷は、少なくとも発射時に入力された未来位置近傍に達した時点で、目標が発する音響を追尾したり、魚雷頭部のアクティブソナーセンサにより目標を追尾したりしていくホーミング式の誘導方式が採用されている。目標Ａ，Ｂのそれぞれの目標運動パラメータを正確に取得することができるということは、目標Ａ，Ｂの未来位置近傍に魚雷が達した時点において、目標Ａ向けに発射された魚雷が目標Ｂを追尾するということや、その逆に目標Ｂ向けに発射された魚雷が目標Ａを追尾するという個別の目標向けに発射された複数魚雷が、結果的に同一目標を追尾するということを防止できる。
また、図５に示すように目標Ｂの状態ベクトル＾ＳB(m)を得る際に制約条件をつけないとすると、図５の符号７０５で示した位置に目標Ｂが現在あると算出してしまう可能性がある。そのときは、目標Ｂ（７０５）の未来位置に対して魚雷を発射しても、魚雷は目標Ｂ（７０５）の未来位置で目標Ｂに対して魚雷頭部のセンサが目標Ｂを捉えることができず、目標を自動追尾できない迷走魚雷となって無駄になってしまう。
【０１００】
また、特許文献１に記載の従来技術では、個別の目標ごとに、目標の距離、速度、針路のうちの少なくとも１つの目標運動パラメータを先験情報に基づいてオペレータ入力するとしているので、オペレータ入力された目標運動パラメータに大きな誤差がある場合、目標移動体の距離、速度及び針路の目標運動パラメータの内の残りの二つの目標運動パラメータの推定精度が低下したり、オペレータ入力された目標運動パラメータの補正を行う工程を伴うために、精度の高い目標移動体の距離、方位、移動速度及び針路の目標運動パラメータを得るのに時間を要したりする問題がある。
本実施形態によれば、複数の目標に対し、短時間に的確な個別目標の未来位置を魚雷に入力し、魚雷を的確に個別の目標に対し最終的な自動追尾誘導状態に投入することができる。
【０１０１】
本実施形態は、自艦が潜水艦であり、魚雷により複数の水上艦艇等に対する襲撃行動を行う場合に限定されるものではなく、潜水艦から水中発射できる対艦ミサイルによる複数の水上艦艇等に対する襲撃行動の場合にも、発射時に目標Ａ，Ｂに対艦ミサイルが達する時点における未来位置を、複数の対艦ミサイルに個別に精度良く入力することができる。対艦ミサイルは、入力された未来位置近傍において、ミサイルの頭部に設けられたシーカー（レーダ等のアクティブセンサの場合もあれば、画像センサ等のパッシブセンサの場合もある）により目標を捕捉して、目標に最終突入するので、発射された複数のミサイルを目標Ａ，Ｂに対して個別に最終突入させることができる。
【０１０２】
本実施形態では、水上艦艇又は船舶である目標Ａに対して一定の相対距離をおいて水上艦艇又は船舶である目標Ｂが存在するとして説明したが、目標Ａの目標Ｂと相対的に反対側に一定の相対距離をおいて目標Ｃが存在する場合にも、つまり、護衛している船団又は艦隊中央の外周部に対潜艦艇が多数対潜陣形で配置されている場合にも、適用でき、オペレータによって選択される目標の数は二つに限定されたものではない。
ちなみに、自艦にとって脅威の大きい対潜艦艇（目標Ａ、目標Ｂ）に対して攻撃を完了後、船団や艦隊の対潜護衛を失った側から、船団又は艦隊の中央部分の船舶や艦艇に対する攻撃に移行しやすいことは当然である。
【０１０３】
（他の実施形態）
本実施形態では、潜水艦による水上艦艇を目標とした複数目標移動体の目標運動パラメータを、パッシブソナーを用いて解析する場合を例に説明したがそれに限定されるものではない。例えば、レーダに対してステルス性の高い艦艇において、複数のアンテナ素子で構成されたパッシブレーダ装置（デジタル式逆探装置）で、目標である水上艦の発する電波の方位及び電波のドップラー効果による受信電波の周波数の変化を捉えることにより、音響を用いた場合と同様に複数の目標間の相対位置の距離情報を制約条件として用いることで、より精度の高い目標運動パラメータを算出することができる。
【０１０４】
更に、ステルス航空機においてパッシブレーダを用いて目標の編隊を捉え、その複数の航空機を目標として長射程の空対空ミサイルを複数発射するために、複数の航空機目標に対する目標運動パラメータの算出にも適用可能である。
【符号の説明】
【０１０５】
１      パッシブセンサ（センサ）
２      センサ信号処理装置
５      航海センサ
６      位置演算部
７      水温センサ
１００  目標運動解析装置
１０１  ＣＰＵ
１０２  記憶部
１０２ａ  プログラム
１０２ｂ  観測データベース
１０２ｃ  自艦運動データベース
１０２ｄ  解析結果データベース
１０２ｅ  海域音響特性データベース
１０３  入力装置
１０４  表示装置
２０１  データ取得部
２０２  制約条件算出部
２０３  解析処理部
２０５  入力制御部
２０７  表示制御部
４０１  目標Ａ（目標移動体）
４０２，７０５  目標Ｂ（目標移動体）
【図１】
2014025786.tif 000031
【図２】
2014025786.tif 000032
【図３】
2014025786.tif 000033
【図４】
2014025786.tif 000034
【図５】
2014025786.tif 000035
【図６】
2014025786.tif 000036
【図７】
2014025786.tif 000037