到来方向推定装置、到来方向推定方法および到来方向推定プログラム
【課題】到来方向の推定を高速かつ高精度に行うこと。
【解決手段】Ｎ個のセンサが受信した目標からの到来信号より、ｖ１〜ｖＮで表されるＮ個のベースバンド信号ベクトルを生成する信号ベクトル生成手段と、１≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１から、（Ｎ−Ｍ）×Ｍ次の行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち、少なくとも１つの行列を生成するＨａｎｋｅｌ行列生成手段と、前記行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、前記行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔによって２つの部分行列Ｒ１およびＲ２に分割し、前記部分行列Ｒ１とＲ２とから前記到来信号の到来方向を推定する推定手段と、を具備する到来方向推定装置、方法およびプログラムである。
【解決手段】Ｎ個のセンサが受信した目標からの到来信号より、ｖ１〜ｖＮで表されるＮ個のベースバンド信号ベクトルを生成する信号ベクトル生成手段と、１≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１から、（Ｎ−Ｍ）×Ｍ次の行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち、少なくとも１つの行列を生成するＨａｎｋｅｌ行列生成手段と、前記行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、前記行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔによって２つの部分行列Ｒ１およびＲ２に分割し、前記部分行列Ｒ１とＲ２とから前記到来信号の到来方向を推定する推定手段と、を具備する到来方向推定装置、方法およびプログラムである。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は到来方向推定装置、到来方向推定方法および到来方向推定プログラムに関し、特に複数のアンテナを用い目標からの到来信号が到来する方向を推定する到来方向推定装置、到来方向推定方法および到来方向推定プログラムに関する。
【背景技術】
【０００２】
複数のセンサを用いた目標（信号送信源または信号反射物体）からの到来信号のＤＯＡ（Direction Of Arrival：到来方向）推定には、ＤＢＦ（デジタルビームフォーミング法）、ＳＳＭ（Sub Space Method：部分空間法）、ＭＬ（Maximum Likelihood Method：最尤推定法）等が用いられている。
【０００３】
特許文献１には、ＤＯＡ推定を高速かつ高精度に行うため、到来信号のベースバンド信号の相関ベクトルを組み合わせた空間平均共役分散行列を用いてＤＯＡを推定する技術が開示されている。
【特許文献１】国際公開第２００６／０６７８６９号パンフレット
【発明の開示】
【発明が解決しようとする課題】
【０００４】
しかしながら、特許文献１の技術においてもＤＯＡ推定の計算速度・ロバスト性には改善の余地がある。本件開示の到来方向推定装置、到来方向推定方法および到来方向推定プログラムは、到来方向の推定を高速かつ高精度に行うことが可能な到来方向推定装置、到来方向推定方法および到来方向推定プログラムを提供することを目的とする。
【課題を解決するための手段】
【０００５】
本到来方向推定装置は、Ｎ個のセンサが受信した目標からの到来信号よりｖ１〜ｖＮで表されるＮ個のベースバンド信号ベクトルを生成する信号ベクトル生成手段と、 １≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１から（Ｎ−Ｍ）×Ｍ次の行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を生成するＨａｎｋｅｌ行列生成手段と、前記行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、前記行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔによって２つの部分行列Ｒ１およびＲ２に分割し、前記部分行列Ｒ１とＲ２とから前記到来信号の到来方向を推定する推定手段と、を具備し、前記行列Ｒｆ１は前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列であり、前記行列Ｒｆ２は前記信号ベクトルの部分ベクトルｖ２〜ｖＮを構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｅｌ行列であり、前記行列Ｒｂ１は、Ｒｆ１＊を前記信号ベクトルの共役複素部分ベクトルｖ１＊〜ｖＮ−１＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ１＝ＪＮ−ＭＲｆ１＊ＪＭであり、前記行列Ｒｆ２＊は前記信号ベクトルの共役複素部分ベクトルｖ２＊〜ｖＮ＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ２＝ＪＮ−ＭＲｆ２＊ＪＭであることを特徴とする。
【０００６】
また、本到来方向推定方法は、Ｎ個のセンサが受信した目標からの到来信号よりｖ１〜ｖＮで表されるＮ個のベースバンド信号ベクトルを生成するステップと、１≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１から（Ｎ−Ｍ）×Ｍ次の行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を生成するステップと、前記行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、前記行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔによって２つの部分行列Ｒ１およびＲ２に分割し、前記部分行列Ｒ１とＲ２とから前記到来信号の到来方向を推定するステップと、を有する方法である。さらに、本到来方向推定プログラムは、コンピュータに上記方法を実行させるプログラムである。
【発明の効果】
【０００７】
本到来方向推定装置、到来方向推定方法および到来方向推定プログラムによれば、到来方向の推定を高速かつ高精度に行うことが可能となる。
【発明を実施するための最良の形態】
【０００８】
まず、特許文献１に係る技術を用いてＤＯＡ推定を行う場合に、計算速度およびロバスト性を改善する余地が残る要因について説明する。なお、以下で、行列やベクトル後に付けた上添え字Ｈは複素共役転置を、Ｔは置換を、＊は各複素共役を取ることを、そして、行列の後に付けた上添え字−１は逆行列を示す。また、行列ＪＰは反対角要素が１でその他の要素が０であるＰ次反対角単位行列を示し、行列ＩｐはＰ次単位行列である。０Ｐ×ＱはＰ×Ｑ次のゼロ行列、０ＰはＰ次正方ゼロ行列を示す。
【０００９】
図１は、本発明の受信用センサ（例としてアンテナ１０とする）を示す図である。図１では、Ｎ個のセンサとして、アンテナＡ１〜ＡＮがＸ軸方向の直線上に等しい間隔ｄを隔てて配置されている。これにより、複数のアンテナ１０は等間隔リニアアレーアンテナ（ＵＬＡ）を構成する。アンテナ１０から距離ｒｍ、角度θｍ（Ｙ軸の正方向を角度の０度とし、時計回りを正方向として測る）なる位置に、独立したＭ個の目標１２（信号送信源または信号反射物体、図１では煩雑さを避ける為、１個で代表させている）が存在しているものとする。
【００１０】
さて、説明を行う便宜上、具体的な装置として波長λなる高周波信号をキャリア信号とするレーダを考える。レーダは、受信用アンテナ１０と空間的に隔てた位置に設置された送信用アンテナから目標に向けて探知信号を送信し、これが目標で反射される事によって生成される、所謂、エコー信号をアンテナ１０で受信し、エコー信号と送信信号とをミキシングしてベースバンド信号に復調し、更に、適切な信号処理を施して各目標の位置や速度情報を推定する。よって図１の状況では、Ｍ（≦Ｎ−１）個の目標からのエコー信号ｘｍ（ｔ）がアンテナ１０に到来する。なお、ｍ＝１，．．．，Ｍであり、ｔは時間である。
【００１１】
このとき、ｎ番目のアンテナＡｎからの入力信号を復調して得られた信号ｖｎ（ｔ）は式１で表される。
【数１】
ただし、
【数２】
である。ここで、ｘｍ（ｔ）はベースバンド信号、ｎｍ（ｔ）は加法性ガウス雑音、λはキャリア信号である。
【００１２】
信号ｖｎ（ｔ）をベースバンド信号ベクトルとして表すと式３となる
【数３】
ただし、
【数４】
【数５】
【数６】
である。ｘ（ｔ）とｎ（ｔ）との間に相関が無いものとして、式３からｖ（ｔ）の共分散行列Ｒｖｖを計算すると、式７となる。
【数７】
ここで、Ｒｘｘはベースバンド信号ベクトルｘ（ｔ）の共分散行列であり、式８で表される。
【数８】
ＲｖｖがＤＯＡ推定する場合の対象となる。
【００１３】
ここではレーダについて考えているので、ＵＬＡで受信された信号は、同じ信号源からの送信された信号が各目標によって反射されてきものであり、コヒーレント性が高い為、Ｎ×Ｎ次行列であるＲｖｖの階数（ランク）は１に縮退する。したがって、このままではＲｖｖからＤＯＡを推定することは難しい。
【００１４】
そこで、Ｒｖｖの主対角線方向にＬ×Ｌ次部分行列（但し、Ｌ＜Ｎ）をＮ−Ｌ＋１個取り、これらの部分行列を足し合わせて平均することにより、Ｒｖｖのランクを回復する、前方空間平均法（forward spatial smoothing：ｆｓｓ）が用いられる。また、ＵＬＡの基準点を反転させｆｓｓと同様の操作を行う後方空間平均法（backward spatial smoothing：ｂｓｓ）も用いられる。通常は双方を組み合わせ、ｆｂｓｓ（forward backward spatial smoothing）として用いられる。
【００１５】
Ｒｖｖにｆｂｓｓを適用して得られるＲｖｖｆｂｓｓと表すことにすれば、例えば、ｃａｐｏｎ法によるＤＯＡ推定には、式９が用いられる。ただし、ａ（θ）はＬ次のモードベクトルであり、式４を構成するベクトル要素ａ（θｍ）と同様の構造を持つ。
【数９】
ｍｕｓｉｃ（Multiple Signal Classification）法では、Ｒｖｖｆｂｓｓを式１０のように固有値分解して行列ＥＮを求め、式１１によってＤＯＡ推定を行う。
【数１０】
【数１１】
すなわち、式９または式１１において、θをパラメータとするモードベクトルa(θ)=[1,exp(j2παsinθ)],…,exp(j2πα(L-1)sinθ)]Tを定義し、θを走査しながら式１１を計算することで行列に含まれる角度情報を調べ、ＤＯＡを推定する。ただし、α＝ｄ／λである。
【００１６】
通常、信号共分散行列には式７の様に雑音成分σ２Ｉが不可避的に含まれるが、その影響を最小限に抑えてＤＯＡの推定精度を向上させる為、特許文献１の技術では以下のような処理を行っている。すなわち、ＵＬＡを構成する各アンテナで受信した信号を復調して得られたベースバンド信号から相関ベクトルを計算し、擬似共分散行列を生成する。次いで、この擬似共分散行列からプロパゲータ行列を生成し、ＤＯＡ推定を行う。以下にこの推定方法について詳しく説明する。
【００１７】
ここでは、簡単のため、ｆｓｓのみについて考える。ベースバンド信号の相関ベクトルｒｖ１を式１２のように計算する。
【数１２】
これ以降、Ｍは、Ｎ個のアンテナで構成されたＵＬＡによって、ＤＯＡを推定する事が可能な最大の目標数を表すものとする。例えば、Ｍは（Ｎ−１）／２以下の最大の自然数である。
【００１８】
次に、ｒｖ１の要素を並べ替えて、Ｒｆ１を式１３のように生成する。
【数１３】
ここで、
【数１４】
であるから、Ｒｆ１はＲｆ１＝ＡＸなる形をしている。すなわち、式１５に示す行列Ｒは、空間平均を適用した共分散行列Ｒｖｖと同じ位相情報を含んでいる。
【数１５】
【００１９】
以下、Ｎ＝５、Ｍ＝２の場合を例に説明する。このとき行列Ｒは式１６となる。
【数１６】
次に、式１７のように、行列ＲをＭ×Ｍ次部分行列Ｒ１および（Ｎ−２Ｍ）×Ｍ次部分行列Ｒ２とに分割する。
【数１７】
部分行列Ｒ１およびＲ２を用いて、Ｍ×（Ｎ−２Ｍ）次行列Γ＝（Ｒ１Ｒ１Ｈ）−１Ｒ１Ｒ２Ｈを生成する。
【００２０】
行列ΓとＮ−２Ｍ次単位行列：ＩＮ−２Ｍから（Ｎ−Ｍ）×（Ｎ−２Ｍ）次プロパゲータ行列Π＝［Γ⊥−ＩＮ−２Ｍ］Ｔを生成する。
【００２１】
さらに、この行列Πを用いて、例えば、（Ｎ−Ｍ）×（Ｎ−Ｍ）次核行列Ωを、Ω＝Π（ΠＨΠ）−１ΠＨによって定義し、これとアレイモードベクトルa(θ)=[1,exp(j2παsinθ)],…,exp(j2πα(N-M-1)sinθ)]Tとを用いて、例えば、角度スペクトラムＰ（θ）を式１８のように定義する。ただし、α＝ｄ／λである。
【数１８】
そして、パラメータθを走査しながら角度スペクトラムＰ（θ）を計算し、そのピーク位置から到来信号の到来角度θｍを得ることができる。または、ｚ＝ｅｘｐ（ｊ２παｓｉｎ（θ））、α＝ｄ／λ、に対し、代数方程式ａ（１／ｚ）ＴΩａ（ｚ）＝０の解から到来信号の到来角度θｍを得ることができる。
【００２２】
以上のように、特許文献１では、式１５および１６のように行列Ｒを求め、行列Ｒからプロパゲータ行列Πを求め、核行列Ωを求め、到来方向を算出している。
【００２３】
説明を簡単にするため、以下式１６に示した様に、Ｎ＝５、Ｍ＝２の場合を例に取る。プロパゲータ行列Πの主要部分であるΓは式１９のようになる。
【数１９】
ＤＯＡ推定を行うためには、Γを構成するベクトルの位相成分があれば十分である。しかしながら、式１９に掛かっている因子１／｜ｖ５｜２は、ＤＯＡ推定に直接関与しないばかりか、計算の高速化を阻害さえする。また実際の計算に於いては、桁落ちや丸め誤差等の原因となり、高精度化を阻害する。
【００２４】
そこで、本到来方向推定装置は、上記因子を含まない行列を用いてＤＯＡ推定を行うことを特徴とする。
【００２５】
１≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、Ｎ−Ｍを列の次元とし、ベースバンド信号ベクトルｖ１〜ｖＮの部分ベクトルｖ１〜ｖＮ−１を構成要素として並べ、式２０のようなＨａｎｋｅｌ行列Ｒｆ１を生成する。
【数２０】
【００２６】
具体例として、Ｎ＝５およびＭ＝２の場合についてＲｆ１を式２１に示す。
【数２１】
【００２７】
次に、例えば、式２２に示す様に行列Ｒを定義し、このＲをＭ×Ｍ次部分行列Ｒ１および（Ｎ−２Ｍ）×Ｍ次部分行列Ｒ２に分割する。
【数２２】
【００２８】
上記部分行列Ｒ１およびＲ２から、Ｍ×（Ｎ−２Ｍ）次行列Γ＝（Ｒ１Ｒ１Ｈ）−１Ｒ１Ｒ２Ｈを計算すると、式２３のようになる。
【数２３】
式２３の行列には、ＤＯＡ推定に用いられる位相情報が保存されている反面、式１９に含まれていた余分な因子：１／｜ｖ５｜２は含まれない。
【００２９】
以下、式１７と同様に、行列ΓとＮ−２Ｍ次元単位行列ＩＮ−２Ｍとからプロパゲータ行列Π＝［Γ⊥−ＩＮ−２Ｍ］Ｔを生成する。
【００３０】
更に、この行列Πを用いて、例えば、（Ｎ−Ｍ）×（Ｎ−Ｍ）次核行列Ωを、Ω＝Π（ΠＨΠ）−１ΠＨによって定義し、これとアレイモードベクトルa(θ)=[1,exp(j2παsinθ)],…,exp(j2πα(N-M-1)sinθ)]Tとを用いて、例えば、角度スペクトラムＰ（θ）を式１８のように定義する。但し、α＝ｄ／λである。
【数２４】
そして、パラメータθを走査しながら角度スペクトラムＰ（θ）を計算し、そのピーク位置から到来信号の到来角度θｍを得ることができる。または、ｚ＝ｅｘｐ（ｊ２παｓｉｎ（θ））、α＝ｄ／λ、として、代数方程式ａ（１／ｚ）ＴΩａ（ｚ）＝０の解から到来信号の到来角度θｍを得ることができる。
【００３１】
本到来方向推定装置によれば、行列Ｒｆ１として、式１３の代わりに式２１のような行列Ｒｆ１を用いる。これにより式２３のΓは、式１９の１／｜ｖ５｜２のような因子を含まない。よって、Γの計算、プロパゲータ行列Πの計算を高速高精度に行うことができる。以下、上記原理に基づく本発明の実施例について説明する。
【実施例１】
【００３２】
図２は実施例１に係る到来方向推定装置のブロック図である。図３は到来方向推定装置が行う処理のフローチャートである。図２および図３を参照して説明を行う。図２に於いて、到来方向推定装置１００は、複数のアンテナ１０（アンテナＡ１〜ＡＮ）、受信部２０および到来方向推定部３０を含んでいる。複数のアンテナ１０は目標からの到来信号を受信する。受信部２０は受信器２２および信号ベクトル生成手段２４を含んでいる。受信器２２は、例えば、キャリア信号の周波数が７６ＧＨｚである到来信号を受信し、変調方式に応じた復調を行い、適切な帯域のベースバンド信号にダウンコンバートする。そして更に、このベースバンド信号をＡ／Ｄ変換してデジタル信号に変換する（以下、単にベースバンド信号と記す）。信号ベクトル生成手段２４は、各アンテナに対応するベースバンド信号をベクトル列として整形する。以上により、信号ベクトル生成手段２４は式３のベースバンド信号ベクトルｖ１〜ｖＮを生成する（図３のステップＳ１０）。
【００３３】
到来方向推定部３０は、到来信号数予測／設定手段３１、Ｈａｎｋｅｌ行列生成手段３２および推定手段４５を含んでおり、ベースバンド信号ベクトルからＤＯＡを推定する。推定手段４５は、プロパゲータ行列生成手段３３、核行列生成手段３４および到来方向算出手段３５を含んでいる。到来方向推定部３０はコンピュータのＣＰＵ（中央演算処理装置）等の演算装置、演算装置に実行させるプログラム、或いはその記憶装置上の実体である。また、受信部２０と到来方向推定部３０との全てが、ソフトウェア無線装置と同様の意味で、演算装置から構成されていても良い。
【００３４】
到来信号数予測／設定手段３１は、到来信号数Ｍを予測、若しくは設定して出力する（図３のステップＳ１２）。到来信号数Ｍを予測する場合は、最尤法に基づくＡＩＣ（Akaike Information Criteria）、ＭＤＬ（Minimum Description Length）、若しくは他の適切な指標を用いることができる。到来信号数予測／設定手段３１は、到来信号数Ｍの予測を行わない場合は、アンテナ数Ｎから自動的に決定される値である［（Ｎ−１）／２］をＭの値としてもよい。なお、この数式において記号［ｘ］は、実数ｘに対して、ｘ以下の最大の自然数を戻す事を意味する。
【００３５】
推定手段４５は、先ず、式３の信号ベクトルのＳＮ比が事前に決定しておいた或る閾値以上かを判定し（ステップＳ１４）、ＳＮ比が閾値より小さい場合には、例えば、相関ベクトルｃｏｎｊ（ｖＮ）×［ｖ１，…，ｖＮ］Ｔを算出し、これを信号ベクトルとして置き直す（ステップＳ１６）。ここでｃｏｎｊ（ｖＮ）は、信号ｖＮの共役複素数を意味する。
【００３６】
Ｈａｎｋａｌ行列生成手段３２は、式２０のようなＨａｎｋｅｌ行列Ｒｆ１を生成する（ステップＳ１６に於いて信号ベクトルを置き換えた場合は、式１６の様な構造となるが、今後は特に触れない）。プロパゲータ行列生成手段３３では、行列Ｒ＝［Ｒｆ１］とする（ステップＳ１８）。
【００３７】
プロパゲータ行列生成手段３３は、行列Ｒを２つの部分行列Ｒ１およびＲ２に分解し、前述のように行列Γ＝（Ｒ１Ｒ１Ｈ）−１Ｒ１Ｒ２Ｈを生成し、さらに行列Γと単位行列ＩＮ−２Ｍとからプロパゲータ行列Π＝［Γ⊥−ＩＮ−２Ｍ］Ｔを生成する（ステップＳ２２）。
【００３８】
核行列生成手段３４は、核行列Ω＝Π（ΠＨΠ）−１ΠＨを生成する（ステップＳ２４）。
【００３９】
到来方向算出手段３５は、パラメータθを含むアレイモードベクトルa(θ)=[1,exp(j2παsinθ)],…,exp(j2πα(N-M-1)sinθ)]Tと核行列Ωとを用いて、式２４のような角度スペクトラムＰ（θ）を生成する（ステップＳ２６）。但し、α＝ｄ／λである。到来方向算出手段３５は、θを走査しながら角度スペクトラムＰ（θ）を計算し、角度スペクトラムＰ（θ）のピーク位置から到来信号の到来角度θｍを計算する。または、ｚ＝ｅｘｐ (ｊ２παｓｉｎ（θ）)、α＝ｄ／λ、として、代数方程式ａ（１／ｚ）ＴΩａ（ｚ）＝０の解から到来信号の到来角度θｍを計算する。到来方向算出手段３５は、到来角度θｍを推定した到来方向として出力する（ステップＳ２８）。
【００４０】
Ｈａｎｋａｌ行列生成手段３２は、信号ベクトルの他の部分ベクトルｖ２〜ｖＮを構成要素とし、式２５のようなＨａｎｋｅｌ行列Ｒｆ２をＲｆ１の代わりに生成し、プロパゲータ行列生成手段３３は、このＲｆ２を用いてプロパゲータ行列Πを生成してもよい。
【数２５】
【００４１】
また、Ｈａｎｋｅｌ行列生成手段３２は、共役複素部分ベクトルｖ１＊〜ｖＮ−１＊を構成要素とするＨａｎｋｅｌ行列、Ｒｂ１＝ＪＮ−ＭＲ＊ｆ１ＪＭをＲｆ１の代わりに生成し、プロパゲータ行列生成手段３３は、このＲｂ１を用いてプロパゲータ行列Πを生成してもよい。なお、Ｒ＊ｆ１は、Ｒｆ１の各要素について共役複素数をとった行列である。
【００４２】
一例として、Ｎ＝５かつＭ＝２の場合を説明すれば、Ｒｂ１は式２６のように表される。
【数２６】
次に、行列Ｒを式２７のように定義し、更に、行列ＲをＭ×Ｍ次部分行列Ｒ１および（Ｎ−２Ｍ）×Ｍ次部分行列Ｒ２に分割する。
【数２７】
【００４３】
上記部分行列Ｒ１およびＲ２から行列Γ＝（Ｒ１Ｒ１Ｈ）−１Ｒ１Ｒ２Ｈを計算すると、式２８のようになる。
【数２８】
式２８においても式２３と同様に、ＤＯＡ推定に用いられる位相情報が保存され、式１９の１／｜ｖ５｜２のような因子は含まない。よって、Γの計算、プロパゲータ行列Πの計算を高速高精度に行うことができる。
【００４４】
さらに、Ｈａｎｋａｌ行列生成手段３２は、共役複素部分ベクトルｖ２＊〜ｖＮ＊を用いてＨａｎｋｅｌ行列Ｒｂ２＝ＪＮ−ＭＲ＊ｆ２ＪＭをＲｆ１の代わりに生成し、プロパゲータ行列生成手段３３は、このＲｂ２を用いてプロパゲータ行列Πを生成してもよい。なお、Ｒ＊ｆ２は、Ｒｆ２の各要素について共役複素数をとった行列である。
【００４５】
プロパゲータ行列生成手段３３は、（Ｎ−Ｍ）×Ｍ次のＨａｎｋｅｌ行列Ｒｆ１、Ｒｆ２、Ｒｂ１、またはＲｂ２をＫ個組み合わせ、（Ｎ−Ｍ）×（Ｋ×Ｍ）の行列として行列Ｒを構成してもよい。例えば、Ｋ＝２の場合Ｒ＝［Ｒｆ１⊥Ｒｆ２］、Ｒ＝［Ｒｂ１⊥Ｒｂ２］、Ｒ＝［Ｒｆ１⊥Ｒｂ２］、またはＲ＝［Ｒｆ２⊥Ｒｆ２］、Ｋ＝４の場合Ｒ＝［Ｒｆ１⊥Ｒｆ２｜Ｒｂ１⊥Ｒｂ２］とすることができる。この場合、Ｋの値に応じて行列ＲをＭ×（Ｋ×Ｍ）の部分行列Ｒ１および（Ｎ−２Ｍ）×（Ｋ×Ｍ）の部分行列Ｒ２に分割し、Ｒ＝［Ｒ１⊥Ｒ２］Ｔとする。
【００４６】
以上のように、プロパゲータ行列生成手段３３は、行列Ｒを行列Ｒｆ１、Ｒｆ２、Ｒｂ１、およびＲｂ２の少なくとも１つを用いて生成し、この行列Ｒを適切な次元の２つの部分行列Ｒ１、Ｒ２に分割する。
【００４７】
実施例１によれば、到来方向推定装置１００は、Ｈａｎｋｅｌ行列生成手段３２がＲｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を生成する。推定手段４５が、行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔなる２つの部分行列Ｒ１およびＲ２に分割し、部分行列Ｒ１とＲ２とを用いて前記行列Γ、Π、そしてΩを計算し、到来信号の到来方向を推定する。これにより、式１９の１／｜ｖ５｜２のような因子を含まないΓ（式２３）を用いる事で、Γの計算、プロパゲータ行列Πの計算を高速高精度に行うことができる。
【実施例２】
【００４８】
図４は、実施例２に係る到来方向推定装置のブロック図である。図２を参照して説明を行えば、到来方向推定装置１００ａは、実施例１の図２と比較し、到来方向推定部３０ａがスケーリング行列生成手段３６を含んでいる。同図の他の構成は、実施例１の図２と同じであるから説明を省略する。
【００４９】
スケーリング行列生成手段３６は、プロパゲータ行列生成手段３３が生成した行列Ｒの部分行列Ｒ１およびＲ２を用い、（Ｎ−２Ｍ）×（Ｎ−２Ｍ）次スケーリング行列Λ＝Ｒ２Ｒ２Ｈ−Ｒ２Ｒ１ＨΓを生成する。核行列生成手段３４は、プロパゲータ行列Πとスケーリング行列Λとから核行列Ω＝Π（Λ）−１ΠＨを生成する。または、Π´＝Π（ΠＨΠ）−１／２とし、核行列をΩ＝Π´（Λ）−１（Π´）Ｈによって生成してもよい。
【００５０】
核行列生成手段３４は、Π（Λ）−１ΠＨの代わりにΠΠＨを核行列Ωとして生成してもよい。この場合、角度スペクトラムは式２９となり、
【数２９】
代数方程式は式３０となる。
【数３０】
ただし、ｚ＝ｅｘｐ(ｊ２παｓｉｎ（θ）)、α＝ｄ／λ、に対し、
【数３１】
である。
【００５１】
核行列生成手段３４は、核行列をΩ＝Π［α（ΠＨΠ）＋（１−α）Λ］−１ΠＨとして生成してもよい。ここで、０≦α≦１である。このように、核行列生成手段３４は、Ｈａｎｋｅｌ行列から導かれる複数の行列を複合させたものとして、核行列を生成してもよい。
【００５２】
ここで、α≠０、１とし、ΠＨΠとΛとが特異点でないとして、Ａ＝ΠＨΠ、Ｂ＝Λ、１−α＝βとすると、式３２となる。
【数３２】
以上より上の核行列Ωは、公知の核行列Π（ΠＨΠ）−１ΠＨおよびΠ（Λ）−１ΠＨの加重平均となる。よって、公知の核行列を単独で用いたのでは推定困難な目標についても角度推定が可能となる。
【００５３】
核行列生成手段３４は、核行列Ωの生成を、実施例１ではプロパゲータ行列Πから行い、また、実施例２ではプロパゲータ行列Πとスケーリング行列Λとから行う。
【００５４】
到来方向算出手段３５は、プロパゲータ行列Πと、Ｄ、（ＢＨＡ−１Ｂ）、（ＢＨＢ）、（Ｄ−ＢＨＢ）、Λ（ΠＨΠ）等、以下に定める行列Ａ，Ｂ，Ｃ、Ｄを適当に組み合わせた代用行列をスケーリング行列Λの代わりに用いて核行列Ωを生成し、上記実施例１と同様にして、角度スペクトラムまたは代数方程式から到来方向を算出してもよい。なお、Ａ＝Ｒ１Ｒ１Ｈ、Ｂ＝Ｒ１Ｒ２Ｈ、Ｃ＝Ｒ２Ｒ１Ｈ、Ｄ＝Ｒ２Ｒ２Ｈである。
【００５５】
例としてＤを用いた角度スペクトラムを式３３、代数方程式を式３４に示す。
【数３３】
【数３４】
ただし、ｚ＝ｅｘｐ(ｊ２παｓｉｎ（θ）)、α＝ｄ／λ、に対し、
【数３５】
である。
【００５６】
さらに、到来方向算出手段３５は、スケーリング行列Λ、その代用行列、またはスケーリング行列Λのノルムを用いて到来方向を算出してもよい。
【００５７】
図５から８は、実施例２をＦＭＣＷ（Frequency Modulated Continuous Wave）レーダに適用し、目標定位を実施した例である（これらの図では、距離の単位はｂｉｎである）。
【００５８】
図５（ａ）から（ｄ）は、目標が２つ存在する場合の角度スペクトラムを角度と距離に対し示した図である。設定した目標の距離は４０ｍ、角度は０°および３°である。図５（ａ）はＦＦＴ−ＤＢＦ法、図５（ｂ）はＦＢＳＳ−ＭＵＳＩＣ法、図５（ｃ）は特許文献１の方法（比較例）、図５（ｄ）は実施例２の方法を用いて計算した角度スペクトラムである。角度スペクトラムが仮定した２つの目標位置に集約している場合、つまり、角度スペクトラムが２つの目標付近で点状に表示されている程、高い精度で目標の到来方向が推定されていることを示している。
【００５９】
図６は、図５（ａ）、（ｃ）および（ｄ）の角度スペクトラムを距離４０ｍの位置で角度軸に平行に切り出したものを重ね描きしたものである。図５および図６より、比較例（特許文献１）および実施例２の方法は、ＦＦＴ−ＤＢＦ法、ＦＢＳＳ−ＭＵＳＩＣ法と比較して高い精度で目標の角度推定を実現している事が分かる。
【００６０】
図７および図８は、空間的に近接した位置に配置された複数の目標に対する定位性能を、５個の目標を仮定して（図中×の点）角度スペクトラムによって評価した例である。モノクロ表示で若干判別しづらいが、等高線の込み合っている箇所が角度スペクトラムのピークである。図７は特許文献１の方法を用いて計算した角度スペクトラムを示し、図８は、実施例２において、核行列をΩ＝Π［α（Π）＋（１−α）Λ］−１ΠＨ（α＝０．５）として計算した角度スペクトラムを示している。図８は図７と比較して、×印の周辺の等高線がより込み合っており、仮定した目標付近に角度スペクトラムのピークがある。このように、核行列として複合された核行列を用いることにより、エコー信号のコヒーレンスが高く、推定が難しい多くの目標を正しく検出することができる。
【実施例３】
【００６１】
実施例３は、自動車等の車両に搭載された到来方向推定装置の例である。実施例３に係る到来方向推定装置１００ｂは、到来方向推定装置１００ｂを含む、例えば、ミリ波レーダを搭載した自車に対して脅威のある目標の角度（到来方向）を出力する。目標としては、例えば自車の前を走行する自動車等の車両である。図９は、実施例３に係る到来方向推定装置１００ｂのブロック図であるが、到来方向推定装置１００ｂは、受信部２０と到来方向推定部３０の間に脅威度評価手段４０を含んでいる。脅威度評価手段４０は、距離／速度／粗到来方向推定手段３７、目標選別手段３８およびビーム形成手段３９を含んでいる。
【００６２】
図１０は、実施例３に係る到来方向推定装置１００ｂの処理を示すフローチャートである。以下、図１０を参照しながら説明を行う。例えば、ミリ波信号にＦＭ変調を掛ける為に用いられる信号の１周期分に渡り、受信アンテナ毎にエコー信号のサンプリングを行い、信号ベクトル生成手段２４が対応する信号ベクトルのデータセットを抽出する（ステップＳ３０）。即ち、アンテナ数をＮ、サンプリング数をＱとすると、このデータはＮ×Ｑ次行列（データセット）となる。距離／速度／粗到来方向推定手段３７は、時間領域に沿ってＱ次元ベクトルがＮ組並んだデータセットをそれぞれのアンテナ毎にＦＦＴ（高速フーリエ変換）処理する（ステップＳ３２） − Ｑ次元ベクトルのＦＦＴをＮ回行い、周波数領域のデータセットに変換されたＮ×Ｑ行列を得る。距離／速度／粗到来方向推定手段３７は、ＦＦＴの結果として周波数領域に現れたピークに対応するｂｉｎ位置（周波数インデックスの値と考えても良い）から、各目標と自車との相対距離ｒｍおよび相対速度ｓｍを推定し、更に、各目標（ｂｉｎ）位置に於けるＮ次元データベクトルを取り出して、例えばＦＦＴ−ＤＢＦ（デジタルビームフォーミング）処理を行って各目標の大雑把な角度位置を把握する（ステップＳ３４）。また、推定した相対距離、相対速度、角度位置を、各目標固有の識別情報としてペアリングする。
【００６３】
目標選別手段３８は、各目標の脅威度を、目標の相対距離ｒｍと相対速度ｓｍとに基づいて評価する（ステップＳ３６）。脅威度はｒｍが小さい程、そして、ｓｍが大きな負の値となる程（即ち、近接した位置に存在し、自車に高速で接近してくる目標程）、大きくなる。脅威度は、例えば式３６の評価関数を用いる。
【数３６】
ただし、
【数３７】
ここで、ｓＲＥＦ、ｒＲＥＦは、車両メーカーの安全基準等に従ってあらかじめ設定された所定値であるが、例えばＤＳＲＣ（Dedicated Short Range Communications）等のインフラ装置からの情報を受け、交通状況に応じて動的に変更する様にしても良い。
【００６４】
目標選別手段３８は、最大の脅威度を有する目標の脅威度が、基準値以上か否かを判断する（ステップＳ３８）。判定がＹｅｓの場合、目標選別手段３８は、先に粗到来方向推定手段に於いて推定した目標の角度を参照して、最大脅威度を有する目標が自車レーン等、自車の動線と交差する位置を走行するか否かを判断する（ステップＳ４０）。
【００６５】
ステップＳ３８またはＳ４０で判定がＮｏとなった場合、脅威度推定手段４０は、周波数軸上で定位される各目標の相対距離ｒｍおよび相対速度ｓｍのデータ（ｂｉｎ）を、後段の到来方向推定部３０に引き渡す。到来方向推定部３０は、指定された各ｂｉｎ位置のＮ次元データベクトルに対して、例えばデジタルビームフォーミング等、精度は比較的低いが、計算量の小さい通常の到来方向推定処理を適用して各目標の角度推定を行う（ステップＳ４２）。もちろん、粗到来方向推定手段による推定値を用いても、実施例１または実施例２の方法を適用しても良い。到来方向推定部３０は、各目標の識別情報である、相対距離、相対速度、および角度推定値のそれぞれを１つ前の測定期間のデータとして更新保存する（ステップＳ４４）。更新されたデータは、次回の測定期間に於いて得られたデータに基づいて、脅威度評価手段４０が各目標の脅威度を再評価する際に参照される。
【００６６】
ステップＳ４０において判定がＹｅｓとなった場合、目標選別手段３８は、最大の脅威度を持つと判定された目標の相対距離ｒｍおよび相対速度ｓｍのデータ（ｂｉｎ）に基づいて、対応するＮ次元データベクトルを抽出してビーム形成手段３９に引き渡す。ビーム形成手段３９は、距離／速度／粗到来方向推定手段３威度を持つ目標の角度情報を参照し、当該目標が存在する方向に、ビームスペース法によるビーム形成を行う（ステップＳ４６）。この様にして、脅威度推定手段４０は、脅威度の高い目標に関するデータベクトルをビームスペースデータに変換し、後段の到来方向推定部３０に引き渡す。到来方向推定部３０は、詳細な角度推定が可能な実施例１または実施例２の方法で到来方向の推定を行う（ステップＳ４８）。装置が、新たにデータ測定を開始した場合、ステップＳ３０に戻って同様な処理を繰り返す。
【００６７】
実施例３の脅威度評価手段４０は、信号ベクトルを用い、目標の自車に対する脅威を示す脅威度を評価する。到来方向推定部３０の推定手段４５は脅威度に基づき、目標の到来方向を推定する、または推定しない。このように、実施例３においては、脅威度評価手段４０を導入する事により、自車の安全に対する目標の重要性を脅威度という指標で評価し、各目標の脅威度に応じて到来方向推定部３０と連携して適切な到来方向推定手法を適用する事で、装置としてのレスポンスタイムを最適化できる。
【実施例４】
【００６８】
実施例４は、信号品質や目標間の相対離角といった追加指標を導入する事で、より、現実的な運用環境に適した処理を実現する装置である。図１１は、実施例４に係る到来方向推定装置１００ｃのブロック図である。同図に於いて、到来方向推定装置１００ｃは、受信部２０と到来方向推定部３０との間に脅威度評価手段５０を含んでおり、ハードウェアとしての基本的な構成は実施例３と殆ど同じであるから説明を省略する。
【００６９】
図１２は、到来方向処理装置１００ｃの処理を示すフローチャートである。脅威度評価手段５０は、信号ベクトルの抽出（ステップＳ５０）から目標毎に相対距離と相対速度の推定（ステップＳ５３）を行うまで、実施例３のそれと同じ処理を行う。但し、実施例３で行っていたＦＦＴ−ＤＢＦによる粗到来方向推定処理は実施しなくても良い。
【００７０】
脅威度評価手段５０の機能として新たに加わったステップは、脅威度の評価（ステップＳ５６）を実施する前の、信号品質の評価（ステップＳ５４）、一つ前の測定期間に推定された目標緒元データの参照（ステップＳ６６）、そして、ステップＳ６６の参照値の相対離角を基準値との比較（ステップＳ６８）、である。
【００７１】
新たに加わった処理の流れは、先ず、ステップＳ５４で、信号ベクトルのＳＮ比が基準値以上か否かを判断し、判定がＹｅｓの場合、ステップＳ５６に進む。一方、判定がＮｏの場合、ステップＳ６６に進み、一つ前の測定期間に推定された目標緒元データのうち到来方向のデータを参照し、ステップＳ６８に於いて参照値の相対離角と基準値との比較を行う。そして、脅威度評価手段５０は、ステップＳ５４の判定がＹｅｓの場合、または、ステップＳ６８の判定がＮｏの場合、各目標について脅威度の評価を行い（ステップＳ５６）、最大の脅威度を有する目標の脅威度が、基準値以上か否かを判断する（ステップＳ５８）。脅威度の評価は、実施例３のステップＳ３６と同じ方法で行うことができる。
【００７２】
脅威度評価手段５０は、ステップＳ５８の判定がＹｅｓである場合、制御を到来方向推定部３０に渡し、詳細な到来方向の推定を行う（ステップＳ６０）。一方、ステップＳ５８の判定がＮｏ、若しくは、ステップＳ６８の判定がＹｅｓである場合、脅威度評価手段５０は、制御を到来方向推定部３０に渡し、粗い到来方向推定を行う（ステップＳ７０）。このように、脅威度評価手段５０は、信号ベクトルを用い目標の自車に対する脅威を示す脅威度を評価する。到来方向推定部３０の推定手段４５は脅威度に基づき、目標の到来方向の推定を粗い角度分解能で行う、または細かい角度分解能で行う。
【００７３】
以下に、ステップＳ６０で行う処理とステップＳ７０で行う処理との違いを示す。行列Ｒの部分行列Ｒ１、Ｒ２を用いて、Ｒ１Ｒ１Ｈ＝Ａ、Ｒ１Ｒ２Ｈ＝Ｂ、Ｒ２Ｒ１Ｈ＝ＢＨ、Ｄ＝Ｒ２Ｒ２Ｈ＝Ｄとし、ＲＲＨに逆行列の定理を適用すると、式３８となり、
【数３８】
これから直ちに、式３９が得られる。
【数３９】
【００７４】
ステップＳ６０では、核行列としてΩ＝（ＲＲＨ）−１−［（Ｒ１Ｒ１Ｈ）−１⊥０Ｍ×Ｎ−２Ｍ；０Ｎ−２Ｍ×Ｍ⊥０Ｎ−２Ｍ］を用い、また、ステップＳ７０では、核行列としてΩ＝（ＲＲＨ）−１を用いて到来方向の推定を行う。ステップＳ７０の行列は信号部分空間成分の成分を含んでいる為、推定精度が劣る。しかしながら、ステップＳ７０の行列を計算する途中で、（Ｒ１Ｒ１Ｈ）−１の計算をしておく事は比較的容易であるから、実施例４の方法は推定精度の切り替えが簡単に実現できる。
【００７５】
以上、実施例４によれば、エコー信号の品質（ＳＮ比）が基準値以上であり、かつ目標の脅威度が基準値以上である場合、ステップＳ６０のように到来方向の詳細な推定を行う。一方、目標の脅威度が小さく詳細な位置の把握が必要が無い場合や、エコー信号のＳＮ比が小さく時間を掛けても精度の良い推定が望めない場合、ステップＳ７０のように到来方向の粗い推定を行う。即ち、実施例４は運用環境に応じてより適切な動作 − レスポンス重視、または精度重視 − を選択できる装置である。
【００７６】
以上、本実施例について詳述したが、本発明は係る特定の実施例に限定されるものではなく、特許請求の範囲に記載された本発明の要旨の範囲内において、種々の変形・変更が可能である。
【図面の簡単な説明】
【００７７】
【図１】図１は受信用センサを示す図である。
【図２】図２は実施例１に係る到来方向推定装置のブロック図である。
【図３】図３は実施例１に係る到来方向推定装置のフローチャートである。
【図４】図４は実施例２に係る到来方向推定装置のブロック図である。
【図５】図５（ａ）から図５（ｄ）は各方法で計算された角度スペクトラムを示す図である。
【図６】図６は距離が一定のときの角度に対する角度スペクトラムを示した図である。
【図７】図７は実施例２と比較例の角度スペクトラムを示した図である。
【図８】図８は実施例２の角度スペクトラムを示した図である。
【図９】図９は実施例３に係る到来方向推定装置のブロック図である。
【図１０】図１０は実施例３に係る到来方向推定装置のフローチャートである。
【図１１】図１１は実施例４に係る到来方向推定装置のブロック図である。
【図１２】図１２は実施例４に係る到来方向推定装置のフローチャートである。
【符号の説明】
【００７８】
１０ アンテナ
１２ 目標
２０ 受信部
３０ 到来方向推定部
３１ 到来信号予想／設定手段
３２ Ｈａｎｋｅｌ行列生成手段
３３ プロパゲータ行列生成手段
３４ 核行列生成手段
３５ 到来方向算出手段
３６ スケーリング行列生成手段
３７ 距離／速度／到来方向推定手段
３８ 目標選別手段
３９ ビーム形成手段
４５ 推定手段
４０、５０ 脅威度評価手段
【技術分野】
【０００１】
本発明は到来方向推定装置、到来方向推定方法および到来方向推定プログラムに関し、特に複数のアンテナを用い目標からの到来信号が到来する方向を推定する到来方向推定装置、到来方向推定方法および到来方向推定プログラムに関する。
【背景技術】
【０００２】
複数のセンサを用いた目標（信号送信源または信号反射物体）からの到来信号のＤＯＡ（Direction Of Arrival：到来方向）推定には、ＤＢＦ（デジタルビームフォーミング法）、ＳＳＭ（Sub Space Method：部分空間法）、ＭＬ（Maximum Likelihood Method：最尤推定法）等が用いられている。
【０００３】
特許文献１には、ＤＯＡ推定を高速かつ高精度に行うため、到来信号のベースバンド信号の相関ベクトルを組み合わせた空間平均共役分散行列を用いてＤＯＡを推定する技術が開示されている。
【特許文献１】国際公開第２００６／０６７８６９号パンフレット
【発明の開示】
【発明が解決しようとする課題】
【０００４】
しかしながら、特許文献１の技術においてもＤＯＡ推定の計算速度・ロバスト性には改善の余地がある。本件開示の到来方向推定装置、到来方向推定方法および到来方向推定プログラムは、到来方向の推定を高速かつ高精度に行うことが可能な到来方向推定装置、到来方向推定方法および到来方向推定プログラムを提供することを目的とする。
【課題を解決するための手段】
【０００５】
本到来方向推定装置は、Ｎ個のセンサが受信した目標からの到来信号よりｖ１〜ｖＮで表されるＮ個のベースバンド信号ベクトルを生成する信号ベクトル生成手段と、 １≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１から（Ｎ−Ｍ）×Ｍ次の行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を生成するＨａｎｋｅｌ行列生成手段と、前記行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、前記行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔによって２つの部分行列Ｒ１およびＲ２に分割し、前記部分行列Ｒ１とＲ２とから前記到来信号の到来方向を推定する推定手段と、を具備し、前記行列Ｒｆ１は前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列であり、前記行列Ｒｆ２は前記信号ベクトルの部分ベクトルｖ２〜ｖＮを構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｅｌ行列であり、前記行列Ｒｂ１は、Ｒｆ１＊を前記信号ベクトルの共役複素部分ベクトルｖ１＊〜ｖＮ−１＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ１＝ＪＮ−ＭＲｆ１＊ＪＭであり、前記行列Ｒｆ２＊は前記信号ベクトルの共役複素部分ベクトルｖ２＊〜ｖＮ＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ２＝ＪＮ−ＭＲｆ２＊ＪＭであることを特徴とする。
【０００６】
また、本到来方向推定方法は、Ｎ個のセンサが受信した目標からの到来信号よりｖ１〜ｖＮで表されるＮ個のベースバンド信号ベクトルを生成するステップと、１≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１から（Ｎ−Ｍ）×Ｍ次の行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を生成するステップと、前記行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、前記行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔによって２つの部分行列Ｒ１およびＲ２に分割し、前記部分行列Ｒ１とＲ２とから前記到来信号の到来方向を推定するステップと、を有する方法である。さらに、本到来方向推定プログラムは、コンピュータに上記方法を実行させるプログラムである。
【発明の効果】
【０００７】
本到来方向推定装置、到来方向推定方法および到来方向推定プログラムによれば、到来方向の推定を高速かつ高精度に行うことが可能となる。
【発明を実施するための最良の形態】
【０００８】
まず、特許文献１に係る技術を用いてＤＯＡ推定を行う場合に、計算速度およびロバスト性を改善する余地が残る要因について説明する。なお、以下で、行列やベクトル後に付けた上添え字Ｈは複素共役転置を、Ｔは置換を、＊は各複素共役を取ることを、そして、行列の後に付けた上添え字−１は逆行列を示す。また、行列ＪＰは反対角要素が１でその他の要素が０であるＰ次反対角単位行列を示し、行列ＩｐはＰ次単位行列である。０Ｐ×ＱはＰ×Ｑ次のゼロ行列、０ＰはＰ次正方ゼロ行列を示す。
【０００９】
図１は、本発明の受信用センサ（例としてアンテナ１０とする）を示す図である。図１では、Ｎ個のセンサとして、アンテナＡ１〜ＡＮがＸ軸方向の直線上に等しい間隔ｄを隔てて配置されている。これにより、複数のアンテナ１０は等間隔リニアアレーアンテナ（ＵＬＡ）を構成する。アンテナ１０から距離ｒｍ、角度θｍ（Ｙ軸の正方向を角度の０度とし、時計回りを正方向として測る）なる位置に、独立したＭ個の目標１２（信号送信源または信号反射物体、図１では煩雑さを避ける為、１個で代表させている）が存在しているものとする。
【００１０】
さて、説明を行う便宜上、具体的な装置として波長λなる高周波信号をキャリア信号とするレーダを考える。レーダは、受信用アンテナ１０と空間的に隔てた位置に設置された送信用アンテナから目標に向けて探知信号を送信し、これが目標で反射される事によって生成される、所謂、エコー信号をアンテナ１０で受信し、エコー信号と送信信号とをミキシングしてベースバンド信号に復調し、更に、適切な信号処理を施して各目標の位置や速度情報を推定する。よって図１の状況では、Ｍ（≦Ｎ−１）個の目標からのエコー信号ｘｍ（ｔ）がアンテナ１０に到来する。なお、ｍ＝１，．．．，Ｍであり、ｔは時間である。
【００１１】
このとき、ｎ番目のアンテナＡｎからの入力信号を復調して得られた信号ｖｎ（ｔ）は式１で表される。
【数１】
ただし、
【数２】
である。ここで、ｘｍ（ｔ）はベースバンド信号、ｎｍ（ｔ）は加法性ガウス雑音、λはキャリア信号である。
【００１２】
信号ｖｎ（ｔ）をベースバンド信号ベクトルとして表すと式３となる
【数３】
ただし、
【数４】
【数５】
【数６】
である。ｘ（ｔ）とｎ（ｔ）との間に相関が無いものとして、式３からｖ（ｔ）の共分散行列Ｒｖｖを計算すると、式７となる。
【数７】
ここで、Ｒｘｘはベースバンド信号ベクトルｘ（ｔ）の共分散行列であり、式８で表される。
【数８】
ＲｖｖがＤＯＡ推定する場合の対象となる。
【００１３】
ここではレーダについて考えているので、ＵＬＡで受信された信号は、同じ信号源からの送信された信号が各目標によって反射されてきものであり、コヒーレント性が高い為、Ｎ×Ｎ次行列であるＲｖｖの階数（ランク）は１に縮退する。したがって、このままではＲｖｖからＤＯＡを推定することは難しい。
【００１４】
そこで、Ｒｖｖの主対角線方向にＬ×Ｌ次部分行列（但し、Ｌ＜Ｎ）をＮ−Ｌ＋１個取り、これらの部分行列を足し合わせて平均することにより、Ｒｖｖのランクを回復する、前方空間平均法（forward spatial smoothing：ｆｓｓ）が用いられる。また、ＵＬＡの基準点を反転させｆｓｓと同様の操作を行う後方空間平均法（backward spatial smoothing：ｂｓｓ）も用いられる。通常は双方を組み合わせ、ｆｂｓｓ（forward backward spatial smoothing）として用いられる。
【００１５】
Ｒｖｖにｆｂｓｓを適用して得られるＲｖｖｆｂｓｓと表すことにすれば、例えば、ｃａｐｏｎ法によるＤＯＡ推定には、式９が用いられる。ただし、ａ（θ）はＬ次のモードベクトルであり、式４を構成するベクトル要素ａ（θｍ）と同様の構造を持つ。
【数９】
ｍｕｓｉｃ（Multiple Signal Classification）法では、Ｒｖｖｆｂｓｓを式１０のように固有値分解して行列ＥＮを求め、式１１によってＤＯＡ推定を行う。
【数１０】
【数１１】
すなわち、式９または式１１において、θをパラメータとするモードベクトルa(θ)=[1,exp(j2παsinθ)],…,exp(j2πα(L-1)sinθ)]Tを定義し、θを走査しながら式１１を計算することで行列に含まれる角度情報を調べ、ＤＯＡを推定する。ただし、α＝ｄ／λである。
【００１６】
通常、信号共分散行列には式７の様に雑音成分σ２Ｉが不可避的に含まれるが、その影響を最小限に抑えてＤＯＡの推定精度を向上させる為、特許文献１の技術では以下のような処理を行っている。すなわち、ＵＬＡを構成する各アンテナで受信した信号を復調して得られたベースバンド信号から相関ベクトルを計算し、擬似共分散行列を生成する。次いで、この擬似共分散行列からプロパゲータ行列を生成し、ＤＯＡ推定を行う。以下にこの推定方法について詳しく説明する。
【００１７】
ここでは、簡単のため、ｆｓｓのみについて考える。ベースバンド信号の相関ベクトルｒｖ１を式１２のように計算する。
【数１２】
これ以降、Ｍは、Ｎ個のアンテナで構成されたＵＬＡによって、ＤＯＡを推定する事が可能な最大の目標数を表すものとする。例えば、Ｍは（Ｎ−１）／２以下の最大の自然数である。
【００１８】
次に、ｒｖ１の要素を並べ替えて、Ｒｆ１を式１３のように生成する。
【数１３】
ここで、
【数１４】
であるから、Ｒｆ１はＲｆ１＝ＡＸなる形をしている。すなわち、式１５に示す行列Ｒは、空間平均を適用した共分散行列Ｒｖｖと同じ位相情報を含んでいる。
【数１５】
【００１９】
以下、Ｎ＝５、Ｍ＝２の場合を例に説明する。このとき行列Ｒは式１６となる。
【数１６】
次に、式１７のように、行列ＲをＭ×Ｍ次部分行列Ｒ１および（Ｎ−２Ｍ）×Ｍ次部分行列Ｒ２とに分割する。
【数１７】
部分行列Ｒ１およびＲ２を用いて、Ｍ×（Ｎ−２Ｍ）次行列Γ＝（Ｒ１Ｒ１Ｈ）−１Ｒ１Ｒ２Ｈを生成する。
【００２０】
行列ΓとＮ−２Ｍ次単位行列：ＩＮ−２Ｍから（Ｎ−Ｍ）×（Ｎ−２Ｍ）次プロパゲータ行列Π＝［Γ⊥−ＩＮ−２Ｍ］Ｔを生成する。
【００２１】
さらに、この行列Πを用いて、例えば、（Ｎ−Ｍ）×（Ｎ−Ｍ）次核行列Ωを、Ω＝Π（ΠＨΠ）−１ΠＨによって定義し、これとアレイモードベクトルa(θ)=[1,exp(j2παsinθ)],…,exp(j2πα(N-M-1)sinθ)]Tとを用いて、例えば、角度スペクトラムＰ（θ）を式１８のように定義する。ただし、α＝ｄ／λである。
【数１８】
そして、パラメータθを走査しながら角度スペクトラムＰ（θ）を計算し、そのピーク位置から到来信号の到来角度θｍを得ることができる。または、ｚ＝ｅｘｐ（ｊ２παｓｉｎ（θ））、α＝ｄ／λ、に対し、代数方程式ａ（１／ｚ）ＴΩａ（ｚ）＝０の解から到来信号の到来角度θｍを得ることができる。
【００２２】
以上のように、特許文献１では、式１５および１６のように行列Ｒを求め、行列Ｒからプロパゲータ行列Πを求め、核行列Ωを求め、到来方向を算出している。
【００２３】
説明を簡単にするため、以下式１６に示した様に、Ｎ＝５、Ｍ＝２の場合を例に取る。プロパゲータ行列Πの主要部分であるΓは式１９のようになる。
【数１９】
ＤＯＡ推定を行うためには、Γを構成するベクトルの位相成分があれば十分である。しかしながら、式１９に掛かっている因子１／｜ｖ５｜２は、ＤＯＡ推定に直接関与しないばかりか、計算の高速化を阻害さえする。また実際の計算に於いては、桁落ちや丸め誤差等の原因となり、高精度化を阻害する。
【００２４】
そこで、本到来方向推定装置は、上記因子を含まない行列を用いてＤＯＡ推定を行うことを特徴とする。
【００２５】
１≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、Ｎ−Ｍを列の次元とし、ベースバンド信号ベクトルｖ１〜ｖＮの部分ベクトルｖ１〜ｖＮ−１を構成要素として並べ、式２０のようなＨａｎｋｅｌ行列Ｒｆ１を生成する。
【数２０】
【００２６】
具体例として、Ｎ＝５およびＭ＝２の場合についてＲｆ１を式２１に示す。
【数２１】
【００２７】
次に、例えば、式２２に示す様に行列Ｒを定義し、このＲをＭ×Ｍ次部分行列Ｒ１および（Ｎ−２Ｍ）×Ｍ次部分行列Ｒ２に分割する。
【数２２】
【００２８】
上記部分行列Ｒ１およびＲ２から、Ｍ×（Ｎ−２Ｍ）次行列Γ＝（Ｒ１Ｒ１Ｈ）−１Ｒ１Ｒ２Ｈを計算すると、式２３のようになる。
【数２３】
式２３の行列には、ＤＯＡ推定に用いられる位相情報が保存されている反面、式１９に含まれていた余分な因子：１／｜ｖ５｜２は含まれない。
【００２９】
以下、式１７と同様に、行列ΓとＮ−２Ｍ次元単位行列ＩＮ−２Ｍとからプロパゲータ行列Π＝［Γ⊥−ＩＮ−２Ｍ］Ｔを生成する。
【００３０】
更に、この行列Πを用いて、例えば、（Ｎ−Ｍ）×（Ｎ−Ｍ）次核行列Ωを、Ω＝Π（ΠＨΠ）−１ΠＨによって定義し、これとアレイモードベクトルa(θ)=[1,exp(j2παsinθ)],…,exp(j2πα(N-M-1)sinθ)]Tとを用いて、例えば、角度スペクトラムＰ（θ）を式１８のように定義する。但し、α＝ｄ／λである。
【数２４】
そして、パラメータθを走査しながら角度スペクトラムＰ（θ）を計算し、そのピーク位置から到来信号の到来角度θｍを得ることができる。または、ｚ＝ｅｘｐ（ｊ２παｓｉｎ（θ））、α＝ｄ／λ、として、代数方程式ａ（１／ｚ）ＴΩａ（ｚ）＝０の解から到来信号の到来角度θｍを得ることができる。
【００３１】
本到来方向推定装置によれば、行列Ｒｆ１として、式１３の代わりに式２１のような行列Ｒｆ１を用いる。これにより式２３のΓは、式１９の１／｜ｖ５｜２のような因子を含まない。よって、Γの計算、プロパゲータ行列Πの計算を高速高精度に行うことができる。以下、上記原理に基づく本発明の実施例について説明する。
【実施例１】
【００３２】
図２は実施例１に係る到来方向推定装置のブロック図である。図３は到来方向推定装置が行う処理のフローチャートである。図２および図３を参照して説明を行う。図２に於いて、到来方向推定装置１００は、複数のアンテナ１０（アンテナＡ１〜ＡＮ）、受信部２０および到来方向推定部３０を含んでいる。複数のアンテナ１０は目標からの到来信号を受信する。受信部２０は受信器２２および信号ベクトル生成手段２４を含んでいる。受信器２２は、例えば、キャリア信号の周波数が７６ＧＨｚである到来信号を受信し、変調方式に応じた復調を行い、適切な帯域のベースバンド信号にダウンコンバートする。そして更に、このベースバンド信号をＡ／Ｄ変換してデジタル信号に変換する（以下、単にベースバンド信号と記す）。信号ベクトル生成手段２４は、各アンテナに対応するベースバンド信号をベクトル列として整形する。以上により、信号ベクトル生成手段２４は式３のベースバンド信号ベクトルｖ１〜ｖＮを生成する（図３のステップＳ１０）。
【００３３】
到来方向推定部３０は、到来信号数予測／設定手段３１、Ｈａｎｋｅｌ行列生成手段３２および推定手段４５を含んでおり、ベースバンド信号ベクトルからＤＯＡを推定する。推定手段４５は、プロパゲータ行列生成手段３３、核行列生成手段３４および到来方向算出手段３５を含んでいる。到来方向推定部３０はコンピュータのＣＰＵ（中央演算処理装置）等の演算装置、演算装置に実行させるプログラム、或いはその記憶装置上の実体である。また、受信部２０と到来方向推定部３０との全てが、ソフトウェア無線装置と同様の意味で、演算装置から構成されていても良い。
【００３４】
到来信号数予測／設定手段３１は、到来信号数Ｍを予測、若しくは設定して出力する（図３のステップＳ１２）。到来信号数Ｍを予測する場合は、最尤法に基づくＡＩＣ（Akaike Information Criteria）、ＭＤＬ（Minimum Description Length）、若しくは他の適切な指標を用いることができる。到来信号数予測／設定手段３１は、到来信号数Ｍの予測を行わない場合は、アンテナ数Ｎから自動的に決定される値である［（Ｎ−１）／２］をＭの値としてもよい。なお、この数式において記号［ｘ］は、実数ｘに対して、ｘ以下の最大の自然数を戻す事を意味する。
【００３５】
推定手段４５は、先ず、式３の信号ベクトルのＳＮ比が事前に決定しておいた或る閾値以上かを判定し（ステップＳ１４）、ＳＮ比が閾値より小さい場合には、例えば、相関ベクトルｃｏｎｊ（ｖＮ）×［ｖ１，…，ｖＮ］Ｔを算出し、これを信号ベクトルとして置き直す（ステップＳ１６）。ここでｃｏｎｊ（ｖＮ）は、信号ｖＮの共役複素数を意味する。
【００３６】
Ｈａｎｋａｌ行列生成手段３２は、式２０のようなＨａｎｋｅｌ行列Ｒｆ１を生成する（ステップＳ１６に於いて信号ベクトルを置き換えた場合は、式１６の様な構造となるが、今後は特に触れない）。プロパゲータ行列生成手段３３では、行列Ｒ＝［Ｒｆ１］とする（ステップＳ１８）。
【００３７】
プロパゲータ行列生成手段３３は、行列Ｒを２つの部分行列Ｒ１およびＲ２に分解し、前述のように行列Γ＝（Ｒ１Ｒ１Ｈ）−１Ｒ１Ｒ２Ｈを生成し、さらに行列Γと単位行列ＩＮ−２Ｍとからプロパゲータ行列Π＝［Γ⊥−ＩＮ−２Ｍ］Ｔを生成する（ステップＳ２２）。
【００３８】
核行列生成手段３４は、核行列Ω＝Π（ΠＨΠ）−１ΠＨを生成する（ステップＳ２４）。
【００３９】
到来方向算出手段３５は、パラメータθを含むアレイモードベクトルa(θ)=[1,exp(j2παsinθ)],…,exp(j2πα(N-M-1)sinθ)]Tと核行列Ωとを用いて、式２４のような角度スペクトラムＰ（θ）を生成する（ステップＳ２６）。但し、α＝ｄ／λである。到来方向算出手段３５は、θを走査しながら角度スペクトラムＰ（θ）を計算し、角度スペクトラムＰ（θ）のピーク位置から到来信号の到来角度θｍを計算する。または、ｚ＝ｅｘｐ (ｊ２παｓｉｎ（θ）)、α＝ｄ／λ、として、代数方程式ａ（１／ｚ）ＴΩａ（ｚ）＝０の解から到来信号の到来角度θｍを計算する。到来方向算出手段３５は、到来角度θｍを推定した到来方向として出力する（ステップＳ２８）。
【００４０】
Ｈａｎｋａｌ行列生成手段３２は、信号ベクトルの他の部分ベクトルｖ２〜ｖＮを構成要素とし、式２５のようなＨａｎｋｅｌ行列Ｒｆ２をＲｆ１の代わりに生成し、プロパゲータ行列生成手段３３は、このＲｆ２を用いてプロパゲータ行列Πを生成してもよい。
【数２５】
【００４１】
また、Ｈａｎｋｅｌ行列生成手段３２は、共役複素部分ベクトルｖ１＊〜ｖＮ−１＊を構成要素とするＨａｎｋｅｌ行列、Ｒｂ１＝ＪＮ−ＭＲ＊ｆ１ＪＭをＲｆ１の代わりに生成し、プロパゲータ行列生成手段３３は、このＲｂ１を用いてプロパゲータ行列Πを生成してもよい。なお、Ｒ＊ｆ１は、Ｒｆ１の各要素について共役複素数をとった行列である。
【００４２】
一例として、Ｎ＝５かつＭ＝２の場合を説明すれば、Ｒｂ１は式２６のように表される。
【数２６】
次に、行列Ｒを式２７のように定義し、更に、行列ＲをＭ×Ｍ次部分行列Ｒ１および（Ｎ−２Ｍ）×Ｍ次部分行列Ｒ２に分割する。
【数２７】
【００４３】
上記部分行列Ｒ１およびＲ２から行列Γ＝（Ｒ１Ｒ１Ｈ）−１Ｒ１Ｒ２Ｈを計算すると、式２８のようになる。
【数２８】
式２８においても式２３と同様に、ＤＯＡ推定に用いられる位相情報が保存され、式１９の１／｜ｖ５｜２のような因子は含まない。よって、Γの計算、プロパゲータ行列Πの計算を高速高精度に行うことができる。
【００４４】
さらに、Ｈａｎｋａｌ行列生成手段３２は、共役複素部分ベクトルｖ２＊〜ｖＮ＊を用いてＨａｎｋｅｌ行列Ｒｂ２＝ＪＮ−ＭＲ＊ｆ２ＪＭをＲｆ１の代わりに生成し、プロパゲータ行列生成手段３３は、このＲｂ２を用いてプロパゲータ行列Πを生成してもよい。なお、Ｒ＊ｆ２は、Ｒｆ２の各要素について共役複素数をとった行列である。
【００４５】
プロパゲータ行列生成手段３３は、（Ｎ−Ｍ）×Ｍ次のＨａｎｋｅｌ行列Ｒｆ１、Ｒｆ２、Ｒｂ１、またはＲｂ２をＫ個組み合わせ、（Ｎ−Ｍ）×（Ｋ×Ｍ）の行列として行列Ｒを構成してもよい。例えば、Ｋ＝２の場合Ｒ＝［Ｒｆ１⊥Ｒｆ２］、Ｒ＝［Ｒｂ１⊥Ｒｂ２］、Ｒ＝［Ｒｆ１⊥Ｒｂ２］、またはＲ＝［Ｒｆ２⊥Ｒｆ２］、Ｋ＝４の場合Ｒ＝［Ｒｆ１⊥Ｒｆ２｜Ｒｂ１⊥Ｒｂ２］とすることができる。この場合、Ｋの値に応じて行列ＲをＭ×（Ｋ×Ｍ）の部分行列Ｒ１および（Ｎ−２Ｍ）×（Ｋ×Ｍ）の部分行列Ｒ２に分割し、Ｒ＝［Ｒ１⊥Ｒ２］Ｔとする。
【００４６】
以上のように、プロパゲータ行列生成手段３３は、行列Ｒを行列Ｒｆ１、Ｒｆ２、Ｒｂ１、およびＲｂ２の少なくとも１つを用いて生成し、この行列Ｒを適切な次元の２つの部分行列Ｒ１、Ｒ２に分割する。
【００４７】
実施例１によれば、到来方向推定装置１００は、Ｈａｎｋｅｌ行列生成手段３２がＲｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を生成する。推定手段４５が、行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔなる２つの部分行列Ｒ１およびＲ２に分割し、部分行列Ｒ１とＲ２とを用いて前記行列Γ、Π、そしてΩを計算し、到来信号の到来方向を推定する。これにより、式１９の１／｜ｖ５｜２のような因子を含まないΓ（式２３）を用いる事で、Γの計算、プロパゲータ行列Πの計算を高速高精度に行うことができる。
【実施例２】
【００４８】
図４は、実施例２に係る到来方向推定装置のブロック図である。図２を参照して説明を行えば、到来方向推定装置１００ａは、実施例１の図２と比較し、到来方向推定部３０ａがスケーリング行列生成手段３６を含んでいる。同図の他の構成は、実施例１の図２と同じであるから説明を省略する。
【００４９】
スケーリング行列生成手段３６は、プロパゲータ行列生成手段３３が生成した行列Ｒの部分行列Ｒ１およびＲ２を用い、（Ｎ−２Ｍ）×（Ｎ−２Ｍ）次スケーリング行列Λ＝Ｒ２Ｒ２Ｈ−Ｒ２Ｒ１ＨΓを生成する。核行列生成手段３４は、プロパゲータ行列Πとスケーリング行列Λとから核行列Ω＝Π（Λ）−１ΠＨを生成する。または、Π´＝Π（ΠＨΠ）−１／２とし、核行列をΩ＝Π´（Λ）−１（Π´）Ｈによって生成してもよい。
【００５０】
核行列生成手段３４は、Π（Λ）−１ΠＨの代わりにΠΠＨを核行列Ωとして生成してもよい。この場合、角度スペクトラムは式２９となり、
【数２９】
代数方程式は式３０となる。
【数３０】
ただし、ｚ＝ｅｘｐ(ｊ２παｓｉｎ（θ）)、α＝ｄ／λ、に対し、
【数３１】
である。
【００５１】
核行列生成手段３４は、核行列をΩ＝Π［α（ΠＨΠ）＋（１−α）Λ］−１ΠＨとして生成してもよい。ここで、０≦α≦１である。このように、核行列生成手段３４は、Ｈａｎｋｅｌ行列から導かれる複数の行列を複合させたものとして、核行列を生成してもよい。
【００５２】
ここで、α≠０、１とし、ΠＨΠとΛとが特異点でないとして、Ａ＝ΠＨΠ、Ｂ＝Λ、１−α＝βとすると、式３２となる。
【数３２】
以上より上の核行列Ωは、公知の核行列Π（ΠＨΠ）−１ΠＨおよびΠ（Λ）−１ΠＨの加重平均となる。よって、公知の核行列を単独で用いたのでは推定困難な目標についても角度推定が可能となる。
【００５３】
核行列生成手段３４は、核行列Ωの生成を、実施例１ではプロパゲータ行列Πから行い、また、実施例２ではプロパゲータ行列Πとスケーリング行列Λとから行う。
【００５４】
到来方向算出手段３５は、プロパゲータ行列Πと、Ｄ、（ＢＨＡ−１Ｂ）、（ＢＨＢ）、（Ｄ−ＢＨＢ）、Λ（ΠＨΠ）等、以下に定める行列Ａ，Ｂ，Ｃ、Ｄを適当に組み合わせた代用行列をスケーリング行列Λの代わりに用いて核行列Ωを生成し、上記実施例１と同様にして、角度スペクトラムまたは代数方程式から到来方向を算出してもよい。なお、Ａ＝Ｒ１Ｒ１Ｈ、Ｂ＝Ｒ１Ｒ２Ｈ、Ｃ＝Ｒ２Ｒ１Ｈ、Ｄ＝Ｒ２Ｒ２Ｈである。
【００５５】
例としてＤを用いた角度スペクトラムを式３３、代数方程式を式３４に示す。
【数３３】
【数３４】
ただし、ｚ＝ｅｘｐ(ｊ２παｓｉｎ（θ）)、α＝ｄ／λ、に対し、
【数３５】
である。
【００５６】
さらに、到来方向算出手段３５は、スケーリング行列Λ、その代用行列、またはスケーリング行列Λのノルムを用いて到来方向を算出してもよい。
【００５７】
図５から８は、実施例２をＦＭＣＷ（Frequency Modulated Continuous Wave）レーダに適用し、目標定位を実施した例である（これらの図では、距離の単位はｂｉｎである）。
【００５８】
図５（ａ）から（ｄ）は、目標が２つ存在する場合の角度スペクトラムを角度と距離に対し示した図である。設定した目標の距離は４０ｍ、角度は０°および３°である。図５（ａ）はＦＦＴ−ＤＢＦ法、図５（ｂ）はＦＢＳＳ−ＭＵＳＩＣ法、図５（ｃ）は特許文献１の方法（比較例）、図５（ｄ）は実施例２の方法を用いて計算した角度スペクトラムである。角度スペクトラムが仮定した２つの目標位置に集約している場合、つまり、角度スペクトラムが２つの目標付近で点状に表示されている程、高い精度で目標の到来方向が推定されていることを示している。
【００５９】
図６は、図５（ａ）、（ｃ）および（ｄ）の角度スペクトラムを距離４０ｍの位置で角度軸に平行に切り出したものを重ね描きしたものである。図５および図６より、比較例（特許文献１）および実施例２の方法は、ＦＦＴ−ＤＢＦ法、ＦＢＳＳ−ＭＵＳＩＣ法と比較して高い精度で目標の角度推定を実現している事が分かる。
【００６０】
図７および図８は、空間的に近接した位置に配置された複数の目標に対する定位性能を、５個の目標を仮定して（図中×の点）角度スペクトラムによって評価した例である。モノクロ表示で若干判別しづらいが、等高線の込み合っている箇所が角度スペクトラムのピークである。図７は特許文献１の方法を用いて計算した角度スペクトラムを示し、図８は、実施例２において、核行列をΩ＝Π［α（Π）＋（１−α）Λ］−１ΠＨ（α＝０．５）として計算した角度スペクトラムを示している。図８は図７と比較して、×印の周辺の等高線がより込み合っており、仮定した目標付近に角度スペクトラムのピークがある。このように、核行列として複合された核行列を用いることにより、エコー信号のコヒーレンスが高く、推定が難しい多くの目標を正しく検出することができる。
【実施例３】
【００６１】
実施例３は、自動車等の車両に搭載された到来方向推定装置の例である。実施例３に係る到来方向推定装置１００ｂは、到来方向推定装置１００ｂを含む、例えば、ミリ波レーダを搭載した自車に対して脅威のある目標の角度（到来方向）を出力する。目標としては、例えば自車の前を走行する自動車等の車両である。図９は、実施例３に係る到来方向推定装置１００ｂのブロック図であるが、到来方向推定装置１００ｂは、受信部２０と到来方向推定部３０の間に脅威度評価手段４０を含んでいる。脅威度評価手段４０は、距離／速度／粗到来方向推定手段３７、目標選別手段３８およびビーム形成手段３９を含んでいる。
【００６２】
図１０は、実施例３に係る到来方向推定装置１００ｂの処理を示すフローチャートである。以下、図１０を参照しながら説明を行う。例えば、ミリ波信号にＦＭ変調を掛ける為に用いられる信号の１周期分に渡り、受信アンテナ毎にエコー信号のサンプリングを行い、信号ベクトル生成手段２４が対応する信号ベクトルのデータセットを抽出する（ステップＳ３０）。即ち、アンテナ数をＮ、サンプリング数をＱとすると、このデータはＮ×Ｑ次行列（データセット）となる。距離／速度／粗到来方向推定手段３７は、時間領域に沿ってＱ次元ベクトルがＮ組並んだデータセットをそれぞれのアンテナ毎にＦＦＴ（高速フーリエ変換）処理する（ステップＳ３２） − Ｑ次元ベクトルのＦＦＴをＮ回行い、周波数領域のデータセットに変換されたＮ×Ｑ行列を得る。距離／速度／粗到来方向推定手段３７は、ＦＦＴの結果として周波数領域に現れたピークに対応するｂｉｎ位置（周波数インデックスの値と考えても良い）から、各目標と自車との相対距離ｒｍおよび相対速度ｓｍを推定し、更に、各目標（ｂｉｎ）位置に於けるＮ次元データベクトルを取り出して、例えばＦＦＴ−ＤＢＦ（デジタルビームフォーミング）処理を行って各目標の大雑把な角度位置を把握する（ステップＳ３４）。また、推定した相対距離、相対速度、角度位置を、各目標固有の識別情報としてペアリングする。
【００６３】
目標選別手段３８は、各目標の脅威度を、目標の相対距離ｒｍと相対速度ｓｍとに基づいて評価する（ステップＳ３６）。脅威度はｒｍが小さい程、そして、ｓｍが大きな負の値となる程（即ち、近接した位置に存在し、自車に高速で接近してくる目標程）、大きくなる。脅威度は、例えば式３６の評価関数を用いる。
【数３６】
ただし、
【数３７】
ここで、ｓＲＥＦ、ｒＲＥＦは、車両メーカーの安全基準等に従ってあらかじめ設定された所定値であるが、例えばＤＳＲＣ（Dedicated Short Range Communications）等のインフラ装置からの情報を受け、交通状況に応じて動的に変更する様にしても良い。
【００６４】
目標選別手段３８は、最大の脅威度を有する目標の脅威度が、基準値以上か否かを判断する（ステップＳ３８）。判定がＹｅｓの場合、目標選別手段３８は、先に粗到来方向推定手段に於いて推定した目標の角度を参照して、最大脅威度を有する目標が自車レーン等、自車の動線と交差する位置を走行するか否かを判断する（ステップＳ４０）。
【００６５】
ステップＳ３８またはＳ４０で判定がＮｏとなった場合、脅威度推定手段４０は、周波数軸上で定位される各目標の相対距離ｒｍおよび相対速度ｓｍのデータ（ｂｉｎ）を、後段の到来方向推定部３０に引き渡す。到来方向推定部３０は、指定された各ｂｉｎ位置のＮ次元データベクトルに対して、例えばデジタルビームフォーミング等、精度は比較的低いが、計算量の小さい通常の到来方向推定処理を適用して各目標の角度推定を行う（ステップＳ４２）。もちろん、粗到来方向推定手段による推定値を用いても、実施例１または実施例２の方法を適用しても良い。到来方向推定部３０は、各目標の識別情報である、相対距離、相対速度、および角度推定値のそれぞれを１つ前の測定期間のデータとして更新保存する（ステップＳ４４）。更新されたデータは、次回の測定期間に於いて得られたデータに基づいて、脅威度評価手段４０が各目標の脅威度を再評価する際に参照される。
【００６６】
ステップＳ４０において判定がＹｅｓとなった場合、目標選別手段３８は、最大の脅威度を持つと判定された目標の相対距離ｒｍおよび相対速度ｓｍのデータ（ｂｉｎ）に基づいて、対応するＮ次元データベクトルを抽出してビーム形成手段３９に引き渡す。ビーム形成手段３９は、距離／速度／粗到来方向推定手段３威度を持つ目標の角度情報を参照し、当該目標が存在する方向に、ビームスペース法によるビーム形成を行う（ステップＳ４６）。この様にして、脅威度推定手段４０は、脅威度の高い目標に関するデータベクトルをビームスペースデータに変換し、後段の到来方向推定部３０に引き渡す。到来方向推定部３０は、詳細な角度推定が可能な実施例１または実施例２の方法で到来方向の推定を行う（ステップＳ４８）。装置が、新たにデータ測定を開始した場合、ステップＳ３０に戻って同様な処理を繰り返す。
【００６７】
実施例３の脅威度評価手段４０は、信号ベクトルを用い、目標の自車に対する脅威を示す脅威度を評価する。到来方向推定部３０の推定手段４５は脅威度に基づき、目標の到来方向を推定する、または推定しない。このように、実施例３においては、脅威度評価手段４０を導入する事により、自車の安全に対する目標の重要性を脅威度という指標で評価し、各目標の脅威度に応じて到来方向推定部３０と連携して適切な到来方向推定手法を適用する事で、装置としてのレスポンスタイムを最適化できる。
【実施例４】
【００６８】
実施例４は、信号品質や目標間の相対離角といった追加指標を導入する事で、より、現実的な運用環境に適した処理を実現する装置である。図１１は、実施例４に係る到来方向推定装置１００ｃのブロック図である。同図に於いて、到来方向推定装置１００ｃは、受信部２０と到来方向推定部３０との間に脅威度評価手段５０を含んでおり、ハードウェアとしての基本的な構成は実施例３と殆ど同じであるから説明を省略する。
【００６９】
図１２は、到来方向処理装置１００ｃの処理を示すフローチャートである。脅威度評価手段５０は、信号ベクトルの抽出（ステップＳ５０）から目標毎に相対距離と相対速度の推定（ステップＳ５３）を行うまで、実施例３のそれと同じ処理を行う。但し、実施例３で行っていたＦＦＴ−ＤＢＦによる粗到来方向推定処理は実施しなくても良い。
【００７０】
脅威度評価手段５０の機能として新たに加わったステップは、脅威度の評価（ステップＳ５６）を実施する前の、信号品質の評価（ステップＳ５４）、一つ前の測定期間に推定された目標緒元データの参照（ステップＳ６６）、そして、ステップＳ６６の参照値の相対離角を基準値との比較（ステップＳ６８）、である。
【００７１】
新たに加わった処理の流れは、先ず、ステップＳ５４で、信号ベクトルのＳＮ比が基準値以上か否かを判断し、判定がＹｅｓの場合、ステップＳ５６に進む。一方、判定がＮｏの場合、ステップＳ６６に進み、一つ前の測定期間に推定された目標緒元データのうち到来方向のデータを参照し、ステップＳ６８に於いて参照値の相対離角と基準値との比較を行う。そして、脅威度評価手段５０は、ステップＳ５４の判定がＹｅｓの場合、または、ステップＳ６８の判定がＮｏの場合、各目標について脅威度の評価を行い（ステップＳ５６）、最大の脅威度を有する目標の脅威度が、基準値以上か否かを判断する（ステップＳ５８）。脅威度の評価は、実施例３のステップＳ３６と同じ方法で行うことができる。
【００７２】
脅威度評価手段５０は、ステップＳ５８の判定がＹｅｓである場合、制御を到来方向推定部３０に渡し、詳細な到来方向の推定を行う（ステップＳ６０）。一方、ステップＳ５８の判定がＮｏ、若しくは、ステップＳ６８の判定がＹｅｓである場合、脅威度評価手段５０は、制御を到来方向推定部３０に渡し、粗い到来方向推定を行う（ステップＳ７０）。このように、脅威度評価手段５０は、信号ベクトルを用い目標の自車に対する脅威を示す脅威度を評価する。到来方向推定部３０の推定手段４５は脅威度に基づき、目標の到来方向の推定を粗い角度分解能で行う、または細かい角度分解能で行う。
【００７３】
以下に、ステップＳ６０で行う処理とステップＳ７０で行う処理との違いを示す。行列Ｒの部分行列Ｒ１、Ｒ２を用いて、Ｒ１Ｒ１Ｈ＝Ａ、Ｒ１Ｒ２Ｈ＝Ｂ、Ｒ２Ｒ１Ｈ＝ＢＨ、Ｄ＝Ｒ２Ｒ２Ｈ＝Ｄとし、ＲＲＨに逆行列の定理を適用すると、式３８となり、
【数３８】
これから直ちに、式３９が得られる。
【数３９】
【００７４】
ステップＳ６０では、核行列としてΩ＝（ＲＲＨ）−１−［（Ｒ１Ｒ１Ｈ）−１⊥０Ｍ×Ｎ−２Ｍ；０Ｎ−２Ｍ×Ｍ⊥０Ｎ−２Ｍ］を用い、また、ステップＳ７０では、核行列としてΩ＝（ＲＲＨ）−１を用いて到来方向の推定を行う。ステップＳ７０の行列は信号部分空間成分の成分を含んでいる為、推定精度が劣る。しかしながら、ステップＳ７０の行列を計算する途中で、（Ｒ１Ｒ１Ｈ）−１の計算をしておく事は比較的容易であるから、実施例４の方法は推定精度の切り替えが簡単に実現できる。
【００７５】
以上、実施例４によれば、エコー信号の品質（ＳＮ比）が基準値以上であり、かつ目標の脅威度が基準値以上である場合、ステップＳ６０のように到来方向の詳細な推定を行う。一方、目標の脅威度が小さく詳細な位置の把握が必要が無い場合や、エコー信号のＳＮ比が小さく時間を掛けても精度の良い推定が望めない場合、ステップＳ７０のように到来方向の粗い推定を行う。即ち、実施例４は運用環境に応じてより適切な動作 − レスポンス重視、または精度重視 − を選択できる装置である。
【００７６】
以上、本実施例について詳述したが、本発明は係る特定の実施例に限定されるものではなく、特許請求の範囲に記載された本発明の要旨の範囲内において、種々の変形・変更が可能である。
【図面の簡単な説明】
【００７７】
【図１】図１は受信用センサを示す図である。
【図２】図２は実施例１に係る到来方向推定装置のブロック図である。
【図３】図３は実施例１に係る到来方向推定装置のフローチャートである。
【図４】図４は実施例２に係る到来方向推定装置のブロック図である。
【図５】図５（ａ）から図５（ｄ）は各方法で計算された角度スペクトラムを示す図である。
【図６】図６は距離が一定のときの角度に対する角度スペクトラムを示した図である。
【図７】図７は実施例２と比較例の角度スペクトラムを示した図である。
【図８】図８は実施例２の角度スペクトラムを示した図である。
【図９】図９は実施例３に係る到来方向推定装置のブロック図である。
【図１０】図１０は実施例３に係る到来方向推定装置のフローチャートである。
【図１１】図１１は実施例４に係る到来方向推定装置のブロック図である。
【図１２】図１２は実施例４に係る到来方向推定装置のフローチャートである。
【符号の説明】
【００７８】
１０ アンテナ
１２ 目標
２０ 受信部
３０ 到来方向推定部
３１ 到来信号予想／設定手段
３２ Ｈａｎｋｅｌ行列生成手段
３３ プロパゲータ行列生成手段
３４ 核行列生成手段
３５ 到来方向算出手段
３６ スケーリング行列生成手段
３７ 距離／速度／到来方向推定手段
３８ 目標選別手段
３９ ビーム形成手段
４５ 推定手段
４０、５０ 脅威度評価手段
【特許請求の範囲】
【請求項１】
Ｎ個のセンサが受信した目標からの到来信号よりｖ１〜ｖＮで表されるＮ個のベースバンド信号ベクトルを生成する信号ベクトル生成手段と、
１≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１から（Ｎ−Ｍ）×Ｍ次の行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を生成するＨａｎｋｅｌ行列生成手段と、
前記行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、前記行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔによって２つの部分行列Ｒ１およびＲ２に分割し、前記部分行列Ｒ１とＲ２とから前記到来信号の到来方向を推定する推定手段と、
を具備し、
前記行列Ｒｆ１は前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列であり、
前記行列Ｒｆ２は前記信号ベクトルの部分ベクトルｖ２〜ｖＮを構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｅｌ行列であり、
前記行列Ｒｂ１は、Ｒｆ１＊を前記信号ベクトルの共役複素部分ベクトルｖ１＊〜ｖＮ−１＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ１＝ＪＮ−ＭＲｆ１＊ＪＭであり、
前記行列Ｒｆ２＊は前記信号ベクトルの共役複素部分ベクトルｖ２＊〜ｖＮ＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ２＝ＪＮ−ＭＲｆ２＊ＪＭであることを特徴とする到来方向推定装置。
【請求項２】
前記推定手段は、行列Γ＝（Ｒ１Ｒ１Ｈ）−１Ｒ１Ｒ２Ｈとし、ＩＮ−２ＭをＮ−２Ｍ次単位行列とし、プロパゲータ行列Π＝［Γ⊥−ＩＮ−２Ｍ］Ｔを生成する、またはプロパゲータ行列Π´＝Π（ΠＨΠ）−１／２を生成し、前記プロパゲータ行列ΠまたはΠ´から前記到来方向を推定することを特徴とする請求項１記載の到来方向推定装置。
【請求項３】
前記推定手段は、前記プロパゲータ行列ΠまたはΠ´を用い、核行列Ω＝Π（ΠＨΠ）−１ΠＨまたはΩ＝Π´（Π´ＨΠ´）−１Π´Ｈを生成し、核行列Ωから前記到来方向を推定することを特徴とする請求項２記載の到来方向推定装置。
【請求項４】
前記推定手段は、スケーリング行列Λ＝Ｒ２Ｒ２Ｈ−Ｒ２Ｒ１ＨΓを生成し、プロパゲータ行列ΠまたはΠ´とスケーリング行列Λとから前記到来方向を推定することを特徴とする請求項２記載の到来方向推定装置。
【請求項５】
前記推定手段は、前記プロパゲータ行列ΠまたはΠ´とスケーリング行列Λとから核行列Ω＝Π（Λ）−１ΠＨまたはΩ＝Π´（Λ）−１Π´Ｈを生成し、核行列Ωから前記到来方向を推定することを特徴とする請求項４記載の到来方向推定装置。
【請求項６】
核行列生成手段は、複数の核行列を複合させ、複合した核行列から前記到来方向を推定することを特徴とする請求項３または５記載の到来方向推定装置。
【請求項７】
前記信号ベクトルを用い、前記目標の自車に対する脅威を示す脅威度を評価する脅威度評価手段を具備し、
前記推定手段は前記脅威度に基づき、前記目標の到来方向を推定する、または推定しないことを特徴とする請求項１から６のいずれか一項記載の到来方向推定装置。
【請求項８】
前記信号ベクトルを用い、前記目標の自車に対する脅威を示す脅威度を評価する脅威度評価手段を具備し、
前記推定手段は前記脅威度に基づき、前記目標の到来方向の推定を粗い角度分解能で行う、または細かい角度分解能で行う、ことを特徴とする請求項１から６のいずれか一項記載の到来方向推定装置。
【請求項９】
Ｎ個のセンサが受信した目標からの到来信号よりｖ１〜ｖＮで表されるＮ個のベースバンド信号ベクトルを生成するステップと、
１≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１から（Ｎ−Ｍ）×Ｍ次の行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を生成するステップと、
前記行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、前記行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔによって２つの部分行列Ｒ１およびＲ２に分割し、前記部分行列Ｒ１とＲ２とから前記到来信号の到来方向を推定するステップと、
を有し、
前記行列Ｒｆ１は前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列であり、
前記行列Ｒｆ２は前記信号ベクトルの部分ベクトルｖ２〜ｖＮを構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｅｌ行列であり、
前記行列Ｒｂ１は、Ｒｆ１＊を前記信号ベクトルの共役複素部分ベクトルｖ１＊〜ｖＮ−１＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ１＝ＪＮ−ＭＲｆ１＊ＪＭであり、
前記行列Ｒｆ２＊は前記信号ベクトルの共役複素部分ベクトルｖ２＊〜ｖＮ＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ２＝ＪＮ−ＭＲｆ２＊ＪＭであることを特徴とする到来方向推定方法。
【請求項１０】
Ｎ個のセンサが受信した目標からの到来信号よりｖ１〜ｖＮで表されるＮ個のベースバンド信号ベクトルを生成するステップと、
１≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１から（Ｎ−Ｍ）×Ｍ次の行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を生成するステップと、
前記行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、前記行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔによって２つの部分行列Ｒ１およびＲ２に分割し、前記部分行列Ｒ１とＲ２とから前記到来信号の到来方向を推定するステップと、
をコンピュータに実行させ、
前記行列Ｒｆ１は前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列であり、
前記行列Ｒｆ２は前記信号ベクトルの部分ベクトルｖ２〜ｖＮを構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｅｌ行列であり、
前記行列Ｒｂ１は、Ｒｆ１＊を前記信号ベクトルの共役複素部分ベクトルｖ１＊〜ｖＮ−１＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ１＝ＪＮ−ＭＲｆ１＊ＪＭであり、
前記行列Ｒｆ２＊は前記信号ベクトルの共役複素部分ベクトルｖ２＊〜ｖＮ＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ２＝ＪＮ−ＭＲｆ２＊ＪＭであることを特徴とする到来方向推定プログラム。
【請求項１】
Ｎ個のセンサが受信した目標からの到来信号よりｖ１〜ｖＮで表されるＮ個のベースバンド信号ベクトルを生成する信号ベクトル生成手段と、
１≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１から（Ｎ−Ｍ）×Ｍ次の行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を生成するＨａｎｋｅｌ行列生成手段と、
前記行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、前記行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔによって２つの部分行列Ｒ１およびＲ２に分割し、前記部分行列Ｒ１とＲ２とから前記到来信号の到来方向を推定する推定手段と、
を具備し、
前記行列Ｒｆ１は前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列であり、
前記行列Ｒｆ２は前記信号ベクトルの部分ベクトルｖ２〜ｖＮを構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｅｌ行列であり、
前記行列Ｒｂ１は、Ｒｆ１＊を前記信号ベクトルの共役複素部分ベクトルｖ１＊〜ｖＮ−１＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ１＝ＪＮ−ＭＲｆ１＊ＪＭであり、
前記行列Ｒｆ２＊は前記信号ベクトルの共役複素部分ベクトルｖ２＊〜ｖＮ＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ２＝ＪＮ−ＭＲｆ２＊ＪＭであることを特徴とする到来方向推定装置。
【請求項２】
前記推定手段は、行列Γ＝（Ｒ１Ｒ１Ｈ）−１Ｒ１Ｒ２Ｈとし、ＩＮ−２ＭをＮ−２Ｍ次単位行列とし、プロパゲータ行列Π＝［Γ⊥−ＩＮ−２Ｍ］Ｔを生成する、またはプロパゲータ行列Π´＝Π（ΠＨΠ）−１／２を生成し、前記プロパゲータ行列ΠまたはΠ´から前記到来方向を推定することを特徴とする請求項１記載の到来方向推定装置。
【請求項３】
前記推定手段は、前記プロパゲータ行列ΠまたはΠ´を用い、核行列Ω＝Π（ΠＨΠ）−１ΠＨまたはΩ＝Π´（Π´ＨΠ´）−１Π´Ｈを生成し、核行列Ωから前記到来方向を推定することを特徴とする請求項２記載の到来方向推定装置。
【請求項４】
前記推定手段は、スケーリング行列Λ＝Ｒ２Ｒ２Ｈ−Ｒ２Ｒ１ＨΓを生成し、プロパゲータ行列ΠまたはΠ´とスケーリング行列Λとから前記到来方向を推定することを特徴とする請求項２記載の到来方向推定装置。
【請求項５】
前記推定手段は、前記プロパゲータ行列ΠまたはΠ´とスケーリング行列Λとから核行列Ω＝Π（Λ）−１ΠＨまたはΩ＝Π´（Λ）−１Π´Ｈを生成し、核行列Ωから前記到来方向を推定することを特徴とする請求項４記載の到来方向推定装置。
【請求項６】
核行列生成手段は、複数の核行列を複合させ、複合した核行列から前記到来方向を推定することを特徴とする請求項３または５記載の到来方向推定装置。
【請求項７】
前記信号ベクトルを用い、前記目標の自車に対する脅威を示す脅威度を評価する脅威度評価手段を具備し、
前記推定手段は前記脅威度に基づき、前記目標の到来方向を推定する、または推定しないことを特徴とする請求項１から６のいずれか一項記載の到来方向推定装置。
【請求項８】
前記信号ベクトルを用い、前記目標の自車に対する脅威を示す脅威度を評価する脅威度評価手段を具備し、
前記推定手段は前記脅威度に基づき、前記目標の到来方向の推定を粗い角度分解能で行う、または細かい角度分解能で行う、ことを特徴とする請求項１から６のいずれか一項記載の到来方向推定装置。
【請求項９】
Ｎ個のセンサが受信した目標からの到来信号よりｖ１〜ｖＮで表されるＮ個のベースバンド信号ベクトルを生成するステップと、
１≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１から（Ｎ−Ｍ）×Ｍ次の行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を生成するステップと、
前記行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、前記行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔによって２つの部分行列Ｒ１およびＲ２に分割し、前記部分行列Ｒ１とＲ２とから前記到来信号の到来方向を推定するステップと、
を有し、
前記行列Ｒｆ１は前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列であり、
前記行列Ｒｆ２は前記信号ベクトルの部分ベクトルｖ２〜ｖＮを構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｅｌ行列であり、
前記行列Ｒｂ１は、Ｒｆ１＊を前記信号ベクトルの共役複素部分ベクトルｖ１＊〜ｖＮ−１＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ１＝ＪＮ−ＭＲｆ１＊ＪＭであり、
前記行列Ｒｆ２＊は前記信号ベクトルの共役複素部分ベクトルｖ２＊〜ｖＮ＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ２＝ＪＮ−ＭＲｆ２＊ＪＭであることを特徴とする到来方向推定方法。
【請求項１０】
Ｎ個のセンサが受信した目標からの到来信号よりｖ１〜ｖＮで表されるＮ個のベースバンド信号ベクトルを生成するステップと、
１≦ＭかつＭ≦（Ｎ−１）／２なる自然数Ｍを行の次元として優先的に定め、前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１から（Ｎ−Ｍ）×Ｍ次の行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を生成するステップと、
前記行列Ｒｆ１、Ｒｆ２、Ｒｂ１およびＲｂ２のうち少なくとも１つの行列を用いて行列Ｒを生成し、前記行列ＲをＲ＝［Ｒ１⊥Ｒ２］Ｔによって２つの部分行列Ｒ１およびＲ２に分割し、前記部分行列Ｒ１とＲ２とから前記到来信号の到来方向を推定するステップと、
をコンピュータに実行させ、
前記行列Ｒｆ１は前記信号ベクトルの部分ベクトルｖ１〜ｖＮ−１を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列であり、
前記行列Ｒｆ２は前記信号ベクトルの部分ベクトルｖ２〜ｖＮを構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｅｌ行列であり、
前記行列Ｒｂ１は、Ｒｆ１＊を前記信号ベクトルの共役複素部分ベクトルｖ１＊〜ｖＮ−１＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ１＝ＪＮ−ＭＲｆ１＊ＪＭであり、
前記行列Ｒｆ２＊は前記信号ベクトルの共役複素部分ベクトルｖ２＊〜ｖＮ＊を構成要素とする（Ｎ−Ｍ）×Ｍ次のＨａｎｋｌｅ行列、ＪＮ−ＭおよびＪＭをそれぞれＮ−Ｍ次およびＭ次反対角単位行列としたとき、Ｒｂ２＝ＪＮ−ＭＲｆ２＊ＪＭであることを特徴とする到来方向推定プログラム。
【図１】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【公開番号】特開２００９−２４３９４７（Ｐ２００９−２４３９４７Ａ）
【公開日】平成２１年１０月２２日（２００９．１０．２２）
【国際特許分類】
【出願番号】特願２００８−８８０９１（Ｐ２００８−８８０９１）
【出願日】平成２０年３月２８日（２００８．３．２８）
【出願人】（０００００５２２３）富士通株式会社 (25,993)
【Ｆターム（参考）】
【公開日】平成２１年１０月２２日（２００９．１０．２２）
【国際特許分類】
【出願日】平成２０年３月２８日（２００８．３．２８）
【出願人】（０００００５２２３）富士通株式会社 (25,993)
【Ｆターム（参考）】
[ Back to top ]