移動軌跡算出装置および移動軌跡算出方法
【課題】少ない計算量かつ高精度で、移動体の軌跡を算出する。
【解決手段】本発明の移動軌跡算出装置は、観測エリアを観測して得た、観測時刻と複数の観測値とを含む複数のセンサデータを記憶する手段と、前記観測エリアを分割した複数の区間における隣接区間同士毎に一方から他方の区間およびこの逆方向の移動コストを記憶する手段と、直前の観測時刻までに特定された第１〜第Ｌの観測値系列毎に、前記観測値系列の最新の観測値を含む区間からの合計移動コストが所定コスト以下である次の時刻の観測値群を第１〜第Ｌの観測値群として選択する手段と、前記第１〜第Ｌの観測値系列に対して第１〜第Ｌの観測値群からそれぞれ異なる１つの第１〜第Ｌの観測値を選択して割り当てる複数のパターンを生成する手段と、前記第１〜第Ｌの観測値系列の前記パターンの各観測値に対する尤度の合計をパターン毎に計算する手段と、尤度の合計が最も高いパターンを選択する手段とを備える。
【解決手段】本発明の移動軌跡算出装置は、観測エリアを観測して得た、観測時刻と複数の観測値とを含む複数のセンサデータを記憶する手段と、前記観測エリアを分割した複数の区間における隣接区間同士毎に一方から他方の区間およびこの逆方向の移動コストを記憶する手段と、直前の観測時刻までに特定された第１〜第Ｌの観測値系列毎に、前記観測値系列の最新の観測値を含む区間からの合計移動コストが所定コスト以下である次の時刻の観測値群を第１〜第Ｌの観測値群として選択する手段と、前記第１〜第Ｌの観測値系列に対して第１〜第Ｌの観測値群からそれぞれ異なる１つの第１〜第Ｌの観測値を選択して割り当てる複数のパターンを生成する手段と、前記第１〜第Ｌの観測値系列の前記パターンの各観測値に対する尤度の合計をパターン毎に計算する手段と、尤度の合計が最も高いパターンを選択する手段とを備える。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、観測エリアにおける複数の移動体の移動軌跡を算出する移動軌跡算出装置および移動軌跡算出方法に関する。
【背景技術】
【０００２】
レーダーによって広域な領域を監視し、移動するターゲットを検知することができるようになった。レーダーによる観測は捕捉される対象の詳細までは把握できないため、ターゲットの数だけの影が写ったスナップショット程度の情報しか得られない。すなわち、走査した時刻、ターゲットの数とそれぞれの位置だけが情報として入手できる。これらを時間順に並べて、隣接するスナップショット間で観測値の対応付けを行うことで各ターゲットの軌跡をつかむことができる。軌跡が得られると次時刻の移動場所を予測するなどの追跡が実現できる。複数移動体の追跡では大別して２種類の方法がこれまでに提案されている。
【０００３】
一つはスナップショット間で観測値の組み合わせをつくり、各組み合わせについて統計的な尤度を計算する方法である。統計的に最適な組み合わせを見つけるためには、スナップショットの枚数分の長さの組み合わせをすべて作り、各組み合わせの尤度計算をして最大のものを見つける必要があり、スナップショットの枚数に対して指数オーダで組み合わせが増加する。最初の数枚のスナップショットで組み合わせを作り、スナップショット１枚を追加するたびに尤度計算を行い、そこまでの尤度が絶対的あるいは相対的に低いものを順次捨てていくことで計算量増加を回避する方法（従来法１）も提案されている（非特許文献１参照）。しかしながら、この方法では、尤度計算のベースとなるカルマンフィルタがターゲットの挙動を良く表していないと途中で最適解が捨てられてしまう可能性がある。
【０００４】
たとえば時刻1で観測値11,12,13が、時刻2で観測値21,22,23が,時刻3で観測値31,32,33が得られているとする。まず３×３×３でできる27通りの系列それぞれについて尤度を計算する。この後、次の時刻に３つの観測値が得られると27×3通りの系列ができてしまい、系列長が長くなると組み合わせの数が指数オーダで増加することから、尤度が低い系列を途中でふるいにかけることにより、保持する系列が一定数を超えないようする。このとき、系列の尤度は、カルマンフィルタを用いて予測される次時刻での出現場所と新たな観測値の位置のずれが大きいほど低い値となるため、予測モデルが現実と合っていないと正しい系列が途中でふるい落とされる可能性がある。
【０００５】
二つ目の方法は計算量増加回避のために確率やシミュレーションをベースとした方法（従来法２）である（非特許文献２、３参照）。この方法は観測値間の対応付けを１対１にせず、確率を与えることによって上述の組み合わせ問題を回避し、計算量を一定値に収めるが、ターゲットの移動軌跡を地図上に表示するなどの、「軌跡を提示する」必要が発生した場合には、結局組み合わせ問題が発生することになる。
【０００６】
従来法２では従来法１のような組み合わせを作らず、各時刻での観測値と移動体の間の関連の強さを定量化した表を作成する。移動体同士が十分離れている場合には、表の行あるいは列ごとの最大値だけに注目すると移動体と観測値に１対１対応を得ることができ、簡単に軌跡を表示することができる。しかし、移動体同士が接近すると、最大値だけに着目しても必ずしも１対１の対応を作ることができなくなる。もっともらしい対応付けを、時刻毎の表を組み合わせて考えることになると、結局従来法１と同様の組み合わせ問題が生じてしまう。
【０００７】
上記の従来技術の他にも、非特許文献４にはＶＳ（Variable Structure）−ＩＭＭ（Interacting Multiple Model）推定技術を利用した予測方法が開示されている。このＶＳ−ＩＭＭ推定がベースとするＩＭＭ推定では、ターゲットの移動予測をカルマンフィルタによって行う場合に、複数のモデル（カルマンフィルタを規定する状態遷移関数、システムノイズ、観測関数、観測ノイズを合わせたもの）を用意し、各モデルの信用度に応じた荷重平均によって予測値を得る。ＶＳ−ＩＭＭ推定では、ＩＭＭ推定において複数のモデルを何らかの条件（たとえば地理的条件）によって動的に切り替える。この方法を用いた場合でも、結局、多くのモデルを使ったＩＭＭ推定を行う必要があるため、計算量が多くなる問題がある。
【先行技術文献】
【非特許文献】
【０００８】
【非特許文献１】Donald B.Reid: An Algorithm for Tracking Multiple Targets, IEEE trans. on Automatic Control, Vol.AC-24,No.6,1979
【非特許文献２】Yaakov Bar-Shalom: Tracking methods in a Multitarget Environment, IEEE trans. on Automatic control vol,ac-23,No.4,1978
【非特許文献３】C.Hue, et al: Tracking Multiple Objects with Particle Filtering, IEEE trans. on Aerospace and Electronic systems, vol.38,No.3,2002.
【非特許文献４】Yaakov Bar-shalom, William Dale Blair 編集：Multitarget-Multisensor Tracking: Applications and Advances Vol.III, p.259-319に記載のChapter 6“Large-Scale Ground Target Tracking With Single and Multiple MTI sensors”, Artech House, Inc., 2000(ISBN 1-58053-091-5)
【発明の概要】
【発明が解決しようとする課題】
【０００９】
本発明は、少ない計算量かつ高精度で、複数の移動体の軌跡を算出できるようにした移動軌跡算出装置および移動軌跡算出方法を提供する。
【課題を解決するための手段】
【００１０】
本発明の一態様としての移動軌跡算出装置は、
観測エリアを移動する複数の移動体の移動軌跡を算出する移動軌跡算出装置であって、
前記観測エリアを位置観測センサにより一定時間毎に観測することによって得られた、それぞれ観測時刻と複数の観測値とを含む複数のセンサデータを記憶するセンサデータ記憶手段と、
前記観測エリアを分割した複数の区間における隣接区間同士毎に一方の区間から他方の区間および前記他方の区間から前記一方の区間への移動コストを記憶する第１の記憶手段と、
直前の観測時刻までに特定された第１〜第Ｌの観測値系列に対してそれぞれ個別に、次の観測時刻における前記センサデータの前記複数の観測値のうち、前記観測値系列の最新の観測値を含む区間から隣接区間を順次移動した場合の合計移動コストが所定コスト以下で移動可能な区間に含まれる観測値群を第１〜第Ｌの観測値群として選択する選択手段と、
前記第１〜第Ｌの観測値系列に対して前記第１〜第Ｌの観測値群からそれぞれ異なる１つの第１〜第Ｌの観測値を選択して割り当てる複数の割り当てパターンを生成する割り当てパターン生成手段と、
前記複数の割り当てパターンのそれぞれ毎に、前記第１〜第Ｌの観測値系列の、前記第１〜第Ｌの観測値に対する尤度を計算する尤度計算手段と、
前記第１〜第Ｌの観測値系列に対する尤度の合計が最も高いまたは閾値以上である１つの割り当てパターンを前記複数の割り当てパターンの中から選択し、選択した割り当てパターンに示される第１〜第Ｌの観測値を前記第１〜第Ｌの観測値系列の末尾に追加する観測値系列生成手段と
を有する軌跡算出手段と、
前記第１〜第Ｌの観測値系列を表示する軌跡表示手段と
を備える。
【発明の効果】
【００１１】
本発明により、少ない計算量かつ高精度で、複数の移動体の軌跡を算出できる。
【図面の簡単な説明】
【００１２】
【図１】本発明の一実施形態としての移動軌跡算出装置の構成を概略的に示すブロック図。
【図２】センサデータの一例を示す図。
【図３】図２のセンサデータを時間軸に沿って表示したグラフを示す図。
【図４】平地移動コストの説明図。
【図５】平地移動コストの別の説明図。
【図６】隣接する８方向への移動コストを計算する様子を示す図。
【図７】移動軌跡算出部の処理の流れを示すフローチャート。
【図８】軌跡の初期化処理の流れを示すフローチャート。
【図９】図７のステップＳ５の詳細な計算手順を示すフローチャート。
【図１０】図７のステップＳ６の計算手順の詳細を示すフローチャート。
【図１１】移動可能範囲を求める具体例を説明する図。
【図１２】求めたL本の軌跡群を軌跡表示部に表示した例を示す図。
【図１３】軌跡算出部の詳細ブロック図。
【図１４】本発明の実施形態に係わるハードウェア構成図。
【発明を実施するための形態】
【００１３】
図１は本発明の一実施形態としての移動軌跡算出装置の構成を概略的に示すブロック図である。
【００１４】
この移動軌跡算出装置は、センサデータ記憶部１０１、地理情報記憶部１０２、地理制約計算部１０３，軌跡算出部１０４および軌跡表示部１０５を備える。軌跡算出部１０４は、図１３に示すように選択手段２０１、割り当てパターン生成手段２０２、尤度計算手段２０３、観測値系列生成手段２０４、内部バッファ２０５、読み出し手段２０６を備える。
【００１５】
図１の処理部１０３〜１０５はハードウェアによって構成しても、プログラムモジュールによって構成してもよい。 処理部１０３〜１０５をプログラムモジュールで構成する場合のハードウェア構成の一例を図１４に示す。処理部１０３〜１０５の各プログラムモジュールはハードディスク３０３、ＣＤ−ＲＯＭ３０４またはＵＳＢメモリ３０５等の記録媒体に格納され、ＣＰＵ３０１等のコンピュータにより、当該記録媒体から読み出され、バス３０６を介してＲＡＭ等のメモリ３０２に展開されて実行される。図１４のハードウェア構成の場合、センサデータ記憶部１０１および地理情報記憶部１０２はたとえば、ハードディスク３０３、ＣＤ−ＲＯＭ３０４、またはＵＳＢメモリ３０５等の記録媒体によって構成されることができる。
【００１６】
センサデータ記憶部１０１は、複数の移動体（ターゲット）を含む観測エリアをレーダーなどの位置観測センサにより一定時間毎に観測して得られた複数のセンサデータを記憶している。ここでは移動体は観測エリアとして地表を移動しているとする。各センサデータは複数の観測値（位置データ）と観測時刻とを含む。各センサデータに、観測値の個数情報をさらに含めてもよい。
【００１７】
センサデータの一例を図２に示す。ここで、ｚａ（ｂ）は、時刻ｂで観測されたａ番目の位置を表す。また図２のセンサデータを観測エリアの座標にマッピングして時間軸に沿って表示したグラフを図３に示す。
【００１８】
図１の構成では位置観測センサをその構成要素として含んでいないが、位置観測センサを本発明の構成要素として含んでもよい。この場合地理情報記憶部１０２は位置観測センサにより一定時間毎に観測されるセンサデータを逐次内部に記憶する。
【００１９】
センサデータは軌跡算出部１０４の読み出し手段２０６によって読み出されることができる。
【００２０】
地理情報記憶部１０２は、観測エリア全体の地理情報を記憶している。地理情報記憶部１０２は、地表を矩形のメッシュ状に区切った複数の矩形（区間）の地理情報を記憶している。地理情報は、たとえば各矩形のそれぞれについて、当該矩形の平均標高、隣接する８個の各矩形の方向（後述する図６参照）に対する斜度などを含む。
【００２１】
例えば、国土交通省が公開している標高1kmメッシュデータでは、地表をメッシュ状に複数の正方形（区間）に分割し、正方形ごとに平均標高、平均斜度とその方向を保持している。この国土交通省のデータを地理情報として利用することも可能である。
【００２２】
矩形サイズ（区間サイズ）は、観測エリアの広さ、移動体の移動速度、位置観測センサのサンプリング間隔・精度、軌跡表示部１０５に表示するエリアサイズによって変えても良い。また矩形サイズ毎に地理情報を用意してもよい。
【００２３】
ここで本実施形態における地理情報は、隣接矩形同士間の平地移動コストを含んでいる。平地移動コストは、地表の勾配を無視した場合に、隣接矩形同士間の移動の難しさ（たとえば時間）を表すものである。ここで隣接矩形同士の一方から他方、他方から一方への平地移動コストはそれぞれ同一であっても、それぞれ別個に異なる値としてもよい。以降の説明では両者間の平地移動コストはそれぞれ同一であるとする。
【００２４】
平地移動コストの例を示すと、隣接矩形同士のうちの一方の矩形の中心座標が河川ポリゴン外に位置し、他方の矩形の中心座標が河川ポリゴン内に位置するなら河川を渡る移動が発生するため高い平地移動コストを与える。この場合の一例を図４に示す。中心の矩形からその右の矩形への移動では河川を渡る移動が発生するため高めの平地移動コストが与えられる。
【００２５】
また、隣接矩形同士間に道路があるなら道路上での移動が可能であると判断され低めの平地移動コストを与える。この場合の一例を図５に示す。中心の矩形からその右の矩形への移動では道路があるため低めの平地移動コストが与えられる。
【００２６】
なお、道路および河川のポリゴン座標、あるいは河川の中心線座標等の必要な座標情報はGISを利用することで得ることができる。
【００２７】
このようにして各種の地理条件を考慮して隣接矩形同士間の平地移動コストをあらかじめ定義して、地理情報に含めておく。
【００２８】
地理情報記憶部１０２は、平地移動コスト、平均標高、斜度を記憶する本発明の第２の記憶手段を含む。
【００２９】
地理制約計算部１０３は、地理情報記憶部１０２における地理情報を用いて、隣接矩形同士毎に、一方の区間から他方の区間、他方の区間から一方の区間への移動コストをそれぞれ計算する。観測エリアを矩形に区切った場合、図６に示すように、各矩形について、隣接する８方向への移動コストを計算する。移動コストは、隣接する矩形への移動の難しさ（たとえば時間）を、傾斜も考慮して表したものである。
【００３０】
矩形Ａから、隣接する矩形Ｂへの移動コスト（ｘ（Ａ、Ｂ））はたとえば以下の式により定義することができる。
ｘ（Ａ，Ｂ）＝矩形Ａ、Ｂ間の平地移動コスト＋α×矩形Ｂの方向に対する矩形Ａの斜度
ここでαは、
（１）矩形Ｂの平均標高−矩形Ａの平均標高＞第１所定値、が満たされるとき（傾斜登り方向に一致するとき）１、
（２）矩形Ａの平均標高−矩形Ｂの平均標高＞第２所定値、が満たされるとき（下り方向に一致するとき）−１、
（３）矩形Ｂの平均標高−矩形Ａの平均標高≦第１所定値、または、
矩形Ａの平均標高−矩形Ｂの平均標高≦第２所定値、が満たされるとき（傾斜と直交する移動等、ほぼ同じ高さでの移動とき）０
となるパラメータである。第１所定値および第２所定値はそれぞれ０以上の任意の実数である。
【００３１】
地理制約計算部１０３は、隣接矩形同士毎に計算した移動コストを内部の情報記憶手段あるいは地理情報記憶部１０２に記憶する。したがって、本発明の第１の記憶手段はたとえば情報記憶手段または地理情報記憶手段に対応する。移動コストのデータは軌跡算出部１０４の読み出し手段２０６によって読み出されることができる。
【００３２】
軌跡算出部１０４は、センサデータ記憶部１０１に記憶された各センサデータを用いて、各センサデータ間で観測値の対応付けを行うことにより複数の移動軌跡（観測値系列）を算出する。そして、軌跡表示部１０５は、算出された複数の移動軌跡をディスプレイ等の表示デバイスに表示する。以下、軌跡算出部１０４および軌跡表示部１０５の動作について詳細に説明する。
【００３３】
図７は、軌跡算出部１０４の処理の流れを示すフローチャートである。説明の簡単のため、各観測時刻における観測値の個数（すなわち各センサデータに含まれる観測値の個数）はちょうどL個で、ターゲット数（移動体の個数）もL個であるとする。すなわち観測時刻毎の観測値の個数と、ターゲット数は一致しているとする。
【００３４】
まず、最初の時刻t1で観測されたL個の観測値を取得する（Ｓ１）。すなわち、t=1での観測値集合{z1(1),z2(1),...,zL(1)}を得る。
【００３５】
次に、図８に示す軌跡の初期化処理を行う（Ｓ２）。すなわち最初の時刻t1における各観測値を初期軌跡
【数１】
として得る（Ｓ２１）。つまり各観測値をそれぞれ軌跡の始点とする。
【数２】
は、時刻t1〜時刻tkまでの観測値からなるl（エル）番目の軌跡を示す。
【００３６】
次に、t←t+1として、時刻を１時刻インクリメントする（Ｓ３）。
【００３７】
次に、インクリメントされた時刻tでの観測値集合{z1(t),z2(t),..., zL(t)}を得る（Ｓ４）。
【００３８】
次に、L本ある既存軌跡θ（t-1,l）(l=1〜L)のそれぞれについて、時刻tの観測値集合のうち、制約条件（所定コスト）を満たす観測値を候補として選択する（残す）処理を行う（Ｓ５）。つまり、L本ある既存軌跡毎に、上記時刻tの観測値集合のうち、当該既存軌跡の最新の観測値から移動可能範囲に含まれる観測値を候補として選択する。L本ある既存軌跡θは軌跡算出部１０４の内部バッファ２０５に一時的に格納されているものとする。本ステップＳ５は軌跡算出部１０４の選択手段２０１により行う。ステップＳ５の詳細は後述する。
【００３９】
ここでθ（t-1,l）は時刻t-1までのｌ（エル）番目の軌跡を表す。１番目の既存軌跡θ（t-1,1）について選択された候補は本発明の第１の観測値群に対応し、２番目の既存軌跡θ（t-1,2）について選択された候補は本発明の第２の観測値群に対応し、・・・、Ｌ番目の既存軌跡θ（t-1,L）について選択された候補は本発明の第Ｌの観測値群に対応する。また、１〜L番目の既存軌跡は本発明の直前の観測時刻までに特定された第１〜第Lの観測値系列に対応する。
【００４０】
次に、ステップＳ６では、L本の既存軌跡毎に選択された候補群からそれぞれ１つの候補を選択して当該L本の既存軌跡に割り当てる。この際、各既存軌跡にそれぞれ異なる候補（観測値）が割り当たるように（すなわち同じ観測値が複数の既存軌跡に重複して割り当たらないように）各既存軌跡に対する候補をそれぞれ選択する。このような割り当てのパターンを複数生成する（割り当てパターン生成手段２０１）。そしてこれら複数の割り当てパターンの中で評価値（後述する合計尤度Ｂ）が最大あるいは閾値以上となる割り当てパターンを選択し、選択した割り当てパターンに示される各候補（第１〜第Lの観測値）をＬ本の既存軌跡の末尾にそれぞれ追加する（尤度計算手段２０３および観測値系列生成手段２０４）。ステップＳ６の詳細は後述する。
【００４１】
次に、各候補（第１〜第Lの観測値）が追加されたL本の軌跡群（すなわち第１〜第Ｌの観測値系列）を内部バッファ２０５に一時的に格納するとともに、軌跡表示部１０５を介して表示する（Ｓ７）。表示の一例を図１２に示す。
【００４２】
この後、処理を継続する場合は（Ｓ８のＮＯ）、ステップＳ３に戻り、処理を継続しない場合は（Ｓ８のＹＥＳ）本処理を終了する。たとえばすべてのセンサデータを処理したら本処理を終了してもよいし、処理終了を指示するユーザ指示が入力されたら本処理を終了してもよい。
【００４３】
＜ステップＳ５（制約チェック）の詳細＞
図９はステップＳ５の詳細な計算手順を示すフローチャートである。
【００４４】
まず各既存軌跡の末尾の観測値、すなわち、各既存軌跡における直前の時刻t-1の観測値を取り出す（Ｓ３１）。
【００４５】
次に、各既存軌跡のそれぞれ毎に、各取り出した観測値を基準に、現時刻tまでに移動可能な範囲を求める（Ｓ３２）。具体的に、取り出した観測値を含む矩形から隣接矩形を順次辿った場合に合計移動コストが所定コスト以下で移動可能な範囲を求める。
【００４６】
次に、既存軌跡毎に、現時刻tにおける観測値集合のうち、当該移動可能範囲に含まる観測値のみを候補として選択する（Ｓ３３）。すなわち観測値集合のうち、当該移動可能範囲に含まれない観測値を候補から外す。
【００４７】
次に、既存軌跡毎に、それぞれ選択された候補群を返す（Ｓ３４）。
【００４８】
図１１はステップＳ３２で移動可能範囲を求める処理の具体例を説明する図である。
【００４９】
図１１(1)のように２次元座標系の地図（観測エリア）をメッシュ分割し、着目する矩形から、各隣接する矩形への移動コストを(2)のように定義する。これを地図全体に適用すると(3)のようになる。例えば、(3)の矩形Sを起点として、ここからコストC（所定コストに対応）で到達可能な範囲を求める手順は次の通りである。
(a) Sを根とする。
(b) Sから直接移動可能な隣接する８つの矩形を子ノードとして並べる（図１１(4)の１〜８番のノード）。すなわち、Sを根とし、８個の子ノードをもつ木構造を生成する。
(c) 根Sから各子ノード１〜８に至る移動コストをノードに記録する。
(d) さらに各子ノードから移動可能な８つの矩形を孫ノードとして並べ、各孫ノードには、その親となる子ノードに記録された移動コストに、その親からその子へ移動する際にかかる移動コストを加えた値（合計移動コスト）を記録する。
(e) (d)を再帰的に繰り返して木を広げていく。途中、ノードのコスト（合計移動コスト）がCを越える場合には子ノードの展開を停止する。
【００５０】
移動コストがすべて正の値であれば上記手順はやがて停止し、木が求まる。その木に現れる全ノードがコストC（所定コスト）で到達可能な範囲（移動可能範囲）となる。
【００５１】
たとえば移動コストが「移動時間」を表し、直前の観測時刻t-1と現時刻tとの差がCとすると、上記方法で求めた各ノードの観測値の集合が、直前の観測値（座標）から時間C内で移動可能な範囲となり、この範囲に含まれない観測値（座標）へは移動不可能ということになる。つまり、現時刻tの観測値集合のうちこの移動可能範囲に含まれるものを候補として残す。
【００５２】
要するに、上述のステップＳ３２およびＳ３３では、直前の観測時刻t-1までに特定された第１〜第Ｌの観測値系列に対してそれぞれ個別に、次の観測時刻tにおけるセンサデータの複数の観測値のうち、当該観測値系列における最新の観測値を含む区間（矩形）から隣接区間を順次移動した場合の合計移動コストが所定コスト以下で移動可能な区間に含まれる観測値群を、第１〜第Ｌの観測値（候補）群として選択する。合計移動コストは、隣接区間同士間では一方の区間から他方の区間への移動コストそのものになる。また直前の観測時刻t-1の観測値と、最新の観測値とが同一区間に含まれる場合の合計移動コストは処理の簡単のためたとえば単純にゼロとする（いったん隣接区間へ移動し再度元の区間へ戻る場合も観測エリアの分割サイズによっては想定されるがここではそのようなケースは考えないものとする）。この段落で述べた処理は軌跡算出部１０４が備える選択手段２０１により行われる。
【００５３】
別の言い方をすれば、上述のステップＳ３２およびＳ３３では、直前の観測時刻t-1までに特定された第１〜第Ｌの観測値系列毎に、末尾の観測値を含む区間から、次の観測時刻tにおける各観測値を含む区間のそれぞれに対して、移動コストが所定コスト以下で行くことが可能であるような経路が存在するかどうかを調べ、そのような経路が存在する観測値については候補として残す。
【００５４】
＜ステップＳ６（軌跡生成）の詳細＞
図１０はステップＳ６の計算手順の詳細を示すフローチャートである。
【００５５】
まず変数Aに十分小さな値を代入する（Ｓ４１）。
【００５６】
次に、L本ある既存軌跡θ（t-1,l）(l=1〜L)のそれぞれ毎に、ステップＳ５で得られた観測値（候補）のうちの１つを、同一の観測値（候補）が複数の既存軌跡に含まれないように選択して割り当てる（Ｓ４２）。これが１つの割り当てパターン（候補の組み合わせ）になる。すなわち、第１〜第Ｌの既存軌跡に対して第１〜第Ｌの観測値群からそれぞれ異なる１つの第１〜第Ｌの観測値（候補）を選択して割り当てる割り当てパターンを複数生成する。この処理は軌跡算出手段１０４が備える割り当てパターン生成手段２０２の処理に対応する。
【００５７】
次に、各既存軌跡（各観測値系列）から、上記割り当てパターンに示されるそれぞれの観測値（候補）への尤度を計算し、各計算した尤度の和B（評価値）を求める（Ｓ４３）。つまり、第１〜第Ｌの観測値系列の、上記割り当てパターンに示される第１〜第Ｌの観測値に対する尤度を計算し、各計算した尤度の和を求める。この尤度の計算は、軌跡算出手段１０４が備える尤度計算手段２０３により行われる。尤度の計算方法の詳細は後述する。
【００５８】
尤度和Bが変数Aより大きいか否かを検査し（Ｓ４４）、大きい場合は（ＹＥＳ）、Bの値を変数Aに代入し、この割り当てパターンを記録する（Ｓ４５）。
【００５９】
次に、生成可能なすべての割り当てパターンを生成したか否かを検査し（Ｓ４６）、まだ生成していない割り当てパターンが存在するときは（ＮＯ）ステップＳ４２に戻り、すべての割り当てパターンを生成済みのときは（ＹＥＳ）、ステップＳ４５で最後に記録された割り当てパターンを返す（Ｓ４７）。すなわち尤度和Bが最も大きいときの割り当てパターンを返す。ここでは尤度和Bが最も大きいときの割り当てパターンを返しているが、閾値以上の尤度和をもつ割り当てパターンのうちの１つ選択して返しても良い。つまり、尤度和が最も高いまたは閾値以上である１つの割り当てパターンを、ステップＳ４２で生成した複数の割り当てパターンの中から選択し、選択した割り当てパターンに示される各観測値（第１〜第Ｌの観測値）を第１〜第Ｌの観測値系列の末尾に追加する。この処理は、軌跡算出手段１０４が備える観測値系列生成手段２０４により行われる。
【００６０】
＜ステップＳ４３（尤度計算）の詳細＞
ここでは尤度の計算方法としてカルマンフィルタを用いる例を示す。
【００６１】
まずカルマンフィルタについて説明する。
【００６２】
ターゲット（移動体）の挙動は例えば線形カルマンフィルタで記述できる。
【００６３】
カルマンフィルタではターゲットの状態変数x(t)を
【数３】
と記述する。Fは状態遷移行列であり、Qはシステムノイズである。状態変数は位置座標および速度からなるベクトル変数である。また、状態変数から観測値を得るために
【数４】
を用いる。Hは観測行列、Rは観測ノイズである。
【００６４】
観測されるものがxとy座標、すなわちターゲットの位置のみであるとすると、観測値は
【数５】
となる。a,bはそれぞれターゲットのｘ座標、ｙ座標である。ターゲットの移動軌跡は事前に分からないため、等速度運動または等加速度運動を仮定し、その大きさおよび方向は観測値から推定する。等速度運動を仮定するなら状態変数は、位置(a,b)と、速度
【数６】
とからなる。
【数７】
となる。
【数８】
である。q1,q2,r1,r2は平均０の白色雑音である。
【００６５】
以下では簡単のため一部、時刻tを省略して式を書く。例えばx(tk)をxkと書く。
【００６６】
カルマンフィルタによる処理は、時刻tkにおける観測値z(tk)が得られる前の「状態変数の事前推定値
【数９】
」と観測値が得られた後の「事後推定値
【数１０】
」を交互に求めながら進められる。
【００６７】
［カルマンフィルタの処理手順］
【数１１】
【００６８】
マルチターゲット追跡ではある時刻tにおいて複数の位置にターゲットが観測される。観測時刻をti(i=1,…,n)、このときに得られた観測座標をzj(ti)(j=1,…,mi)とする。簡単のために各時刻で観測されるターゲット数が一定数ｍであるとすると、時刻t1〜tnの間に合計でｍ×n個の観測値が得られることになる。時刻t1における観測値集合から適当に一つ選びzl1(t1)とする。同様に時刻t2における観測値集合から一つ選びzl2(t2)とする。これを時刻tkまで繰り返すと、zl1(t1),zl2(t2),zl3(t3),....,zlk(tk)という軌跡が作られることになる。このような軌跡をL個作ることにし、各軌跡を
【数１２】
と書くことにすると、
【数１３】
は次式となる。
【数１４】
【００６９】
（式２）の尤度の対数に負号をつけたものを
【数１５】
とする。尤度Λは値が大きいほど「尤もらしい（確度が高い）」ことを表しているので、λは値が小さいほど尤もらしいことを表していることになる。λは
【数１６】
と書ける。
【００７０】
つまり、時刻t1から時刻tk-1までの観測値から作られた軌跡
【数１７】
に、新たに時刻tkでの観測値を１つ加えて、軌跡を１つ伸ばした時のλは、伸ばす前のλに、新たに加えた分の観測値から計算される量を足し加えればよいことが分かる。
【００７１】
したがって、新たに得られた観測値集合から観測値を一つ選び出して軌跡に加える際には、（式３）右辺の第２項が最小のものを選ぶのが良いことになる。ただし、マルチターゲット追跡における複数（ｌ＝１〜L）の軌跡にとって（式３）右辺第２項が最小となる観測値がかち合う場合もあり、ターゲットが互いに接近すると必ずしも尤度が最適となる観測値を割り当てることが正解とはならない。したがって、観測値が得られるたびに、複数の割り当てパターンを考慮して観測値の割り当てを決定する必要がある（図１０のＳ４２参照）。それ故、観測値が得られるたびに、全体における観測値の組み合わせの数が指数的に増えることになるが、本実施形態では、この点、移動可能範囲に含まれる観測値群のみを割り当ての可能性のある候補（観測値）群として扱うため（図７のステップＳ５の説明を参照）、結果として、全体における観測値の組み合わせの数を少なく抑えることが可能である。
【００７２】
以上のように、本実施形態によれば、隣接矩形同士間に地形的な要因に基づく移動コストを定義し、ある時刻の観測値から次の時刻までに移動可能な範囲を移動コストに基づき求めて当該移動可能な範囲に含まれる観測値のみを次の時刻の移動先候補として選択することにより、尤度計算において各既存軌跡（観測値系列）に対する割り当てパターンの数（組み合わせパターン数）を減少させることができるとともに、最も良い割り当てパターンを選択できる。
【００７３】
従来であれば、たとえば１時刻前までに得られた既存軌跡がL個（θ１〜θL）あり、現時刻で観測された観測値がL個（ｚ１〜zＬ）あるとすると、１：１の条件だけを満たすL！通りの組み合わせ（割り当てパターン）
[(θ1,z1) (θ2,z2) (θ3,z3),..., (θL,zL)]、
[(θ1,z2) (θ2,z1) (θ3,z3),..., (θL,zL)]、
・・・
[(θ1,zL) (θ2,zL-1) (θ3,zL-2),..., (θL,z1)]
全部について尤度和を計算し、最も良い尤度和をもつ割り当てパターンを選ぶことになる。したがって計算量が膨大になり、また地形（移動コスト）を考慮していないため実際には不可能な移動を含む割り当てパターンが選択されることもあり得る。
【００７４】
これに対し本実施形態では、地理条件（移動コスト）を元に、対応可能なθとｚの組み合わせを制限しているため、制約を満たさぬ組み合わせ（割り当てパターン）を除いた残りについて尤度和の計算を行い、最も良い組み合わせ（割り当てパターン）を選ぶことができる。
【符号の説明】
【００７５】
１０１：センサデータ記憶部
１０２：地理情報記憶部（第１の記憶手段、第２の記憶手段）
１０３：地理制約計算部（第２の記憶手段）
１０４：軌跡算出部（選択手段、割り当てパターン生成手段、尤度計算手段、観測値系列生成手段、読み出し手段）
１０５：軌跡表示部
【技術分野】
【０００１】
本発明は、観測エリアにおける複数の移動体の移動軌跡を算出する移動軌跡算出装置および移動軌跡算出方法に関する。
【背景技術】
【０００２】
レーダーによって広域な領域を監視し、移動するターゲットを検知することができるようになった。レーダーによる観測は捕捉される対象の詳細までは把握できないため、ターゲットの数だけの影が写ったスナップショット程度の情報しか得られない。すなわち、走査した時刻、ターゲットの数とそれぞれの位置だけが情報として入手できる。これらを時間順に並べて、隣接するスナップショット間で観測値の対応付けを行うことで各ターゲットの軌跡をつかむことができる。軌跡が得られると次時刻の移動場所を予測するなどの追跡が実現できる。複数移動体の追跡では大別して２種類の方法がこれまでに提案されている。
【０００３】
一つはスナップショット間で観測値の組み合わせをつくり、各組み合わせについて統計的な尤度を計算する方法である。統計的に最適な組み合わせを見つけるためには、スナップショットの枚数分の長さの組み合わせをすべて作り、各組み合わせの尤度計算をして最大のものを見つける必要があり、スナップショットの枚数に対して指数オーダで組み合わせが増加する。最初の数枚のスナップショットで組み合わせを作り、スナップショット１枚を追加するたびに尤度計算を行い、そこまでの尤度が絶対的あるいは相対的に低いものを順次捨てていくことで計算量増加を回避する方法（従来法１）も提案されている（非特許文献１参照）。しかしながら、この方法では、尤度計算のベースとなるカルマンフィルタがターゲットの挙動を良く表していないと途中で最適解が捨てられてしまう可能性がある。
【０００４】
たとえば時刻1で観測値11,12,13が、時刻2で観測値21,22,23が,時刻3で観測値31,32,33が得られているとする。まず３×３×３でできる27通りの系列それぞれについて尤度を計算する。この後、次の時刻に３つの観測値が得られると27×3通りの系列ができてしまい、系列長が長くなると組み合わせの数が指数オーダで増加することから、尤度が低い系列を途中でふるいにかけることにより、保持する系列が一定数を超えないようする。このとき、系列の尤度は、カルマンフィルタを用いて予測される次時刻での出現場所と新たな観測値の位置のずれが大きいほど低い値となるため、予測モデルが現実と合っていないと正しい系列が途中でふるい落とされる可能性がある。
【０００５】
二つ目の方法は計算量増加回避のために確率やシミュレーションをベースとした方法（従来法２）である（非特許文献２、３参照）。この方法は観測値間の対応付けを１対１にせず、確率を与えることによって上述の組み合わせ問題を回避し、計算量を一定値に収めるが、ターゲットの移動軌跡を地図上に表示するなどの、「軌跡を提示する」必要が発生した場合には、結局組み合わせ問題が発生することになる。
【０００６】
従来法２では従来法１のような組み合わせを作らず、各時刻での観測値と移動体の間の関連の強さを定量化した表を作成する。移動体同士が十分離れている場合には、表の行あるいは列ごとの最大値だけに注目すると移動体と観測値に１対１対応を得ることができ、簡単に軌跡を表示することができる。しかし、移動体同士が接近すると、最大値だけに着目しても必ずしも１対１の対応を作ることができなくなる。もっともらしい対応付けを、時刻毎の表を組み合わせて考えることになると、結局従来法１と同様の組み合わせ問題が生じてしまう。
【０００７】
上記の従来技術の他にも、非特許文献４にはＶＳ（Variable Structure）−ＩＭＭ（Interacting Multiple Model）推定技術を利用した予測方法が開示されている。このＶＳ−ＩＭＭ推定がベースとするＩＭＭ推定では、ターゲットの移動予測をカルマンフィルタによって行う場合に、複数のモデル（カルマンフィルタを規定する状態遷移関数、システムノイズ、観測関数、観測ノイズを合わせたもの）を用意し、各モデルの信用度に応じた荷重平均によって予測値を得る。ＶＳ−ＩＭＭ推定では、ＩＭＭ推定において複数のモデルを何らかの条件（たとえば地理的条件）によって動的に切り替える。この方法を用いた場合でも、結局、多くのモデルを使ったＩＭＭ推定を行う必要があるため、計算量が多くなる問題がある。
【先行技術文献】
【非特許文献】
【０００８】
【非特許文献１】Donald B.Reid: An Algorithm for Tracking Multiple Targets, IEEE trans. on Automatic Control, Vol.AC-24,No.6,1979
【非特許文献２】Yaakov Bar-Shalom: Tracking methods in a Multitarget Environment, IEEE trans. on Automatic control vol,ac-23,No.4,1978
【非特許文献３】C.Hue, et al: Tracking Multiple Objects with Particle Filtering, IEEE trans. on Aerospace and Electronic systems, vol.38,No.3,2002.
【非特許文献４】Yaakov Bar-shalom, William Dale Blair 編集：Multitarget-Multisensor Tracking: Applications and Advances Vol.III, p.259-319に記載のChapter 6“Large-Scale Ground Target Tracking With Single and Multiple MTI sensors”, Artech House, Inc., 2000(ISBN 1-58053-091-5)
【発明の概要】
【発明が解決しようとする課題】
【０００９】
本発明は、少ない計算量かつ高精度で、複数の移動体の軌跡を算出できるようにした移動軌跡算出装置および移動軌跡算出方法を提供する。
【課題を解決するための手段】
【００１０】
本発明の一態様としての移動軌跡算出装置は、
観測エリアを移動する複数の移動体の移動軌跡を算出する移動軌跡算出装置であって、
前記観測エリアを位置観測センサにより一定時間毎に観測することによって得られた、それぞれ観測時刻と複数の観測値とを含む複数のセンサデータを記憶するセンサデータ記憶手段と、
前記観測エリアを分割した複数の区間における隣接区間同士毎に一方の区間から他方の区間および前記他方の区間から前記一方の区間への移動コストを記憶する第１の記憶手段と、
直前の観測時刻までに特定された第１〜第Ｌの観測値系列に対してそれぞれ個別に、次の観測時刻における前記センサデータの前記複数の観測値のうち、前記観測値系列の最新の観測値を含む区間から隣接区間を順次移動した場合の合計移動コストが所定コスト以下で移動可能な区間に含まれる観測値群を第１〜第Ｌの観測値群として選択する選択手段と、
前記第１〜第Ｌの観測値系列に対して前記第１〜第Ｌの観測値群からそれぞれ異なる１つの第１〜第Ｌの観測値を選択して割り当てる複数の割り当てパターンを生成する割り当てパターン生成手段と、
前記複数の割り当てパターンのそれぞれ毎に、前記第１〜第Ｌの観測値系列の、前記第１〜第Ｌの観測値に対する尤度を計算する尤度計算手段と、
前記第１〜第Ｌの観測値系列に対する尤度の合計が最も高いまたは閾値以上である１つの割り当てパターンを前記複数の割り当てパターンの中から選択し、選択した割り当てパターンに示される第１〜第Ｌの観測値を前記第１〜第Ｌの観測値系列の末尾に追加する観測値系列生成手段と
を有する軌跡算出手段と、
前記第１〜第Ｌの観測値系列を表示する軌跡表示手段と
を備える。
【発明の効果】
【００１１】
本発明により、少ない計算量かつ高精度で、複数の移動体の軌跡を算出できる。
【図面の簡単な説明】
【００１２】
【図１】本発明の一実施形態としての移動軌跡算出装置の構成を概略的に示すブロック図。
【図２】センサデータの一例を示す図。
【図３】図２のセンサデータを時間軸に沿って表示したグラフを示す図。
【図４】平地移動コストの説明図。
【図５】平地移動コストの別の説明図。
【図６】隣接する８方向への移動コストを計算する様子を示す図。
【図７】移動軌跡算出部の処理の流れを示すフローチャート。
【図８】軌跡の初期化処理の流れを示すフローチャート。
【図９】図７のステップＳ５の詳細な計算手順を示すフローチャート。
【図１０】図７のステップＳ６の計算手順の詳細を示すフローチャート。
【図１１】移動可能範囲を求める具体例を説明する図。
【図１２】求めたL本の軌跡群を軌跡表示部に表示した例を示す図。
【図１３】軌跡算出部の詳細ブロック図。
【図１４】本発明の実施形態に係わるハードウェア構成図。
【発明を実施するための形態】
【００１３】
図１は本発明の一実施形態としての移動軌跡算出装置の構成を概略的に示すブロック図である。
【００１４】
この移動軌跡算出装置は、センサデータ記憶部１０１、地理情報記憶部１０２、地理制約計算部１０３，軌跡算出部１０４および軌跡表示部１０５を備える。軌跡算出部１０４は、図１３に示すように選択手段２０１、割り当てパターン生成手段２０２、尤度計算手段２０３、観測値系列生成手段２０４、内部バッファ２０５、読み出し手段２０６を備える。
【００１５】
図１の処理部１０３〜１０５はハードウェアによって構成しても、プログラムモジュールによって構成してもよい。 処理部１０３〜１０５をプログラムモジュールで構成する場合のハードウェア構成の一例を図１４に示す。処理部１０３〜１０５の各プログラムモジュールはハードディスク３０３、ＣＤ−ＲＯＭ３０４またはＵＳＢメモリ３０５等の記録媒体に格納され、ＣＰＵ３０１等のコンピュータにより、当該記録媒体から読み出され、バス３０６を介してＲＡＭ等のメモリ３０２に展開されて実行される。図１４のハードウェア構成の場合、センサデータ記憶部１０１および地理情報記憶部１０２はたとえば、ハードディスク３０３、ＣＤ−ＲＯＭ３０４、またはＵＳＢメモリ３０５等の記録媒体によって構成されることができる。
【００１６】
センサデータ記憶部１０１は、複数の移動体（ターゲット）を含む観測エリアをレーダーなどの位置観測センサにより一定時間毎に観測して得られた複数のセンサデータを記憶している。ここでは移動体は観測エリアとして地表を移動しているとする。各センサデータは複数の観測値（位置データ）と観測時刻とを含む。各センサデータに、観測値の個数情報をさらに含めてもよい。
【００１７】
センサデータの一例を図２に示す。ここで、ｚａ（ｂ）は、時刻ｂで観測されたａ番目の位置を表す。また図２のセンサデータを観測エリアの座標にマッピングして時間軸に沿って表示したグラフを図３に示す。
【００１８】
図１の構成では位置観測センサをその構成要素として含んでいないが、位置観測センサを本発明の構成要素として含んでもよい。この場合地理情報記憶部１０２は位置観測センサにより一定時間毎に観測されるセンサデータを逐次内部に記憶する。
【００１９】
センサデータは軌跡算出部１０４の読み出し手段２０６によって読み出されることができる。
【００２０】
地理情報記憶部１０２は、観測エリア全体の地理情報を記憶している。地理情報記憶部１０２は、地表を矩形のメッシュ状に区切った複数の矩形（区間）の地理情報を記憶している。地理情報は、たとえば各矩形のそれぞれについて、当該矩形の平均標高、隣接する８個の各矩形の方向（後述する図６参照）に対する斜度などを含む。
【００２１】
例えば、国土交通省が公開している標高1kmメッシュデータでは、地表をメッシュ状に複数の正方形（区間）に分割し、正方形ごとに平均標高、平均斜度とその方向を保持している。この国土交通省のデータを地理情報として利用することも可能である。
【００２２】
矩形サイズ（区間サイズ）は、観測エリアの広さ、移動体の移動速度、位置観測センサのサンプリング間隔・精度、軌跡表示部１０５に表示するエリアサイズによって変えても良い。また矩形サイズ毎に地理情報を用意してもよい。
【００２３】
ここで本実施形態における地理情報は、隣接矩形同士間の平地移動コストを含んでいる。平地移動コストは、地表の勾配を無視した場合に、隣接矩形同士間の移動の難しさ（たとえば時間）を表すものである。ここで隣接矩形同士の一方から他方、他方から一方への平地移動コストはそれぞれ同一であっても、それぞれ別個に異なる値としてもよい。以降の説明では両者間の平地移動コストはそれぞれ同一であるとする。
【００２４】
平地移動コストの例を示すと、隣接矩形同士のうちの一方の矩形の中心座標が河川ポリゴン外に位置し、他方の矩形の中心座標が河川ポリゴン内に位置するなら河川を渡る移動が発生するため高い平地移動コストを与える。この場合の一例を図４に示す。中心の矩形からその右の矩形への移動では河川を渡る移動が発生するため高めの平地移動コストが与えられる。
【００２５】
また、隣接矩形同士間に道路があるなら道路上での移動が可能であると判断され低めの平地移動コストを与える。この場合の一例を図５に示す。中心の矩形からその右の矩形への移動では道路があるため低めの平地移動コストが与えられる。
【００２６】
なお、道路および河川のポリゴン座標、あるいは河川の中心線座標等の必要な座標情報はGISを利用することで得ることができる。
【００２７】
このようにして各種の地理条件を考慮して隣接矩形同士間の平地移動コストをあらかじめ定義して、地理情報に含めておく。
【００２８】
地理情報記憶部１０２は、平地移動コスト、平均標高、斜度を記憶する本発明の第２の記憶手段を含む。
【００２９】
地理制約計算部１０３は、地理情報記憶部１０２における地理情報を用いて、隣接矩形同士毎に、一方の区間から他方の区間、他方の区間から一方の区間への移動コストをそれぞれ計算する。観測エリアを矩形に区切った場合、図６に示すように、各矩形について、隣接する８方向への移動コストを計算する。移動コストは、隣接する矩形への移動の難しさ（たとえば時間）を、傾斜も考慮して表したものである。
【００３０】
矩形Ａから、隣接する矩形Ｂへの移動コスト（ｘ（Ａ、Ｂ））はたとえば以下の式により定義することができる。
ｘ（Ａ，Ｂ）＝矩形Ａ、Ｂ間の平地移動コスト＋α×矩形Ｂの方向に対する矩形Ａの斜度
ここでαは、
（１）矩形Ｂの平均標高−矩形Ａの平均標高＞第１所定値、が満たされるとき（傾斜登り方向に一致するとき）１、
（２）矩形Ａの平均標高−矩形Ｂの平均標高＞第２所定値、が満たされるとき（下り方向に一致するとき）−１、
（３）矩形Ｂの平均標高−矩形Ａの平均標高≦第１所定値、または、
矩形Ａの平均標高−矩形Ｂの平均標高≦第２所定値、が満たされるとき（傾斜と直交する移動等、ほぼ同じ高さでの移動とき）０
となるパラメータである。第１所定値および第２所定値はそれぞれ０以上の任意の実数である。
【００３１】
地理制約計算部１０３は、隣接矩形同士毎に計算した移動コストを内部の情報記憶手段あるいは地理情報記憶部１０２に記憶する。したがって、本発明の第１の記憶手段はたとえば情報記憶手段または地理情報記憶手段に対応する。移動コストのデータは軌跡算出部１０４の読み出し手段２０６によって読み出されることができる。
【００３２】
軌跡算出部１０４は、センサデータ記憶部１０１に記憶された各センサデータを用いて、各センサデータ間で観測値の対応付けを行うことにより複数の移動軌跡（観測値系列）を算出する。そして、軌跡表示部１０５は、算出された複数の移動軌跡をディスプレイ等の表示デバイスに表示する。以下、軌跡算出部１０４および軌跡表示部１０５の動作について詳細に説明する。
【００３３】
図７は、軌跡算出部１０４の処理の流れを示すフローチャートである。説明の簡単のため、各観測時刻における観測値の個数（すなわち各センサデータに含まれる観測値の個数）はちょうどL個で、ターゲット数（移動体の個数）もL個であるとする。すなわち観測時刻毎の観測値の個数と、ターゲット数は一致しているとする。
【００３４】
まず、最初の時刻t1で観測されたL個の観測値を取得する（Ｓ１）。すなわち、t=1での観測値集合{z1(1),z2(1),...,zL(1)}を得る。
【００３５】
次に、図８に示す軌跡の初期化処理を行う（Ｓ２）。すなわち最初の時刻t1における各観測値を初期軌跡
【数１】
として得る（Ｓ２１）。つまり各観測値をそれぞれ軌跡の始点とする。
【数２】
は、時刻t1〜時刻tkまでの観測値からなるl（エル）番目の軌跡を示す。
【００３６】
次に、t←t+1として、時刻を１時刻インクリメントする（Ｓ３）。
【００３７】
次に、インクリメントされた時刻tでの観測値集合{z1(t),z2(t),..., zL(t)}を得る（Ｓ４）。
【００３８】
次に、L本ある既存軌跡θ（t-1,l）(l=1〜L)のそれぞれについて、時刻tの観測値集合のうち、制約条件（所定コスト）を満たす観測値を候補として選択する（残す）処理を行う（Ｓ５）。つまり、L本ある既存軌跡毎に、上記時刻tの観測値集合のうち、当該既存軌跡の最新の観測値から移動可能範囲に含まれる観測値を候補として選択する。L本ある既存軌跡θは軌跡算出部１０４の内部バッファ２０５に一時的に格納されているものとする。本ステップＳ５は軌跡算出部１０４の選択手段２０１により行う。ステップＳ５の詳細は後述する。
【００３９】
ここでθ（t-1,l）は時刻t-1までのｌ（エル）番目の軌跡を表す。１番目の既存軌跡θ（t-1,1）について選択された候補は本発明の第１の観測値群に対応し、２番目の既存軌跡θ（t-1,2）について選択された候補は本発明の第２の観測値群に対応し、・・・、Ｌ番目の既存軌跡θ（t-1,L）について選択された候補は本発明の第Ｌの観測値群に対応する。また、１〜L番目の既存軌跡は本発明の直前の観測時刻までに特定された第１〜第Lの観測値系列に対応する。
【００４０】
次に、ステップＳ６では、L本の既存軌跡毎に選択された候補群からそれぞれ１つの候補を選択して当該L本の既存軌跡に割り当てる。この際、各既存軌跡にそれぞれ異なる候補（観測値）が割り当たるように（すなわち同じ観測値が複数の既存軌跡に重複して割り当たらないように）各既存軌跡に対する候補をそれぞれ選択する。このような割り当てのパターンを複数生成する（割り当てパターン生成手段２０１）。そしてこれら複数の割り当てパターンの中で評価値（後述する合計尤度Ｂ）が最大あるいは閾値以上となる割り当てパターンを選択し、選択した割り当てパターンに示される各候補（第１〜第Lの観測値）をＬ本の既存軌跡の末尾にそれぞれ追加する（尤度計算手段２０３および観測値系列生成手段２０４）。ステップＳ６の詳細は後述する。
【００４１】
次に、各候補（第１〜第Lの観測値）が追加されたL本の軌跡群（すなわち第１〜第Ｌの観測値系列）を内部バッファ２０５に一時的に格納するとともに、軌跡表示部１０５を介して表示する（Ｓ７）。表示の一例を図１２に示す。
【００４２】
この後、処理を継続する場合は（Ｓ８のＮＯ）、ステップＳ３に戻り、処理を継続しない場合は（Ｓ８のＹＥＳ）本処理を終了する。たとえばすべてのセンサデータを処理したら本処理を終了してもよいし、処理終了を指示するユーザ指示が入力されたら本処理を終了してもよい。
【００４３】
＜ステップＳ５（制約チェック）の詳細＞
図９はステップＳ５の詳細な計算手順を示すフローチャートである。
【００４４】
まず各既存軌跡の末尾の観測値、すなわち、各既存軌跡における直前の時刻t-1の観測値を取り出す（Ｓ３１）。
【００４５】
次に、各既存軌跡のそれぞれ毎に、各取り出した観測値を基準に、現時刻tまでに移動可能な範囲を求める（Ｓ３２）。具体的に、取り出した観測値を含む矩形から隣接矩形を順次辿った場合に合計移動コストが所定コスト以下で移動可能な範囲を求める。
【００４６】
次に、既存軌跡毎に、現時刻tにおける観測値集合のうち、当該移動可能範囲に含まる観測値のみを候補として選択する（Ｓ３３）。すなわち観測値集合のうち、当該移動可能範囲に含まれない観測値を候補から外す。
【００４７】
次に、既存軌跡毎に、それぞれ選択された候補群を返す（Ｓ３４）。
【００４８】
図１１はステップＳ３２で移動可能範囲を求める処理の具体例を説明する図である。
【００４９】
図１１(1)のように２次元座標系の地図（観測エリア）をメッシュ分割し、着目する矩形から、各隣接する矩形への移動コストを(2)のように定義する。これを地図全体に適用すると(3)のようになる。例えば、(3)の矩形Sを起点として、ここからコストC（所定コストに対応）で到達可能な範囲を求める手順は次の通りである。
(a) Sを根とする。
(b) Sから直接移動可能な隣接する８つの矩形を子ノードとして並べる（図１１(4)の１〜８番のノード）。すなわち、Sを根とし、８個の子ノードをもつ木構造を生成する。
(c) 根Sから各子ノード１〜８に至る移動コストをノードに記録する。
(d) さらに各子ノードから移動可能な８つの矩形を孫ノードとして並べ、各孫ノードには、その親となる子ノードに記録された移動コストに、その親からその子へ移動する際にかかる移動コストを加えた値（合計移動コスト）を記録する。
(e) (d)を再帰的に繰り返して木を広げていく。途中、ノードのコスト（合計移動コスト）がCを越える場合には子ノードの展開を停止する。
【００５０】
移動コストがすべて正の値であれば上記手順はやがて停止し、木が求まる。その木に現れる全ノードがコストC（所定コスト）で到達可能な範囲（移動可能範囲）となる。
【００５１】
たとえば移動コストが「移動時間」を表し、直前の観測時刻t-1と現時刻tとの差がCとすると、上記方法で求めた各ノードの観測値の集合が、直前の観測値（座標）から時間C内で移動可能な範囲となり、この範囲に含まれない観測値（座標）へは移動不可能ということになる。つまり、現時刻tの観測値集合のうちこの移動可能範囲に含まれるものを候補として残す。
【００５２】
要するに、上述のステップＳ３２およびＳ３３では、直前の観測時刻t-1までに特定された第１〜第Ｌの観測値系列に対してそれぞれ個別に、次の観測時刻tにおけるセンサデータの複数の観測値のうち、当該観測値系列における最新の観測値を含む区間（矩形）から隣接区間を順次移動した場合の合計移動コストが所定コスト以下で移動可能な区間に含まれる観測値群を、第１〜第Ｌの観測値（候補）群として選択する。合計移動コストは、隣接区間同士間では一方の区間から他方の区間への移動コストそのものになる。また直前の観測時刻t-1の観測値と、最新の観測値とが同一区間に含まれる場合の合計移動コストは処理の簡単のためたとえば単純にゼロとする（いったん隣接区間へ移動し再度元の区間へ戻る場合も観測エリアの分割サイズによっては想定されるがここではそのようなケースは考えないものとする）。この段落で述べた処理は軌跡算出部１０４が備える選択手段２０１により行われる。
【００５３】
別の言い方をすれば、上述のステップＳ３２およびＳ３３では、直前の観測時刻t-1までに特定された第１〜第Ｌの観測値系列毎に、末尾の観測値を含む区間から、次の観測時刻tにおける各観測値を含む区間のそれぞれに対して、移動コストが所定コスト以下で行くことが可能であるような経路が存在するかどうかを調べ、そのような経路が存在する観測値については候補として残す。
【００５４】
＜ステップＳ６（軌跡生成）の詳細＞
図１０はステップＳ６の計算手順の詳細を示すフローチャートである。
【００５５】
まず変数Aに十分小さな値を代入する（Ｓ４１）。
【００５６】
次に、L本ある既存軌跡θ（t-1,l）(l=1〜L)のそれぞれ毎に、ステップＳ５で得られた観測値（候補）のうちの１つを、同一の観測値（候補）が複数の既存軌跡に含まれないように選択して割り当てる（Ｓ４２）。これが１つの割り当てパターン（候補の組み合わせ）になる。すなわち、第１〜第Ｌの既存軌跡に対して第１〜第Ｌの観測値群からそれぞれ異なる１つの第１〜第Ｌの観測値（候補）を選択して割り当てる割り当てパターンを複数生成する。この処理は軌跡算出手段１０４が備える割り当てパターン生成手段２０２の処理に対応する。
【００５７】
次に、各既存軌跡（各観測値系列）から、上記割り当てパターンに示されるそれぞれの観測値（候補）への尤度を計算し、各計算した尤度の和B（評価値）を求める（Ｓ４３）。つまり、第１〜第Ｌの観測値系列の、上記割り当てパターンに示される第１〜第Ｌの観測値に対する尤度を計算し、各計算した尤度の和を求める。この尤度の計算は、軌跡算出手段１０４が備える尤度計算手段２０３により行われる。尤度の計算方法の詳細は後述する。
【００５８】
尤度和Bが変数Aより大きいか否かを検査し（Ｓ４４）、大きい場合は（ＹＥＳ）、Bの値を変数Aに代入し、この割り当てパターンを記録する（Ｓ４５）。
【００５９】
次に、生成可能なすべての割り当てパターンを生成したか否かを検査し（Ｓ４６）、まだ生成していない割り当てパターンが存在するときは（ＮＯ）ステップＳ４２に戻り、すべての割り当てパターンを生成済みのときは（ＹＥＳ）、ステップＳ４５で最後に記録された割り当てパターンを返す（Ｓ４７）。すなわち尤度和Bが最も大きいときの割り当てパターンを返す。ここでは尤度和Bが最も大きいときの割り当てパターンを返しているが、閾値以上の尤度和をもつ割り当てパターンのうちの１つ選択して返しても良い。つまり、尤度和が最も高いまたは閾値以上である１つの割り当てパターンを、ステップＳ４２で生成した複数の割り当てパターンの中から選択し、選択した割り当てパターンに示される各観測値（第１〜第Ｌの観測値）を第１〜第Ｌの観測値系列の末尾に追加する。この処理は、軌跡算出手段１０４が備える観測値系列生成手段２０４により行われる。
【００６０】
＜ステップＳ４３（尤度計算）の詳細＞
ここでは尤度の計算方法としてカルマンフィルタを用いる例を示す。
【００６１】
まずカルマンフィルタについて説明する。
【００６２】
ターゲット（移動体）の挙動は例えば線形カルマンフィルタで記述できる。
【００６３】
カルマンフィルタではターゲットの状態変数x(t)を
【数３】
と記述する。Fは状態遷移行列であり、Qはシステムノイズである。状態変数は位置座標および速度からなるベクトル変数である。また、状態変数から観測値を得るために
【数４】
を用いる。Hは観測行列、Rは観測ノイズである。
【００６４】
観測されるものがxとy座標、すなわちターゲットの位置のみであるとすると、観測値は
【数５】
となる。a,bはそれぞれターゲットのｘ座標、ｙ座標である。ターゲットの移動軌跡は事前に分からないため、等速度運動または等加速度運動を仮定し、その大きさおよび方向は観測値から推定する。等速度運動を仮定するなら状態変数は、位置(a,b)と、速度
【数６】
とからなる。
【数７】
となる。
【数８】
である。q1,q2,r1,r2は平均０の白色雑音である。
【００６５】
以下では簡単のため一部、時刻tを省略して式を書く。例えばx(tk)をxkと書く。
【００６６】
カルマンフィルタによる処理は、時刻tkにおける観測値z(tk)が得られる前の「状態変数の事前推定値
【数９】
」と観測値が得られた後の「事後推定値
【数１０】
」を交互に求めながら進められる。
【００６７】
［カルマンフィルタの処理手順］
【数１１】
【００６８】
マルチターゲット追跡ではある時刻tにおいて複数の位置にターゲットが観測される。観測時刻をti(i=1,…,n)、このときに得られた観測座標をzj(ti)(j=1,…,mi)とする。簡単のために各時刻で観測されるターゲット数が一定数ｍであるとすると、時刻t1〜tnの間に合計でｍ×n個の観測値が得られることになる。時刻t1における観測値集合から適当に一つ選びzl1(t1)とする。同様に時刻t2における観測値集合から一つ選びzl2(t2)とする。これを時刻tkまで繰り返すと、zl1(t1),zl2(t2),zl3(t3),....,zlk(tk)という軌跡が作られることになる。このような軌跡をL個作ることにし、各軌跡を
【数１２】
と書くことにすると、
【数１３】
は次式となる。
【数１４】
【００６９】
（式２）の尤度の対数に負号をつけたものを
【数１５】
とする。尤度Λは値が大きいほど「尤もらしい（確度が高い）」ことを表しているので、λは値が小さいほど尤もらしいことを表していることになる。λは
【数１６】
と書ける。
【００７０】
つまり、時刻t1から時刻tk-1までの観測値から作られた軌跡
【数１７】
に、新たに時刻tkでの観測値を１つ加えて、軌跡を１つ伸ばした時のλは、伸ばす前のλに、新たに加えた分の観測値から計算される量を足し加えればよいことが分かる。
【００７１】
したがって、新たに得られた観測値集合から観測値を一つ選び出して軌跡に加える際には、（式３）右辺の第２項が最小のものを選ぶのが良いことになる。ただし、マルチターゲット追跡における複数（ｌ＝１〜L）の軌跡にとって（式３）右辺第２項が最小となる観測値がかち合う場合もあり、ターゲットが互いに接近すると必ずしも尤度が最適となる観測値を割り当てることが正解とはならない。したがって、観測値が得られるたびに、複数の割り当てパターンを考慮して観測値の割り当てを決定する必要がある（図１０のＳ４２参照）。それ故、観測値が得られるたびに、全体における観測値の組み合わせの数が指数的に増えることになるが、本実施形態では、この点、移動可能範囲に含まれる観測値群のみを割り当ての可能性のある候補（観測値）群として扱うため（図７のステップＳ５の説明を参照）、結果として、全体における観測値の組み合わせの数を少なく抑えることが可能である。
【００７２】
以上のように、本実施形態によれば、隣接矩形同士間に地形的な要因に基づく移動コストを定義し、ある時刻の観測値から次の時刻までに移動可能な範囲を移動コストに基づき求めて当該移動可能な範囲に含まれる観測値のみを次の時刻の移動先候補として選択することにより、尤度計算において各既存軌跡（観測値系列）に対する割り当てパターンの数（組み合わせパターン数）を減少させることができるとともに、最も良い割り当てパターンを選択できる。
【００７３】
従来であれば、たとえば１時刻前までに得られた既存軌跡がL個（θ１〜θL）あり、現時刻で観測された観測値がL個（ｚ１〜zＬ）あるとすると、１：１の条件だけを満たすL！通りの組み合わせ（割り当てパターン）
[(θ1,z1) (θ2,z2) (θ3,z3),..., (θL,zL)]、
[(θ1,z2) (θ2,z1) (θ3,z3),..., (θL,zL)]、
・・・
[(θ1,zL) (θ2,zL-1) (θ3,zL-2),..., (θL,z1)]
全部について尤度和を計算し、最も良い尤度和をもつ割り当てパターンを選ぶことになる。したがって計算量が膨大になり、また地形（移動コスト）を考慮していないため実際には不可能な移動を含む割り当てパターンが選択されることもあり得る。
【００７４】
これに対し本実施形態では、地理条件（移動コスト）を元に、対応可能なθとｚの組み合わせを制限しているため、制約を満たさぬ組み合わせ（割り当てパターン）を除いた残りについて尤度和の計算を行い、最も良い組み合わせ（割り当てパターン）を選ぶことができる。
【符号の説明】
【００７５】
１０１：センサデータ記憶部
１０２：地理情報記憶部（第１の記憶手段、第２の記憶手段）
１０３：地理制約計算部（第２の記憶手段）
１０４：軌跡算出部（選択手段、割り当てパターン生成手段、尤度計算手段、観測値系列生成手段、読み出し手段）
１０５：軌跡表示部
【特許請求の範囲】
【請求項１】
観測エリアを移動する複数の移動体の移動軌跡を算出する移動軌跡算出装置であって、
前記観測エリアを位置観測センサにより一定時間毎に観測することによって得られた、それぞれ観測時刻と複数の観測値とを含む複数のセンサデータを記憶するセンサデータ記憶手段と、
前記観測エリアを分割した複数の区間における隣接区間同士毎に一方の区間から他方の区間および前記他方の区間から前記一方の区間への移動コストを記憶する第１の記憶手段と、
直前の観測時刻までに特定された第１〜第Ｌ(Ｌは２以上の整数)の観測値系列に対してそれぞれ個別に、次の観測時刻における前記センサデータの前記複数の観測値のうち、前記観測値系列の最新の観測値を含む区間から隣接区間を順次移動した場合の合計移動コストが所定コスト以下で移動可能な区間に含まれる観測値群を第１〜第Ｌの観測値群として選択する選択手段と、
前記第１〜第Ｌの観測値系列に対して前記第１〜第Ｌの観測値群からそれぞれ異なる１つの第１〜第Ｌの観測値を選択して割り当てる複数の割り当てパターンを生成する割り当てパターン生成手段と、
前記複数の割り当てパターンのそれぞれ毎に、前記第１〜第Ｌの観測値系列の、前記第１〜第Ｌの観測値に対する尤度を計算する尤度計算手段と、
前記第１〜第Ｌの観測値系列に対する尤度の合計が最も高いまたは閾値以上である１つの割り当てパターンを前記複数の割り当てパターンの中から選択し、選択した割り当てパターンに示される第１〜第Ｌの観測値を前記第１〜第Ｌの観測値系列の末尾に追加する観測値系列生成手段と、
を有する軌跡算出手段と、
前記第１〜第Ｌの観測値系列を表示する軌跡表示手段と、
を備えた移動軌跡算出装置。
【請求項２】
前記観測エリアを分割した前記複数の区間のそれぞれについて、前記区間の平均標高と、前記区間の各隣接する区間の方向への斜度とを記憶し、
また前記観測エリアを分割した前記複数の区間における隣接区間同士毎に一方の区間から他方の区間および前記他方の区間から前記一方の区間への平地移動コストを記憶する、
第２の記憶手段と、
前記一方の区間である区間Ａから、前記他方の区間である区間Ｂへの移動コストを、前記区間Ａから前記区間Ｂへの平地移動コストに、前記区間Ｂの方向に対する前記区間Ａの斜度にパラメータαを乗じた値を加算することによって計算する地理制約計算手段と、
を備え、
前記αは、
前記区間Ｂの平均標高が前記区間Ａの平均標高より大きくかつその第１の差が第１所定値より大きいとき正の値、
前記区間Ａの平均標高が前記区間Ｂの平均標高より大きくかつその第２の差が第２所定値より大きいとき負の値、
前記第１の差が前記第１所定値以下のとき、または前記第２の差が前記第２所定値以下のときはゼロであり、
前記第１所定値および第２所定値はゼロ以上の任意の実数である、
ことを特徴とする請求項１に記載の移動軌跡算出装置。
【請求項３】
前記区間の分割サイズは、前記位置観測センサの精度、各前記移動体の移動速度、前記観測エリアのサイズ、の少なくとも一つに基づいて決められた
ことを特徴とする請求項１ないし２のいずれか一項に記載の移動軌跡算出装置。
【請求項４】
観測エリアを移動する複数の移動体の移動軌跡を算出する移動軌跡算出方法であって、
読み出し手段が、前記観測エリアを位置観測センサにより一定時間毎に観測することによって得られた、それぞれ観測時刻と複数の観測値とを含む複数のセンサデータを記憶するセンサデータ記憶手段からデータを読み出し、
前記読み出し手段が、前記観測エリアを分割した複数の区間における隣接区間同士毎に一方の区間から他方の区間および前記他方の区間から前記一方の区間への移動コストを記憶する第１の記憶手段からデータを読み出し、
軌跡算出手段が、直前の観測時刻までに特定された第１〜第Ｌの観測値系列に対してそれぞれ個別に、次の観測時刻における前記センサデータの前記複数の観測値のうち、前記観測値系列の最新の観測値を含む区間から隣接区間を順次移動した場合の合計移動コストが所定コスト以下で移動可能な区間に含まれる観測値群を第１〜第Ｌの観測値群として選択し、
前記軌跡算出手段が、前記第１〜第Ｌの観測値系列に対して前記第１〜第Ｌの観測値群からそれぞれ異なる１つの第１〜第Ｌの観測値を選択して割り当てる複数の割り当てパターンを生成し、
前記軌跡算出手段が、前記複数の割り当てパターンのそれぞれ毎に、前記第１〜第Ｌの観測値系列の、前記第１〜第Ｌの観測値に対する尤度を計算し、
前記軌跡算出手段が、前記第１〜第Ｌの観測値系列に対する尤度の合計が最も高いまたは閾値以上である１つの割り当てパターンを前記複数の割り当てパターンの中から選択し、選択した割り当てパターンに示される第１〜第Ｌの観測値を前記第１〜第Ｌの観測値系列の末尾に追加し、
軌跡表示手段が、前記第１〜第Ｌの観測値系列を表示する
ことを特徴とする移動軌跡算出方法。
【請求項１】
観測エリアを移動する複数の移動体の移動軌跡を算出する移動軌跡算出装置であって、
前記観測エリアを位置観測センサにより一定時間毎に観測することによって得られた、それぞれ観測時刻と複数の観測値とを含む複数のセンサデータを記憶するセンサデータ記憶手段と、
前記観測エリアを分割した複数の区間における隣接区間同士毎に一方の区間から他方の区間および前記他方の区間から前記一方の区間への移動コストを記憶する第１の記憶手段と、
直前の観測時刻までに特定された第１〜第Ｌ(Ｌは２以上の整数)の観測値系列に対してそれぞれ個別に、次の観測時刻における前記センサデータの前記複数の観測値のうち、前記観測値系列の最新の観測値を含む区間から隣接区間を順次移動した場合の合計移動コストが所定コスト以下で移動可能な区間に含まれる観測値群を第１〜第Ｌの観測値群として選択する選択手段と、
前記第１〜第Ｌの観測値系列に対して前記第１〜第Ｌの観測値群からそれぞれ異なる１つの第１〜第Ｌの観測値を選択して割り当てる複数の割り当てパターンを生成する割り当てパターン生成手段と、
前記複数の割り当てパターンのそれぞれ毎に、前記第１〜第Ｌの観測値系列の、前記第１〜第Ｌの観測値に対する尤度を計算する尤度計算手段と、
前記第１〜第Ｌの観測値系列に対する尤度の合計が最も高いまたは閾値以上である１つの割り当てパターンを前記複数の割り当てパターンの中から選択し、選択した割り当てパターンに示される第１〜第Ｌの観測値を前記第１〜第Ｌの観測値系列の末尾に追加する観測値系列生成手段と、
を有する軌跡算出手段と、
前記第１〜第Ｌの観測値系列を表示する軌跡表示手段と、
を備えた移動軌跡算出装置。
【請求項２】
前記観測エリアを分割した前記複数の区間のそれぞれについて、前記区間の平均標高と、前記区間の各隣接する区間の方向への斜度とを記憶し、
また前記観測エリアを分割した前記複数の区間における隣接区間同士毎に一方の区間から他方の区間および前記他方の区間から前記一方の区間への平地移動コストを記憶する、
第２の記憶手段と、
前記一方の区間である区間Ａから、前記他方の区間である区間Ｂへの移動コストを、前記区間Ａから前記区間Ｂへの平地移動コストに、前記区間Ｂの方向に対する前記区間Ａの斜度にパラメータαを乗じた値を加算することによって計算する地理制約計算手段と、
を備え、
前記αは、
前記区間Ｂの平均標高が前記区間Ａの平均標高より大きくかつその第１の差が第１所定値より大きいとき正の値、
前記区間Ａの平均標高が前記区間Ｂの平均標高より大きくかつその第２の差が第２所定値より大きいとき負の値、
前記第１の差が前記第１所定値以下のとき、または前記第２の差が前記第２所定値以下のときはゼロであり、
前記第１所定値および第２所定値はゼロ以上の任意の実数である、
ことを特徴とする請求項１に記載の移動軌跡算出装置。
【請求項３】
前記区間の分割サイズは、前記位置観測センサの精度、各前記移動体の移動速度、前記観測エリアのサイズ、の少なくとも一つに基づいて決められた
ことを特徴とする請求項１ないし２のいずれか一項に記載の移動軌跡算出装置。
【請求項４】
観測エリアを移動する複数の移動体の移動軌跡を算出する移動軌跡算出方法であって、
読み出し手段が、前記観測エリアを位置観測センサにより一定時間毎に観測することによって得られた、それぞれ観測時刻と複数の観測値とを含む複数のセンサデータを記憶するセンサデータ記憶手段からデータを読み出し、
前記読み出し手段が、前記観測エリアを分割した複数の区間における隣接区間同士毎に一方の区間から他方の区間および前記他方の区間から前記一方の区間への移動コストを記憶する第１の記憶手段からデータを読み出し、
軌跡算出手段が、直前の観測時刻までに特定された第１〜第Ｌの観測値系列に対してそれぞれ個別に、次の観測時刻における前記センサデータの前記複数の観測値のうち、前記観測値系列の最新の観測値を含む区間から隣接区間を順次移動した場合の合計移動コストが所定コスト以下で移動可能な区間に含まれる観測値群を第１〜第Ｌの観測値群として選択し、
前記軌跡算出手段が、前記第１〜第Ｌの観測値系列に対して前記第１〜第Ｌの観測値群からそれぞれ異なる１つの第１〜第Ｌの観測値を選択して割り当てる複数の割り当てパターンを生成し、
前記軌跡算出手段が、前記複数の割り当てパターンのそれぞれ毎に、前記第１〜第Ｌの観測値系列の、前記第１〜第Ｌの観測値に対する尤度を計算し、
前記軌跡算出手段が、前記第１〜第Ｌの観測値系列に対する尤度の合計が最も高いまたは閾値以上である１つの割り当てパターンを前記複数の割り当てパターンの中から選択し、選択した割り当てパターンに示される第１〜第Ｌの観測値を前記第１〜第Ｌの観測値系列の末尾に追加し、
軌跡表示手段が、前記第１〜第Ｌの観測値系列を表示する
ことを特徴とする移動軌跡算出方法。
【図１】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【図１３】
【図１４】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【図１３】
【図１４】
【公開番号】特開２０１０−２２３８１５（Ｐ２０１０−２２３８１５Ａ）
【公開日】平成２２年１０月７日（２０１０．１０．７）
【国際特許分類】
【出願番号】特願２００９−７２３６１（Ｐ２００９−７２３６１）
【出願日】平成２１年３月２４日（２００９．３．２４）
【出願人】（０００００３０７８）株式会社東芝 (54,554)
【Ｆターム（参考）】
【公開日】平成２２年１０月７日（２０１０．１０．７）
【国際特許分類】
【出願日】平成２１年３月２４日（２００９．３．２４）
【出願人】（０００００３０７８）株式会社東芝 (54,554)
【Ｆターム（参考）】
[ Back to top ]