(57)【要約】【課題】  現実に沿わない仮定、条件等に起因する音源位置の推定誤差をなくし、実際の条件等に基づいて音源の位置を推定できる方法等を得る。【解決手段】  推定海域において複数の仮想音源を配置したときの、海域の環境に応じた、仮想音源から受波地点までの複数の経路による音波の各到来時刻の理論値に基づく値を各仮想音源に対して到来時間差Ｔａｂｌｅ作成器３が算出し、観測した音波に基づいて、マルチパス波到来時間差測定器２が音源からの音波に関する値を算出し、音源からの音波に関する値と理論値に基づく値からコスト関数算出器４が算出したコストにより、音源位置推定器５が音源からの音波による到来時刻と一致又は最も近い到来時刻の仮想音源を決定し、仮想音源の配置位置を音源の位置と推定する。あるいは、決定した仮想音源及びその近傍の仮想音源の配置位置とこれらのコストから音源位置を推定する。【選択図】        図１
2006194627.tif 000002
【特許請求の範囲】
【請求項１】
音源の位置を推定する領域において複数の仮想音源を仮想的に配置したときの、前記領域の環境に応じた、仮想音源から受波地点までの複数の伝搬経路による音波の各到来時刻の理論値に基づく値を各仮想音源に対して算出する工程と、
前記受波地点において観測した音波に基づいて、推定対象の音源からの音波に関する値を算出する工程と、
前記推定対象の音源からの音波に関する値と前記理論値に基づく値から、前記推定対象の音源からの音波による到来時刻と一致又は最も近い到来時刻の仮想音源を決定し、該仮想音源の配置位置を前記推定対象の音源の位置と推定する工程と
を有することを特徴とする音源位置推定方法。
【請求項２】
音源の位置を推定する領域において複数の仮想音源を仮想的に配置したときの、前記領域の環境に応じた、仮想音源から受波地点までの複数の伝搬経路による音波の到来時刻の理論値に基づいて、最初の音波が到来する時刻と他の伝搬経路による音波が到来する時刻との差の理論値を各仮想音源に対して算出する工程と、
前記受波地点において観測した音波に基づいて、推定対象の音源からの最初の音波が到来した時刻と他の経路により到来した時刻との差の観測値を算出する工程と、
前記推定対象の音源による前記差の観測値と前記各仮想音源による前記差の理論値とにより算出した一致度に基づいて、最も高い一致度が算出された前記仮想音源の配置位置を前記推定対象の音源の位置と推定する工程と
を有することを特徴とする音源位置推定方法。
【請求項３】
前記一致度を算出するために、前記推定対象の音源による前記差の観測値と前記各仮想音源による前記差の理論値の組み合わせを決定する際に、
前記差の値が互いに最も近いものを１対１で対応させて組み合わせを決定することを特徴とする請求項１記載の音源位置推定方法。
前記差の観測値と前記差の理論値とが互いに最も近いものを１対１で対応させて組み合わせを決定することを特徴とする請求項２記載の音源位置推定方法。
【請求項４】
前記一致度を算出するために、前記推定対象の音源による前記差の観測値と前記各仮想音源による前記差の理論値の組み合わせを決定する際に、
前記実際の音源による前記差の観測値に対して、前記仮想音源による前記差の理論値で最も近い値をそれぞれ対応させ、
前記差の理論値において組み合わせが決定しなかった値について、前記差の観測値の中で最も近い値を対応させて、組み合わせを決定することを特徴とする請求項２記載の音源位置推定方法。
【請求項５】
音源の位置を推定する領域において複数の仮想音源を仮想的に配置したときの、前記領域の環境に応じた、各仮想音源から受波地点までの音波の伝搬経路及びそれぞれの前記伝搬経路による音波の到来時刻の理論値を求め、最初の音波が到来する時刻と他の経路による音波が到来する時刻との差の理論値を各仮想音源に対して算出する工程と、
観測した音波に基づいて、推定対象の音源から発せられる音波に基づく値を時系列に算出する工程と、
前記推定対象の音源から発せられる音波が最初に到来する時刻を判断し、該時刻と前記仮想音源における最初の音波が到来する時刻とを一致させたときの、前記差の理論値における前記推定対象の音源から発せられる音波に基づく値から、前記仮想音源と前記推定対象の音源との一致度を算出し、最も高い一致度が算出された前記仮想音源の位置を前記推定対象の音源の位置として推定する工程と
を有することを特徴とする音源位置推定方法。
【請求項６】
前記推定対象の音源から発せられる音波が最初に到来する時刻としてとり得る時間幅を設定し、前記時間幅内の複数の時刻と前記仮想音源に対して最初の音波が到来する時刻とを一致させたときの前記一致度を前記各仮想音源について算出することを特徴とする請求項５記載の音源位置推定方法。
【請求項７】
前記推定対象の音源から発せられる音波に基づく値は、推定対象の音源から発せられる音波とあらかじめ設定されたレプリカ用音波との相関関係を表す値とし、前記各仮想音源の前記差の理論値に対応する時刻における前記相関関係を表す値の総和又は平均に基づく値を、前記各仮想音源の前記一致度とすることを特徴とする請求項５又は６記載の音源位置推定方法。
【請求項８】
複数の音響センサにより観測された、前記推定対象の音源から発せられる音波に基づいて前記推定対象の音源の位置を推定することを特徴とする請求項１〜７のいずれかに記載の音源位置推定方法。
【請求項９】
音源の位置を推定する領域において複数の仮想音源を仮想的に配置したときの、前記領域の環境に応じた、仮想音源から受波地点までの複数の伝搬経路による音波の到来時刻の理論値に基づいて、最初の音波が到来する時刻と他の伝搬経路による音波が到来する時刻との差の理論値を各仮想音源に対して算出し、データとして記憶する到来時間差理論値算出手段と、
前記受波地点において観測した音波のデータに基づいて、推定対象の音源からの最初の音波が到来した時刻と他の経路により到来した時刻との差の観測値を算出する到来時間差観測値算出手段と、
前記推定対象の音源による前記差の観測値と前記各仮想音源による前記差の理論値とにより算出した一致度をコストとして算出するコスト関数算出手段と、
該コスト関数算出手段が算出した前記コストに基づいて最も高い一致度が算出された前記仮想音源の配置位置を前記推定対象の音源の位置と推定する音源位置推定手段と
を備えたことを特徴とする音源位置推定装置。
【請求項１０】
音源の位置を推定する領域において複数の仮想音源を仮想的に配置したときの、仮想音源から受波地点までの複数の伝搬経路による音波の到来時刻の理論値に基づいて、最初の音波が到来する時刻と他の伝搬経路による音波が到来する時刻との差の理論値を各仮想音源に対して算出し、データとして記憶する到来時間差理論値算出手段と、
観測した音波に基づいて、推定対象の音源から発せられる音波に基づく値のデータを時系列に算出するデータ算出手段と、
前記推定対象の音源から発せられる音波が最初に到来する時刻を推定する到来時刻推定手段と、
前記推定対象の音源から発せられる音波が最初に到来する時刻と各前記仮想音源に対して最初の音波が到来する時刻とに基づいて、前記差の理論値における前記推定対象の音源から発せられる音波に基づく値を判断する処理手段と、
該処理手段が判断した前記推定対象の音源から発せられる音波に基づく値に基づいて、前記推定対象の音源と各前記仮想音源との一致度をコストとして算出するコスト関数算出手段と、
該コスト関数算出手段が算出した前記コストに基づいて最も高い一致度が算出された前記仮想音源の位置を前記推定対象の音源の位置と推定する音源位置推定手段と
を備えたことを特徴とする音源位置推定装置。
【請求項１１】
音源位置推定手段は、最も高い一致度が算出された前記仮想音源の位置を前記推定対象の音源の位置と推定する代わりに、前記最も高い一致度が算出された仮想音源及び該仮想音源の近傍に位置する仮想音源について算出したコストに基づいて、前記推定対象の音源の位置を推定することを特徴とする請求項９又は１０記載の音源位置推定装置。
【請求項１２】
観測した音波に基づく数値のデータを含む信号に変換する音響センサと、
請求項９、１０又は１１に記載の音源位置推定装置とを
備えたことを特徴とするソーナー。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、例えば、水中でパルスを放射する音源の位置を推定する音源位置推定方法及び装置、逆探ソーナー等を含むソーナーに関するものである。
【背景技術】
【０００２】
例えば、海中で発せられる音の音源を推定するための方法がある。この方法を行うために逆探ソーナーが用いられる。逆探ソーナーとは、音源から発せられる音（例えばアクティブソーナーが送信する音。以下、探信音という）を信号として受信し、音源の位置を推定するためのパッシブソーナーのことである（非特許文献１参照）。ここで逆探ソーナーは、探信音を受信し、サンプリング等を行って、時系列にディジタルデータ化する音響センサシステムとそのデジタルデータを処理し、音源位置推定を行う音源位置推定装置を備えているものとする。
【０００３】
海中において、音源からの探信音（以下、パルスと呼ぶ）の波（音波）は、例えば海面、海底等での反射等により、直接波以外にも複数の伝搬経路（マルチパスと呼ぶ）を伝搬して、時間差を持って受波センサに到来する。例えば、逆探ソーナーにおいて音源を位置推定する場合、それぞれの経路を伝搬して到来したパルス（マルチパス波と呼ぶ）の到来時間差を観測し、下記の２つの仮定及び１つの条件の下に、音源位置を解析的に推定するという方法がある。
［仮定］
（１）音波は直進する
（２）最初に到来するパルスを直接波とし、２番目に到来するパルスを海底１回反射波とする
［条件］
（１）音源深度及び受波点の深度は既知である
【非特許文献１】海洋音響用語辞典、海洋音響学会編、ｐ．４９、（１９９９）
【非特許文献２】海洋音響−基礎と応用−、海洋音響学会編、ｐ．１２１−１２２、（１９８４）
【発明の開示】
【発明が解決しようとする課題】
【０００４】
しかしながら、海水温、水圧等の影響で、海中内において音速は空間的に変化するため、音波はスネルの法則に従って屈折しながら伝搬し、受波センサに到来する。つまり、海中では音波は直進しない。さらに、音波が屈折しながら伝搬することにより、直接波が受波点（受波センサ）に到達しない場合がある（つまり、受波センサに最初に到来するパルスが直接波でない可能性がある）。これらの問題があるにも関わらず、従来は、上記の仮定（１）、（２）及び条件（１）の下、パルスの到来時間差を経路長に換算し、音源位置を推定していた。
【０００５】
この場合、実際には音波は直進しないので、仮定（１）、（２）が成り立たず、音源の位置推定に誤差が生じる。また、直接波が到達しない場合、最初に到来した海底Ｍ回反射波を直接波とし、次の反射波を海底１回反射波として音源位置を推定してしまうため、音源位置の推定に生じる誤差はさらに大きくなる。さらに、音源深度が未知の場合には、条件（１）を満たさないため、音源位置の推定が困難となってしまう。
【０００６】
そこで、現実に沿わない仮定、条件等に起因する音源位置の推定誤差をなくし、より実際の条件等に基づいて音源の位置を推定できる方法、装置等の実現が望まれていた。
【課題を解決するための手段】
【０００７】
本発明に係る音源位置推定方法は、音源の位置を推定する海域において複数の仮想音源を仮想的に配置したときの、海域の環境に応じた、仮想音源から受波地点までの複数の伝搬経路による音波の各到来時刻の理論値に基づく値を各仮想音源に対して算出する工程と、受波地点において観測した音波に基づいて、推定対象の音源からの音波に関する値を算出する工程と、推定対象の音源からの音波に関する値と理論値に基づく値から、推定対象の音源からの音波による到来時刻と一致又は最も近い到来時刻の仮想音源を決定し、仮想音源の配置位置を推定対象の音源の位置と推定する工程とを有するものである。
本発明においては、音源の位置を推定する海域において複数の仮想音源を仮想的に設定、配置をして、海域の環境に応じた複数の音波の伝搬経路、その経路に基づく受波地点に音波が到来する時刻の理論値に基づく値を算出する。一方、観測によって得られる推定対象の音源からの音波に関する値（例えば、各経路の到来時刻、あらかじめ定めたレプリカの音波との相関関係に基づく値、時系列での音波の大きさ等）を算出する。そして、推定対象の音源からの音波に関する値と理論値に基づく値から、推定対象の音源からの音波による到来時刻と一致又は最も近い到来時刻の仮想音源を決定する。この仮想音源の配置位置を推定対象の音源の位置と推定することにより、あらかじめ設定した仮定、条件の下での音源位置推定ではなく、より現実に沿った音源位置推定を行って、精度を高める。
【０００８】
また、本発明に係る音源位置推定装置は、音源の位置を推定する海域において複数の仮想音源を仮想的に配置したときの、海域の環境に応じた、仮想音源から受波地点までの複数の伝搬経路による音波の到来時刻の理論値に基づいて、最初の音波が到来する時刻と他の伝搬経路による音波が到来する時刻との差の理論値を各仮想音源に対して算出し、データとして記憶する到来時間差理論値算出手段と、受波地点において観測した音波のデータに基づいて、推定対象の音源からの最初の音波が到来した時刻と他の経路により到来した時刻との差の観測値を算出する到来時間差観測値算出手段と、推定対象の音源による差の観測値と各仮想音源による差の理論値とにより算出した一致度をコストとして算出するコスト関数算出手段と、コスト関数算出手段が算出したコストに基づいて最も高い一致度が算出された仮想音源の配置位置を推定対象の音源の位置と推定する音源位置推定手段とを備えたものである。
本発明においては、到来時間差理論値算出手段が、音源の位置を推定する海域において複数の仮想音源を仮想的に配置したときの、海域の環境に応じた、仮想音源から受波地点までの複数の伝搬経路による音波の到来時刻の理論値に基づいて、最初の音波が到来する時刻と他の伝搬経路による音波が到来する時刻との差の理論値を各仮想音源に対して算出し、データとして記憶する。一方、到来時間差観測値算出手段が、推定対象の音源からの最初の音波が到来した時刻と他の経路により到来した時刻との差の観測値を算出する。そして、コスト関数算出手段が、差の観測値と各仮想音源による差の理論値とに基づいて、各仮想音源と推定対象の音源との一致度をコストとして算出し、音源位置推定手段がその一致度が最も高い仮想音源の配置位置を推定対象の音源の位置と推定することで、より現実に沿った音源位置推定を行い、精度を高める。
【発明の効果】
【０００９】
本発明によれば、音源の位置を推定する海域において複数の仮想音源を仮想的に配置して、環境に応じた仮想音源による音波の到来時刻の理論値に基づく値を算出しておき、推定対象の音源からの音波に基づいて得られる値と理論値に基づく値との関係から、複数の仮想音源のうち、推定対象の音源と一致又は最も近い仮想音源を判断し、その仮想音源の配置位置を推定対象の音源の位置と推定するようにしたので、現実的でない仮定にとらわれず、また、音源深度が既知であるという条件を満たさない環境下でも、より現実の環境条件に近く、精度の高い音源位置の推定を行うことができる。
【００１０】
また、本発明によれば、到来時間差理論値算出手段が、音源の位置を推定する海域において配置した各仮想音源に対して差の理論値を算出し、一方で、到来時間差観測値算出手段が差の観測値を算出し、差の理論値及び差の観測値に基づいてコスト関数算出手段が一致度をコストとして算出し、音源位置推定手段がコストに基づいて、一致度が最も高い仮想音源の配置位置を推定対象の音源の位置と推定するようにしたので、より現実に沿った精度の高い音源位置推定を行うことができる装置を得ることができる。
【発明を実施するための最良の形態】
【００１１】
実施の形態１．
図１は本発明の第１の実施の形態に係る音源位置推定装置のブロック図である。本実施の形態及び以下の実施の形態において、音源位置推定装置は音響センサ（図示せず）等と共に、ソーナー（例えば逆探ソーナー。アクティブソーナー又はパッシブソーナーは問わない）に備えられている。図１において、入力端子１には音響センサからディジタルデータを含む信号が入力される。マルチパス波到来時間差測定器２は、ディジタルデータに基づいてマルチパス波の到来時間差観測値をデータとして算出する。到来時間差Ｔａｂｌｅ作成器３は設定した各仮想音源に対する到来時間差Ｔａｂｌｅ値（理論値）を算出し、到来時間差Ｔａｂｌｅ（テーブル形式のデータ）として記憶する。
【００１２】
また、コスト関数算出器４は、到来時間差観測値と仮想音源に基づく到来時間差Ｔａｂｌｅ値とに基づいて、それらの一致度をコストとして表すためのコスト関数の値をデータとして算出する。音源位置推定器５はコスト関数の値に基づいて音源位置を推定する。そして、出力端子６は推定した音源位置のデータを含む信号を出力するための端子である。
【００１３】
図２はコスト関数算出器４の詳細構成ブロック図である。図２において、入力端子９には到来時間差観測値のデータを含む信号が入力される。一方、入力端子１０には到来時間差Ｔａｂｌｅ値のデータを含む信号が入力される。Ｔａｂｌｅ値−観測値差算出器１１は到来時間差観測値と到来時間差Ｔａｂｌｅ値との差を算出する。自乗和算出器１２は到来時間差観測値と到来時間差Ｔａｂｌｅ値との差に基づいてコスト関数の値をデータとして算出する。そして、接続端子１３はコスト関数の値のデータを含む信号を音源位置推定器５に送信するための端子である。
【００１４】
ここで、ディジタルデータの処理を行うマルチパス波到来時間差測定器２、到来時間差Ｔａｂｌｅ作成器３、コスト関数算出器４（Ｔａｂｌｅ値−観測値差算出器１１、自乗和算出器１２）及び音源位置推定器５は、それぞれコンピュータ等のデータ処理手段（装置）で構成されているものとする。それぞれの手段（又は装置内のある記憶手段）には、後述する処理手順のプログラムが記憶されており、各手段はそのプログラムに基づいて処理を実行する。ここでは、複数の手段がそれぞれの処理を行っているが、また、１つのデータ処理手段が各手段のプログラムを実行して処理を行うようにしてもよい。これは後述する実施の形態における各算出器等においても同様である。
【００１５】
本実施の形態及び以下の実施の形態は、マルチパスにより各経路を経由して到来したパルスの時刻を理論値として算出した値（又はその値に基づく値）を有する複数の仮想音源を音源の位置を推定する海域に配置し、その値と実際に観測された値とに基づいて、推定対象の音源と一致又は最も近い仮想音源を判断することで、その仮想音源の配置位置を推定対象の音源と推定するようにしたものである。
【００１６】
次に本発明の実施の形態１に係る音源位置推定装置の動作について図１、２のブロック図を用い、図３〜７に沿って説明する。まず音源位置推定装置が位置推定処理を実際に行う前の到来時間差Ｔａｂｌｅの作成手順について説明する。到来時間差Ｔａｂｌｅの作成に際し、あらかじめ電気伝導度水温水深計ＣＴＤ（Conductivity-Temperature-Depth profiler ）やＸＢＴ（Expendable-Bathy-Thermometer）等により、音源位置を推定対象の海域における水温及び塩分濃度を測定し、その測定結果に基づいて音速プロファイル（深度と音速との関係）を算出し、これをデータ化しておく。音速プロファイルの作成については、既存の方法を用いることができる。ここでは、到来時間差Ｔａｂｌｅ作成器３は外部の装置によって作成された音速プロファイルのデータを用いるが、到来時間差Ｔａｂｌｅ作成器３が測定結果に基づいて音速プロファイルのデータを作成するようにしてもよい。
【００１７】
図３は仮想音源の設定を表す図である。ここでは、図３に示すように、推定対象の海域を仮想的にメッシュ（格子）状に区切り、各交点には仮想音源Ｓij（ｉ＝１，２，…，ｐ、ｊ＝１，２，…，ｑ）が存在するものと想定する。本実施の形態では、メッシュ幅（仮想音源の間隔）は特に限定しない。所望する音源位置の推定精度に応じて幅を設定し、仮想的に仮想音源を配置することができる。また、例えば、幅が狭くなるほど、データ量が多くなり、演算量も多くなるので、例えばリアルタイムでの処理を行う必要があれば、幅を広くして演算量を減らす等、状況に応じて設定を行うことができる。ここで、本実施の形態では各仮想音源を等間隔となるように想定しているが、各仮想音源の間隔を一律に想定する必要はない。さらに、本実施の形態では、すべての交点について仮想音源を想定するが、例えば任意の位置（交点）について想定するようにしてもよい。
【００１８】
図４は仮想音源の伝搬経路を示す図である。到来時間差Ｔａｂｌｅ作成器３では、まず、音速プロファイルのデータを用いて、図４に示すように、音線理論に基づいて、仮想音源Ｓijから発せられ、音響センサが受信するまでのマルチパス（音波の伝搬経路）を判断する。これを、想定したすべての仮想音源Ｓijについて判断する。ここで、音線とは音のエネルギ流れの経路（音の空間において伝搬する様子）を示す線である。
【００１９】
図５は到来時間差Ｔａｂｌｅを視覚的に表す図である。次に、判断した各仮想音源Ｓijのマルチパスにおいて、音響センサに到来する時刻が早い順に、第１パルス、第２パルス、…、第Ｎパルスとし、それらのパルスに基づいて、仮想音源Ｓijについて、伝搬時間理論値ｔ01（ｉ，ｊ），ｔ02（ｉ，ｊ），…，ｔ0N（ｉ，ｊ）を算出する。さらに、第１パルスと第２パルスの到来時間差τ120 （ｉ，ｊ）、第１パルスと第３パルスの到来時間差τ130 （ｉ，ｊ）、…、第１パルスと第Ｎパルスの到来時間差τ1N0 （ｉ，ｊ）を、次式（１）に基づいて算出する。これを想定したすべての仮想音源Ｓijについて行い、仮想音源Ｓijについて到来時間差Ｔａｂｌｅを作成する。作成した到来時間差Ｔａｂｌｅは、到来時間差Ｔａｂｌｅ作成器３が有する記憶手段（図示せず）にテーブル形式のデータとして記憶される。
【００２０】
【数１】
2006194627.tif 000003
【００２１】
次に、到来（観測）したパルスに基づいて音源位置推定装置が音源の位置を推定する動作について説明する。音源から到来したパルスを音響センサがディジタルデータ化し、そのデータを含む信号が入力端子１を介してマルチパス波到来時間差測定器２に入力される。マルチパス波到来時間差測定器２では、到来したパルスが既知の場合、音響センサからの信号が入力端子１を介して入力されると、既知のパルスのデータＳ１（ｋ）（ｋは時刻インデクスを表し、ｋ＝１，２，…，Ｋとする）と音響センサからのデータＳ２（ａ）（ａは時刻インデクスを表し、ａ＝１，２，…，Ａとする。Ｓ２（ａ）は複数のパルスを含むデータである）とを用い、次式（２）に基づいて相関結果Ｃ（ｄｔ）を算出する。さらに、次式（３）に示すように、相関結果Ｃ（ｄｔ）の自乗値をＣＰ（ｄｔ）として算出する。
【００２２】
【数２】
2006194627.tif 000004
【００２３】
そして、ＣＰ（ｄｔ）の値（ピーク点）を結んでいき、包絡線を求める。包絡線の極大となる箇所を判断し、判断した箇所の時刻ｔ1 ，ｔ2 ，…，ｔn が、各パルスがそれぞれの経路を伝搬して到来した時刻とする。極大となる箇所を判断する方法については、例えば、包絡線を求めた際に算出した各ピーク点の間の傾きに基づいて、その傾きが正から負に変わるピーク点を極大となる箇所とする。
【００２４】
一方、音源からのパルスが未知で、上記のような相関処理を行おうとする場合、入力端子１を介して音響センサから送信された信号に含まれるデータに基づいて、周波数、パルス長等を判断し、音源からのパルスと同様のパルスをレプリカ信号Ｓ１（ｋ）（ｋ＝１，２，…，Ｋ）として作成する。そして、上述した（２）式及び（３）式に基づいてＣＰ（ｄｔ）を算出し、包絡線を求めて、各パルスの到来時刻を決定する。
【００２５】
また、音源からのパルスが未知で、上記のような相関処理を行わずに各パルスの到来時刻を決定する場合、音響センサからのデータＳ２（ａ）の自乗値の包絡線を算出する。そして、該包絡線の立ち上がりにおいて、あらかじめ定めた閾値を越える箇所の時刻ｔ1 ，ｔ2 ，…，ｔn が、各パルスがそれぞれの経路を伝搬して到来した時刻とする。また、包絡線における立ち上がり時刻又は立ち下がり時刻を判断し、その中心の時刻をパルスが到来した時刻としてもよい。
【００２６】
図６は到来時間差観測値を表す図である。本実施の形態では、決定した各時刻に基づいて、マルチパス波到来時間差算出器２が、第１パルスと第２パルスの到来時間差τ12、第１パルスと第３パルスの到来時間差τ13、…、第１パルスと第ｎパルスの到来時間差τ1nを次式（４）に基づいて、到来時間差観測値のデータとして算出する（ｂ＝２，３，…，ｎ）。そして、到来時間差観測値のデータ（ｎ−１データ分。これらをベクトルτob＝［τ12，τ13，…，τ1n］とした一群のデータとして扱ってもよい）を含む信号をコスト関数算出器４に送信する。
τ1b＝ｔb −ｔ1                                                     …（４）
【００２７】
図７は差を算出する到来時間差観測値と到来時間差Ｔａｂｌｅ値との組み合わせを表す図である。コスト関数算出器４のＴａｂｌｅ値−観測値差算出器１１には、到来時間差観測値のデータτ12，τ13，…，τ1nを含む信号が接続端子９を介して入力される。また、到来時間差Ｔａｂｌｅ作成器３に記憶された、仮想音源Ｓijにおけるパルスの到来時間差Ｔａｂｌｅ値τ120 （ｉ，ｊ），τ130 （ｉ，ｊ），…，τ1N0 （ｉ，ｊ）のデータを含む信号が接続端子１０を介して入力される。
【００２８】
そして、次式（５）に基づいて、仮想音源Ｓijにおける到来時間差Ｔａｂｌｅ値と到来時間差観測値との差Δτ12（ｉ，ｊ），Δτ13（ｉ，ｊ），…，Δτ1n（ｉ，ｊ）を算出する（ここで、ｈ＝２，…，ｎである）。これを各仮想音源Ｓij（ｉ＝１，２，…，ｐ、ｊ＝１，２，…，ｑ）について行う。本実施の形態では、図７に示すように、τ12とτ120 （ｉ，ｊ）、τ13とτ130 （ｉ，ｊ）、…τ1nとτ1n0 （ｉ，ｊ）の組み合わせで演算を行うものとする（本実施の形態ではｎ＋１以上の到来時間差Ｔａｂｌｅ値は用いないものとする）。
Δτ1h（ｉ，ｊ）＝τ1h0 （ｉ，ｊ）−τ1h                            …（５）
【００２９】
自乗和算出器１２は、Ｔａｂｌｅ値−観測値差算出器１１が算出した仮想音源Ｓijにおける到来時間差Ｔａｂｌｅ値と到来時間差観測値との差のデータから、次式（６）に基づいてコスト関数Ｊijを算出する。そして、各仮想音源Ｓijについて算出したコスト関数Ｊijのデータを含む信号を、接続端子１３を介して音源位置推定器５に出力する。
【００３０】
【数３】
2006194627.tif 000005
【００３１】
音源位置推定器５は、各仮想音源Ｓijのコスト関数Ｊijに基づいて、その中で最小値のコスト関数Ｊijを決定する。そして、最小のコスト関数Ｊijが算出された仮想音源Ｓijについて設定した位置を、音源の位置として推定する。ここで、各コスト関数Ｊijは、到来時間差観測値と各仮想音源の到来時間差Ｔａｂｌｅ値との一致度を表す。一致度合いが高いほど、到来時間差観測値と到来時間差Ｔａｂｌｅ値との差が小さくなり、コスト関数Ｊijも小さくなる。したがってコスト関数Ｊijが最小の仮想音源Ｓijが、観測したパルスを発する音源に最も近いものと推定できる。
【００３２】
以上のように第１の実施の形態によれば、海中の音速プロファイルを考慮し、音波の屈折現象を反映した音線理論に基づいて仮想的に設定した仮想音源Ｓijの到来時間差Ｔａｂｌｅ値を到来時間差Ｔａｂｌｅ作成器３が算出して記憶しておき、その到来時間差Ｔａｂｌｅ値とマルチパス波到来時間差測定器２が算出した実際の到来時間差観測値との差を算出し、コスト関数算出器４が第１パルスと第２パルス以降の各パルスとの時間差の一致度となるコスト関数Ｊijを算出し、さらに、音源位置推定器５が最も一致度が高い仮想音源に設定された位置を、音源位置として判断するようにしたので、上述した仮定を満足しない条件下においても、その仮定にとらわれることなく、精度良く音源位置を推定することができる。また、音源深度が既知である等の条件を設定しなくてもよいので、音源深度が未知（条件を充たさない）である場合でも音源位置を推定することができる。
【００３３】
実施の形態２．
図８は差を算出する到来時間差Ｔａｂｌｅ値と到来時間差観測値との関係を表す図である。図８（ａ）は、到来時間差観測値τ12，τ13，…，τ1nと到来時間差Ｔａｂｌｅ値τ120 （ｉ，ｊ），τ130 （ｉ，ｊ），…，τ1N0 （ｉ，ｊ）を示す。また、図８（ｂ）は、到来時間差観測値τ12，τ13，…，τ1n及び到来時間差Ｔａｂｌｅ値τ120 （ｉ，ｊ），τ130 （ｉ，ｊ），…，τ1N0 （ｉ，ｊ）の組み合わせを時間軸を用いて表したものである。第２の実施の形態における装置は、Ｔａｂｌｅ値−観測値差算出器１１において、到来時間差Ｔａｂｌｅ値と到来時間差観測値との組み合わせを決定するための処理が異なるものの、その装置構成は第１の実施の形態で説明したものと同じである。そこで、本実施の形態では、装置構成については図１及び図２を用い、処理については図８及び９を参照しながら、本実施の形態における音源位置推定装置について説明することにする。
【００３４】
例えば受信信号のＳ／Ｎが低い場合、（ア）、（イ）のような誤検出が発生することがあり、到来時間差観測値と到来時間差Ｔａｂｌｅ値を単純に時間差に沿って選択すると、差を算出するための組み合わせが食い違ってしまうことがある。
（ア）雑音を誤ってパルスとして検出する
（イ）存在するはずのパルスを検出できない
差をとる組み合わせに食い違いが生じた結果、例えば真の音源位置に最も近い仮想音源Ｓijの到来時間差Ｔａｂｌｅ値と到来時間差観測値の間でも、伝搬経路（海面及び海底の反射回数）が異なるものどうしの差を算出することになるため、コスト関数が最小とならない可能性がある。本実施の形態は、Ｔａｂｌｅ値−観測値差算出器１１において、到来時間差観測値と到来時間差Ｔａｂｌｅ値との間で、図８（ｂ）の矢印によって示すように、互いに最も近い値どうしを１対１で対応させ、差を算出するための組み合わせを決定するようにしたものである。そして、組み合わせた到来時間差観測値と到来時間差Ｔａｂｌｅ値の差Δτ（ｘ）（ｘ＝１，２，…，Ｘ）を算出する。Ｘは到来時間差観測値と到来時間差Ｔａｂｌｅ値の組み合わせの数を表す。
【００３５】
図９は、Ｔａｂｌｅ値−観測値差算出器１１による、組み合わせ決定を含めた処理の一例を表すフローチャートである。ここでは、ある仮想音源Ｓijの到来時間差Ｔａｂｌｅ値と到来時間差観測値との組み合わせを決定する手順についてさらに詳しく説明する。ここで図９においては仮想音源の位置に関する添字（ｉ，ｊ）を省略している。
【００３６】
まず、コスト関数算出器４のＴａｂｌｅ値−観測値差算出器１１には、到来時間差観測値のデータτ12，τ13，…，τ1nを含む信号が接続端子９を介して入力される（ここでは、これらのデータの一群をベクトルτob＝［τ12，τ13，…，τ1n］として表す）。また、到来時間差Ｔａｂｌｅ作成器３に記憶された、仮想音源Ｓijにおけるパルスの到来時間差Ｔａｂｌｅ値τ120 （ｉ，ｊ），τ130 （ｉ，ｊ），…，τ1N0 （ｉ，ｊ）のデータを含む信号が接続端子１０を介して入力される（ここでは、これらのデータの一群をベクトルτTable ＝［τ120 ，τ130 ，…，τ1N0 ］として表す）。
【００３７】
Ｔａｂｌｅ値−観測値差算出器１１では、まず、τob（１）（＝τ12）について、最も近い到来時間差Ｔａｂｌｅ値をτTable の要素の中から決定する。つまり、τob（１）との差の絶対値が最も小さくなるτTable の要素を決定する。そして、決定した要素をτnearとし、決定したτTable の要素のインデクス（例えば要素を表す名前、番号等）をインデクスＩ（１）にデータとして代入する。さらに、τob（１）と最も近い到来時間差Ｔａｂｌｅ値との差の絶対値を算出し、Δτ（１）に代入する。以上を式で表すと次式（７）のようになる（ブロックａ）。
｜τnear−τob（１）｜＝Δτ（１）                                  …（７）
【００３８】
次に、インデクスを指定するためのポインタｘに２を代入する（ブロックｂ）。そして、τob（ｘ）（ここではｘ＝２）について、前述したように差の絶対値に基づいて、最も近い到来時間差Ｔａｂｌｅ値をτTable の要素の中から決定する。そして、決定した要素をτnearとし、決定したτTable の要素のインデクスをＩ（ｘ）にデータとして代入する。さらに、到来時間差Ｔａｂｌｅ値と到来時間差観測値の差の絶対値を算出し、Δτ（ｘ）に代入する（ブロックｃ）。以上を式で表すと次式（８）のようになる。
｜τnear−τob（ｘ）｜＝Δτ（ｘ）                                  …（８）
【００３９】
次に、インデクスＩ（ｘ−１）のデータとＩ（ｘ）のデータとの関係について、Ｉ（ｘ−１）＝Ｉ（ｘ）であるかどうか、つまり、τob（ｘ−１）とτob（ｘ）とについて、最も近いτTable の要素が共通しているかどうかを判断する（分岐ｄ）。
【００４０】
Ｉ（ｘ−１）＝Ｉ（ｘ）であると判断すると、さらにΔτ（ｘ−１）≧Δτ（ｘ）であるかどうか、つまりτob（ｘ−１）とτob（ｘ）とのどちらがよりτnearに近いかどうかを判断する（分岐ｅ）。Δτ（ｘ−１）≧Δτ（ｘ）であれば（τob（ｘ−１）とτob（ｘ）とが同じ又はτob（ｘ）の方がτnearに近ければ）、τob（ｘ）とτnearとを差をとる組み合わせとして採用することとし、Δτ（ｘ−１）に０を代入する（ブロックｆ）。一方、Δτ（ｘ−１）≧Δτ（ｘ）でなければ（τob（ｘ−１）の方がτnearに近ければ）、τob（ｘ−１）とτnearとを差をとる組み合わせとして採用することとし、Δτ（ｘ）に０を代入する（ブロックｇ）。そして、ｘにｘ＋１を代入してｘを１増加させる（ブロックｈ）。
【００４１】
ｘ＞ｎ−１であるかどうか（つまり、ベクトルτobの最後（ｎ−１番目）の要素（到来時間差観測値τ1n）について処理を行ったかどうか、すべての到来時間差観測値について処理を行ったかどうか）を判断する（分岐ｉ）。すべての要素について処理が終了していないと判断すると、次の要素について処理を続ける。すべての要素について処理が終了していると判断すると、ある仮想音源Ｓijにおける処理を終了する。
【００４２】
以上のように決定した組み合わせによりＴａｂｌｅ値−観測値差算出器１１が算出したΔτ（ｘ）（ｘ＝１，２，…，ｎ−１）に基づいて、第１の実施の形態と同様に自乗和算出器１２がコスト関数Ｊijを算出する。以上の処理を、コスト関数算出器４（Ｔａｂｌｅ値−観測値差算出器１１、自乗和算出器１２）は、各仮想音源Ｓijについて行い、算出したコスト関数Ｊijのデータを含む信号を、接続端子１３を介して音源位置推定器５に出力する。音源位置推定器５での音源位置の推定処理は、第１の実施の形態で説明した処理と同様であるので説明を省略する。
【００４３】
以上のように第２の実施の形態によれば、Ｔａｂｌｅ値−観測値差算出器１１が、互いに最も近い到来時間差Ｔａｂｌｅ値と到来時間差観測値との組み合わせを１対１で決定し、その組み合わせに基づいてコスト関数Ｊij算出のための差を算出するようにしたので、真の音源の位置の近傍において、海面及び海底における反射回数による経路が同じもの同士の到来時間差Ｔａｂｌｅ値と到来時間差観測値との組み合わせの差がコストとして計上されることとなり、誤検出に基づく余分な到来時間差観測値や到来時間差Ｔａｂｌｅ値はコスト関数Ｊijの算出に影響しないため、誤検出に対して、ロバストな音源位置を推定することができる。なお、本実施の形態において説明した図９のフローチャートの処理手順は一例であり、本実施の形態のように、互いに最も近い到来時間差観測値と到来時間差Ｔａｂｌｅ値との組み合わせを決定してΔτ（ｘ）が算出できるのであれば処理手順は限定しない。
【００４４】
実施の形態３．
図１０は差を算出するＴａｂｌｅ値と到来時間差観測値との関係を表す図である。図１０（ａ）は、到来時間差観測値τ12，τ13，…，τ1nと到来時間差Ｔａｂｌｅ値τ120 （ｉ，ｊ），τ130 （ｉ，ｊ），…，τ1N0 （ｉ，ｊ）を示す。また、図１０（ｂ）は、到来時間差観測値τ12，τ13，…，τ1n及び到来時間差Ｔａｂｌｅ値τ120 （ｉ，ｊ），τ130 （ｉ，ｊ），…，τ1N0 （ｉ，ｊ）の組み合わせを時間軸を用いて表したものである。第３の実施の形態における装置についても、Ｔａｂｌｅ値−観測値差算出器１１において、到来時間差Ｔａｂｌｅ値と到来時間差観測値との組み合わせを決定する処理が異なるものの、その構成は第１の実施の形態で説明したものと同じであるので、装置については図１及び図２を用いて説明する。
【００４５】
例えば、マルチパスの多い環境において、Ｓ／Ｎが低く、マルチパス波の到来時刻推定精度が高くない場合、真の音源位置に近い仮想音源とそうでない仮想音源とのコスト関数の値において差が強調されず、グローバルミニマムとなる可能性がある。そこで、本実施の形態は、図１０（ｂ）に示す矢印Ａのように、到来時間差観測値に対し、到来時間差Ｔａｂｌｅ値の近い値を組み合わせていき、組み合わせができなかった到来時間差Ｔａｂｌｅ値についても、矢印Ｂのように、さらに最も近い到来時間差観測値と組み合わせることで、到来時間差Ｔａｂｌｅ値を有効に利用しようとするものである（図１０（ｂ）ではτ13とτ130 （ｉ，ｊ），τ140 （ｉ，ｊ））。したがって、本実施の形態は必ずしも１対１の対応になるとは限らない。ここで、第１及び第２の実施の形態と同様に、到来時間差観測値のうち、最も時間差が大きいτ1nと最も近い到来時間差Ｔａｂｌｅ値（図１０（ｂ）ではτ1(n-1)0 （ｉ，ｊ））よりも大きな時間差を有する到来時間差Ｔａｂｌｅ値は採用しないものとする。
【００４６】
図１１はＴａｂｌｅ値−観測値差算出器１１による、組み合わせ決定を含めた処理の一例を表すフローチャートである。ここでは、ある仮想音源Ｓijの到来時間差Ｔａｂｌｅ値と到来時間差観測値との組み合わせを決定する手順について説明する（図１１では仮想音源の位置に関する添字（ｉ，ｊ）を省略している）。
【００４７】
コスト関数算出器４のＴａｂｌｅ値−観測値差算出器１１には、到来時間差観測値のデータτ12，τ13，…，τ1nを含む信号が接続端子９を介して入力される（ここでは、これらのデータの一群をベクトルτob＝［τ12，τ13，…，τ1n］として表す場合もある）。また、到来時間差Ｔａｂｌｅ作成器３に記憶された、仮想音源Ｓijにおけるパルスの到来時間差Ｔａｂｌｅ値τ120 （ｉ，ｊ），τ130 （ｉ，ｊ），…，τ1N0 （ｉ，ｊ）のデータを含む信号が接続端子１０を介して入力される。
【００４８】
Ｔａｂｌｅ値−観測値差算出器１１は、まず、到来時間差観測値のインデクスを指定するポインタｘに１を代入する（ブロックｊ）。そして、τob（ｘ）（ここではｘ＝１）について、最も近い到来時間差Ｔａｂｌｅ値をτTable の要素の中から決定する。そして、決定した要素をτnearとし、決定したτTable のインデクスをＩ（ｘ）に代入する。さらに、到来時間差Ｔａｂｌｅ値と到来時間差観測値の差を表すベクトルΔτ（ｘ）を前述した（８）式に基づいて算出する（ブロックｋ）。
【００４９】
次にｘにｘ＋１を代入してｘを１増加させる（ブロックｌ）。そして、ベクトルτobを構成するｎ−１の要素（到来時間差観測値τ12，τ13，…，τ1n）について処理を行ったかどうかを判断する（分岐ｍ）。すべての要素について処理が終了していないと判断すると、次の要素について処理を続ける。
【００５０】
分岐ｍにおいて、すべての要素について処理が終了していると判断すると、τTable のインデクスの中での最大値Ｉ（ｎ−１）の値をＩmax とする（ブロックｎ）。図１０（ｂ）では、仮想音源Ｓijの到来時間差Ｔａｂｌｅ値におけるτ1(n-1)0 （ｉ，ｊ）のインデクスとなる。そして、τTable のインデクスを指定するポインタｙに１を代入する（ブロックｏ）。さらにそのポインタｙの値がベクトルＩの要素となっているかどうかを判断する（分岐ｐ）。
【００５１】
ポインタｙの値がベクトルＩの要素となっていない、つまり、τTable （ｙ）がベクトルτobとの間で、差を算出する組み合わせに採用されていないと判断すると、最も近い到来時間差観測値をベクトルτobの要素の中から決定し、その要素をτnearとする。そして、到来時間差Ｔａｂｌｅ値と到来時間差観測値の差を表すベクトルΔτ（ｘ）を算出する（ブロックｑ）。以上を式で表すと、次式（９）のようになる。その後、ｘにｘ＋１を代入してｘを１増加させる（ブロックｒ）。
｜τnear−τTable（ｙ）｜＝Δτ（ｘ）                               …（９）
【００５２】
そして、ｙにｙ＋１を代入してｙを１増加させる（ブロックｓ）。そして、分岐ｐ〜ｓまでの処理をｙがＩmax よりも大きくなったものと判断するまで繰り返す（分岐ｔ）。ｙがＩmax よりも大きくなったものと判断すると、ある仮想音源Ｓijにおいて、到来時間差観測値と到来時間差Ｔａｂｌｅ値との組み合わせを決定する処理を終了する。
【００５３】
以上のように決定した組み合わせによりＴａｂｌｅ値−観測値差算出器１１が算出したΔτ（ｘ）（ｘ＝１，２，…，Ｘ、Ｘは組み合わせの数）に基づいて、第１の実施の形態と同様に自乗和算出器１２がコスト関数Ｊijを算出する。以上の処理をコスト関数算出器４（Ｔａｂｌｅ値−観測値差算出器１１、自乗和算出器１２）は仮想音源Ｓij毎に行い、各仮想音源Ｓijについて算出したコスト関数Ｊijのデータを含む信号を、接続端子１３を介して音源位置推定器５に出力する。音源位置推定器５での音源位置の推定処理は、第１の実施の形態で説明した処理と同様であるので説明を省略する。
【００５４】
以上のように第３の実施の形態によれば、実際に到来したパルスに基づく到来時間差観測値及び到来時間差観測値の最大値に最も近い到来時間差Ｔａｂｌｅ値以下の値を有する到来時間差Ｔａｂｌｅ値を余すことなく使用し、コスト関数Ｊijの値を算出し、音源位置の推定を行うようにしたので、一致度を表すコスト関数Ｊijの値について、実際の音源位置と近い仮想音源とそうでない仮想音源との差を強調し、真の音源位置近傍の仮想音源Ｓijにおいて算出したコスト関数Ｊijの勾配を急峻にすることで、精度のよい推定を行うことができる。さらに、真の音源位置から離れた仮想音源においては、互いにずれが大きく対応しない到来時間差観測値、到来時間差Ｔａｂｌｅ値の差に基づく値がコスト関数Ｊijの値に計上することになり、ローカルミニマムの値を大きくすることができる。これにより、マルチパスが多く、到来時刻推定精度が高くない条件の下でもロバストな音源位置の推定を行うことができる。なお、本実施の形態において説明した図１１のフローチャートの処理手順は一例であり、本実施の形態のように、到来時間差観測値及び到来時間差Ｔａｂｌｅ値を有効に利用してΔτ（ｘ）が算出できるのであれば処理手順は限定しない。
【００５５】
実施の形態４．
図１２は本発明の第４の実施の形態に係る音源位置推定装置のブロック図である。図１２において図１と同じ符号を付しているものは、上述した実施の形態と同様の動作を行うので説明を省略する。レプリカ相関器１６は、観測した音波に基づく値を算出するデータ算出手段であり、レプリカ信号に含まれる既知のパルスのデータと、入力端子１を介して音響センサから送信された信号に含まれるデータとの間で、相関結果であるＣ（ｄｔ）を算出する。図１２では、装置外部（例えば記憶手段等）から送信されたレプリカ信号に基づいて相関結果を算出しているが、例えば、上述したように、パルスが未知の場合、音響センサシステムからの受波したデータに基づいてレプリカ信号とするためのデータを作成することもある。また、一度送信されたレプリカ信号に含まれるデータを記憶手段（図示せず）に記憶させておくこともできる。包絡線算出器１７は、観測した音波に基づく値を処理する処理手段であり、相関結果の自乗値ＣＰ（ｄｔ）を算出し、そのピーク点を結ぶことにより、ＣＰ（ｄｔ）の包絡線ＣＰＥを求める。第１パルス到来時刻推定器１８は、包絡線の中で最初に閾値を越え、かつ極大となる時刻を求める。これを第１パルス到来時刻推定値ｔ１とする。また、コスト関数算出器１９は、コスト関数Ｊijの値をデータとして算出するが、その算出手順が、後述するようにコスト関数算出器４とは異なる。
【００５６】
図１３は包絡線ＣＰＥの一例を表す図である。求めた包絡線においてはそれぞれ１の到来時間差観測値が算出される。ここで複数のパルスが重なっているため、通常分離するはずの包絡線が分離できていない場合、上述の第１〜３の実施の形態では、１の包絡線として扱われ、１つの到来時間差観測値しか算出されないことになる。したがって、真の音源位置又はその近傍の仮想音源の到来時間差Ｔａｂｌｅ値と到来時間差観測値との間でも対応関係が異なってしまい、算出したコスト関数Ｊijの値が大きくなり、最小値にならない場合がある。このような場合、本来とは異なる仮想音源の位置が音源位置として推定されてしまうため、正確な音源位置が推定できていないことになってしまう。そこで、本実施の形態は、各到来時間差Ｔａｂｌｅ値に対応する相関結果の自乗値の和をコスト関数Ｊijの値として算出することにより、分離できていない包絡線においても複数の到来時間差Ｔａｂｌｅ値との対応を考慮したコスト関数を算出することができる。
【００５７】
次に本実施の形態に係る音源位置推定装置の動作について、図１２の構成に基づいて、図１４、１５を参照しながら説明する。レプリカ相関器１６は、入力端子１を介して音響センサから送信された信号に含まれるデータ及び既知のパルスのデータから、第１の実施の形態で説明した（２）式に基づいてレプリカ相関処理を行い、相関結果Ｃ（ｄｔ）のデータを含む信号を包絡線算出器１７に出力する。包絡線算出器１７は、（３）式に基づいて相関結果Ｃ（ｄｔ）の自乗値をＣＰ（ｄｔ）として算出する。そして、ＣＰ（ｄｔ）におけるピーク点を結ぶことにより、包絡線ＣＰＥを求め、そのデータを含む信号を第１パルス到来時刻推定器１８及びコスト関数算出器１９に出力する。ここで、相関出力の自乗値ＣＰ（ｄｔ）の全てのピーク点を、さらに２次関数フィッティング等を行って結び、包絡線ＣＰＥを求めるようにしてもよい。この場合、包絡線ＣＰＥを精度よく求められるため、音源位置推定精度を向上させることができる。
【００５８】
図１４は包絡線ＣＰＥと第１パルス到来時刻推定値ｔ１との関係を表す図である。第１パルス到来時刻推定器１８は、包絡線ＣＰＥのデータに基づいて、包絡線ＣＰＥの中であらかじめ設定した閾値を越えた部分があるかどうかを判断し、越えた部分に基づいて最初の極大点を検出する。そして、検出した極大点における時刻を第１パルスの到来時刻ｔ１として、そのデータを含む信号をコスト関数算出器１９に出力する。
【００５９】
図１５は包絡線ＣＰＥと到来時間差Ｔａｂｌｅ値との関係を表す図である。コスト関数算出器１９には、包絡線ＣＰＥのデータ及び第１パルスの到来時刻推定値ｔ１を含む信号に加え、到来時間差Ｔａｂｌｅ作成器３に記憶された、仮想音源Ｓijにおけるパルスの到来時間差Ｔａｂｌｅ値τ120 （ｉ，ｊ），τ130 （ｉ，ｊ），…，τ1N0 （ｉ，ｊ）のデータを含む信号が入力される。そして、包絡線ＣＰＥのデータ、第１パルスの到来時刻推定値ｔ１及びパルスの到来時間差Ｔａｂｌｅ値のデータに基づいて、図１５に示すように、時刻ｔ１＋τ120 （ｉ，ｊ），ｔ１＋τ130 （ｉ，ｊ），…，ｔ１＋τ1N0 （ｉ，ｊ）における包絡線ＣＰＥの値を求める。そして、その総和値の逆数を、仮想音源Ｓijにおけるコスト関数Ｊijの値とし、そのデータを含む信号を音源位置推定器５に出力する。以上を式で表すと、次式（１０）となる。（１０）式における包絡線ＣＰＥの値ＣＰＥ（ｔ１＋τ1z0 （ｉ，ｊ））（ｚ＝１，２，…，Ｎ）は、線形補間等の補間処理によって求める。また、ｔ１＋τ1z0 （ｉ，ｊ）に最も近い時刻に算出された包絡線による値を用いるようにしてもよい。ここで、τ110 （ｉ，ｊ）＝０とする。音源位置推定器５での音源位置の推定処理は、第１の実施の形態で説明した処理と同様であるので説明を省略する。
【００６０】
【数４】
2006194627.tif 000006
【００６１】
以上のように第４の実施の形態によれば、コスト関数算出器１９が、パルスの到来時間差Ｔａｂｌｅ値に基づく仮想音源Ｓijにおけるコスト関数の値Ｊijを相関処理がなされた結果である包絡線ＣＰＥから直接的に導き出すようにしたので、パルスが重なって受信されたために、本来分離されているべき包絡線が分離されていない場合でも、それを１つの包絡線として扱ってしまうことがなくなるため、真の音源位置及びその近傍に対応する仮想音源Ｓijのコスト関数の値Ｊijが大きくならず、最小となる。そのため、重なって受信したパルスを分離できないことによる誤検出が生じた場合に対してロバストな音源位置を推定することができる。
【００６２】
実施の形態５．
図１６は本発明の第５の実施の形態に係る音源位置推定装置のブロック図である。図１６において図１及び図１２と同じ符号を付しているものは、上述した実施の形態と同様の動作を行うので説明を省略する。第１パルス到来時刻幅算出器２０は、第１パルス到来時刻に幅をもたせるための開始時刻ｔstart と終了時刻ｔend を決定する。
【００６３】
図１７はコスト関数算出器２１の詳細構成ブロック図である。図１７において、入力端子２２には包絡線算出器１７から包絡線のデータを含む信号が入力される。また、入力端子２３には第１パルス到来時刻幅算出器２０から開始時刻ｔstart と終了時刻ｔend のデータを含む信号が入力される。さらに、入力端子２５には到来時間差Ｔａｂｌｅ値のデータを含む信号が入力される。コストベクトル算出器２４は、第１パルス到来時刻の幅内でのそれぞれのコスト値を算出し、コストベクトルとする。コスト算出器２６は、コストベクトルに基づいて仮想音源Ｓijにおけるコスト関数の値Ｊijを算出する。接続端子２７はコスト関数の値Ｊijのデータを含む信号を音源位置推定器５に送信するための端子である。
【００６４】
例えば、Ｓ／Ｎ比が低く、第１パルスの推定精度が高くない場合、基準となる第１パルスの到来時刻推定値を固定してコスト関数の値Ｊijを算出すると、第１パルスの到来時刻にずれが生じた場合には、求める包絡線ＣＰＥの値ＣＰＥ（ｔ１＋τ1z0 （ｉ，ｊ））もずれるため、コスト関数の値Ｊijも異なり、音源位置推定精度も劣化する可能性がある。そこで、本実施の形態では、第１パルスの到来時刻に幅を持たせて、その幅で第１パルス到来時刻推定値を変化させながら、コスト関数の値Ｊijを算出し、音源位置推定の精度を高めようとするものである。
【００６５】
図１８は包絡線ＣＰＥと第１パルス到来時刻幅との関係を表す図である。次に、第５の実施の形態における音源位置推定装置の動作について、図１６、１７のブロック図に沿って、図１８、１９を参照しながら説明する。第１パルス到来時刻幅算出器２０は、包絡線ＣＰＥのデータに基づいて、包絡線ＣＰＥの中であらかじめ設定した閾値を越えた部分をもつ包絡線を第１パルスを含む包絡線であると判断し、越えた部分に基づいて、その閾値を越える開始時刻ｔstart 及び終了時刻ｔend を求める。求めた開始時刻ｔstart 及び終了時刻ｔend のデータを含む信号をコスト関数算出器２１に出力する。
【００６６】
図１９は包絡線ＣＰＥと到来時間差Ｔａｂｌｅ値との関係を表す図である。コスト関数算出器２１のコストベクトル算出器２４には、包絡線ＣＰＥのデータ並びに開始時刻ｔstart 及び終了時刻ｔend のデータを含む信号が、それぞれ入力端子２２、２３を介して入力される。また、到来時間差Ｔａｂｌｅ作成器３に記憶された、仮想音源Ｓijにおけるパルスの到来時間差Ｔａｂｌｅ値τ120 （ｉ，ｊ），τ130 （ｉ，ｊ），…，τ1N0 （ｉ，ｊ）のデータを含む信号が接続端子２５を介して入力される。
【００６７】
コストベクトル算出器２４は、まず、開始時刻ｔstart 及び終了時刻ｔend のデータに基づいて、次式（１１）に示すΔｔ（ｓ）を算出する。ここで、ｆｓはサンプリング周波数であり、ｓ＝１，２，…，Ｗ、Ｗ＝（ｔend −ｔstart ）ｆｓ＋１である。
Δｔ（ｓ）＝（ｓ−１）／ｆｓ                                      …（１１）
【００６８】
そして、包絡線ＣＰＥのデータ、第１パルスの到来時刻推定値ｔ１及びパルスの到来時間差Ｔａｂｌｅ値のデータに基づいて、図１９に示すように、時刻ｔstart ＋Δｔ（ｓ）＋τ120 （ｉ，ｊ），ｔstart ＋Δｔ（ｓ）＋τ130 （ｉ，ｊ），…，ｔstart ＋Δｔ（ｓ）＋τ1N0 （ｉ，ｊ）における包絡線ＣＰＥの値を求める（ｚ＝２，３，…，Ｎ）。あるΔｔ（ｓ）における包絡線ＣＰＥの値の総和に対し、その逆数の値をコストベクトルＪＶＥＣijの第ｓ要素とする。以上を式で表すと、次式（１２）となる。これにより、仮想音源Ｓijについてのコストベクトルは、ＪＶＥＣij＝［Ｊ１ij，Ｊ２ij，…，Ｊｓij，…ＪＷij］となる。
【００６９】
【数５】
2006194627.tif 000007
【００７０】
コスト算出器２６は、次式（１３）に示すように、コストベクトルＪＶＥＣijに基づいて、各要素の値の中で最小のものを仮想音源Ｓijにおけるコスト関数の値Ｊijとして決定する。そして、各仮想音源Ｓijにおけるコスト関数の値のデータを含む信号を、接続端子２７を介して音源位置推定器５に出力する。音源位置推定器５での音源位置の推定処理は、第１の実施の形態で説明した処理と同様であるので説明を省略する。
Ｊij＝MIN［ＪＶＥＣij］
ここで、ｉ＝１，２，…，ｐ、ｊ＝１，２，…，ｑであり、MIN［］ は、［］内の要素の最小値をとることを表す。
【００７１】
以上のように第５の実施の形態によれば、第１パルス到来時刻幅算出器２０が開始時刻ｔstart及び終了時刻ｔendを求め、第１パルスの到来時刻推定値にｔstart〜ｔendまでの幅を持たせ、コストベクトル算出器２４がその幅の間で第１パルス到来時刻推定値を変化させて、コスト算出器２６が仮想音源Ｓijのコスト関数の値Ｊijをデータとして算出し、各仮想音源Ｓijのコスト関数の値Ｊij中の最小値に基づいて音源位置推定器６が音源位置を推定するようにしたので、精度の高い第１パルスの到来時刻推定値に基づいて音源位置を推定することができる。そのため、第１パルスの到来時刻推定誤差に対してロバストな音源位置を推定することができる。
【００７２】
実施の形態６．
上述の第５の実施の形態では、コストベクトルＪＶＥＣijの要素の中の最小値を仮想音源Ｓijにおけるコスト関数の値Ｊijと決定したが、コストベクトルＪＶＥＣijの要素の総和の値、平均値等をコスト関数の値Ｊijとして決定するようにしてもよい。
【００７３】
また、上述の第４及び第５の実施の形態では、レプリカ信号と音響センサによる受信信号との相関処理結果の自乗値の包絡線ＣＰＥ及び到来時間差Ｔａｂｌｅ値に基づいてコスト関数の値Ｊijを算出するようにしたが、例えば、受波した音波自体の自乗値の包絡線を用いてコスト関数の値Ｊijを算出することもできる。本発明は、１つの音響センサが受波した音波による信号を用いて行ってもよいし、複数の音響センサで受波した音波を整相処理した信号を用いて行ってもよい。
【００７４】
さらに、第１〜第５の実施の形態で示した方法を、それぞれ並列に実行し、算出した音源位置の推定結果の中から尤もらしい結果に基づいて、最終的な音源位置を推定するようにすることもできる。
【００７５】
上述の実施の形態では、音源位置推定器５において、推定対象と仮想音源との一致度をコストとして算出し、最も一致度が高い仮想音源の位置を推定対象の位置とした。ただ、本発明はこれに限定されるものではない。例えば、最も一致度が高い仮想音源のコスト関数の値と、その近傍（例えば仮想音源がマトリクス状に配置されている場合、隣接する８つの仮想音源）のコスト関数との値に基づいて、仮想音源間の位置に推定対象の位置を求めることにより、さらに精度の高い位置推定を行うようにしてもよい。
【００７６】
実施の形態７．
上述の実施の形態では、音波が伝搬する媒体を海水としたが、本発明は、液体をはじめとする他の媒体を伝搬する音波を発する音源の位置推定についても適用することができる。
【図面の簡単な説明】
【００７７】
【図１】本発明の第１の実施の形態に係る音源位置推定装置のブロック図である。
【図２】コスト関数算出器４の詳細構成ブロック図である。
【図３】仮想音源の設定を表す図である。
【図４】仮想音源の伝搬経路を示す図である。
【図５】到来時間差Ｔａｂｌｅを視覚的に表す図である。
【図６】到来時間差観測値を表す図である。
【図７】差を算出する到来時間差観測値と到来時間差Ｔａｂｌｅ値との組み合わせを表す図である。
【図８】差を算出する到来時間差Ｔａｂｌｅ値と到来時間差観測値との関係を表す図である。
【図９】Ｔａｂｌｅ値−観測値差算出器１１による、組み合わせ決定を含めた処理を表すフローチャートである。
【図１０】差を算出するＴａｂｌｅ値と到来時間差観測値との関係を表す図である。
【図１１】Ｔａｂｌｅ値−観測値差算出器１１による、組み合わせ決定を含めた処理の一例を表すフローチャートである。
【図１２】本発明の第４の実施の形態に係る音源位置推定装置のブロック図である。
【図１３】包絡線ＣＰＥの一例を表す図である。
【図１４】包絡線ＣＰＥと第１パルス到来時刻推定値ｔ１との関係を表す図である。
【図１５】包絡線ＣＰＥと到来時間差Ｔａｂｌｅ値との関係を表す図である。
【図１６】本発明の第５の実施の形態に係る音源位置推定装置のブロック図である。
【図１７】コスト関数算出器２１の詳細構成ブロック図である。
【図１８】包絡線ＣＰＥと第１パルス到来時刻幅との関係を表す図である。
【図１９】包絡線ＣＰＥと到来時間差Ｔａｂｌｅ値との関係を表す図である。
【符号の説明】
【００７８】
１  入力端子
２  マルチパス波到来時間差測定器
３  到来時間差Ｔａｂｌｅ作成器
４  コスト関数算出器
５  音源位置推定器
６  出力端子
９、１０  入力端子
１１  Ｔａｂｌｅ値−観測値差算出器
１２  自乗和算出器
１３  接続端子
１６  レプリカ相関器
１７  包絡線算出器
１８  第１パルス到来時刻推定器
１９  コスト関数算出器
２０  第１パルス到来時刻幅算出器
２１  コスト関数算出器
２２、２３、２５  入力端子
２４  コストベクトル算出器
２６  コスト算出器
２７  接続端子
【図１】
2006194627.tif 000008
【図２】
2006194627.tif 000009
【図３】
2006194627.tif 000010
【図４】
2006194627.tif 000011
【図５】
2006194627.tif 000012
【図６】
2006194627.tif 000013
【図７】
2006194627.tif 000014
【図８】
2006194627.tif 000015
【図９】
2006194627.tif 000016
【図１０】
2006194627.tif 000017
【図１１】
2006194627.tif 000018
【図１２】
2006194627.tif 000019
【図１３】
2006194627.tif 000020
【図１４】
2006194627.tif 000021
【図１５】
2006194627.tif 000022
【図１６】
2006194627.tif 000023
【図１７】
2006194627.tif 000024
【図１８】
2006194627.tif 000025
【図１９】
2006194627.tif 000026