受信波の波源位置推定装置及び波源方向推定装置、並びに受信波の波源位置推定方法及び波源方向推定方法
【課題】近傍波源の受信波の場合に、計算負荷を低減しつつ、推定精度の劣化が抑制された波源位置推定装置を提供する。
【解決手段】波源位置推定装置は、アレーアンテナ２４ａにより電波又は音波を受信して受信波の波源の位置を推定する。波源位置推定装置は、各受信素子が受信した受信波の入力ベクトルを算出し、代数的手法を用いた波源の位置算出の際に近傍波源の位置特定を可能とする変換行列を算出し、前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の位置を算出する。
【解決手段】波源位置推定装置は、アレーアンテナ２４ａにより電波又は音波を受信して受信波の波源の位置を推定する。波源位置推定装置は、各受信素子が受信した受信波の入力ベクトルを算出し、代数的手法を用いた波源の位置算出の際に近傍波源の位置特定を可能とする変換行列を算出し、前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の位置を算出する。
【発明の詳細な説明】
【技術分野】
【０００１】
この発明は、近傍波源である受信波の波源位置推定装置及び波源方向推定装置、並びに受信波の波源位置推定方法及び波源方向推定方法に関する。
【背景技術】
【０００２】
近年、地上デジタル放送機器、携帯電話、電子キー等の様々な無線機器の使用が日常化しており、希望の相手の位置を把握するため、また希望の相手からの電波を受信するために伝播特性測定の必要性が高まっている。
【０００３】
例えば、車両においては、無線通信により車両に固有のＩＤコードを発信可能な電子キーを車両のキーとして使用する電子キーシステムが搭載されている。この電子キーシステムでは、電子キーから発信されたＩＤコードを車両が受信すると、車両が電子キーのＩＤコードと自身に予め登録されているＩＤコードとを照合するＩＤ照合を行い、これらＩＤコードが一致してＩＤ照合が成立すれば、車載機器の制御、例えば車両のドアロックの施解錠やエンジンの始動が許可又は実行される。
【０００４】
この電子キーシステムにおいて、電子キーの位置が車内と車外とのいずれにあるかで車載機器の制御を切り分けている場合には、車両からのＩＤコードの返信要求を車内と車外とで異なるタイミングで発信して対応している。そこで、車両が電子キーから発信された電波の波源の位置を推定できれば、ＩＤコードの返信要求を異なるタイミングで発信する動作を車両に実行させる必要がなくなるので、この利点を享受するために、波源位置推定によりキー位置を特定する技術が求められている。
【０００５】
波源位置推定を行う方法としては、ＭＵＳＩＣ（ＭＵｌｔｉｐｌｅ ＳＩｇｎａｌ Ｃｌａｓｓｉｆｉｃａｔｉｏｎ）法やＥＳＰＲＩＴ（Ｅｓｔｉｍａｔｉｏｎ ｏｆ Ｓｉｇｎａｌ Ｐａｒａｍｅｔｅｒｓ ｖｉｓａ Ｒｏｔａｔｉｏｎａｌ Ｉｎｖａｒｉａｎｃｅ Ｔｅｃｈｎｉｑｕｅｓ）法等が挙げられるが、これらの方法では波源が十分遠方であり、到来波が平面波として近似できることを前提としている（例えば、非特許文献１参照）。
【０００６】
そこで、モードベクトルを近傍波源に対応したものに変更することで近傍波源位置推定を可能とした、球面波モードベクトルＭＵＳＩＣ法が提案されている。しかし、ＭＵＳＩＣ法を近傍波源に対して適用する場合、２次元のピークサーチが必要となるため、その計算負荷が懸念される。
【先行技術文献】
【非特許文献】
【０００７】
【非特許文献１】菊間信良著「アレーアンテナによる適応信号処理」科学技術出版
【発明の概要】
【発明が解決しようとする課題】
【０００８】
一方、ＥＳＰＲＩＴ法に基づいた近傍波源位置推定法も提案されている。これは、ＥＳＰＲＩＴ法に類似したＤＯＡ（Ｄｉｒｅｃｔｉｏｎ Ｏｆ Ａｒｒｉｖａｌ）−Ｍａｔｒｉｘ法を近傍波源に対して適用可能としたものであり、代数的に位置を算出することができる。したがって、ピークサーチが不要であり、球面波モードベクトルＭＵＳＩＣ法で懸念される計算負荷の低減に有効である。しかし、この方法ではアレー形状が等間隔リニアアレーであることが要求され、また受信波の位相項を２次の級数近似として扱うため、アレー素子間隔が１／４波長以下であることが要求される。更に、アレー開口の増大に伴い近似誤差も増大するという課題がある。そこで、アレー形状に制約が無く、波源が近傍である受信波の位置や方向の推定において計算負荷が低減されつつ、推定精度の劣化が抑制された推定装置が求められていた。
【０００９】
この発明は、こうした実情に鑑みてなされたものであり、その目的は、近傍波源の受信波の場合に、計算負荷を低減しつつ、推定精度の劣化が抑制された波源位置推定装置及び波源方向推定装置、並びに受信波の波源位置推定方法及び波源方向推定方法を提供することにある。
【課題を解決するための手段】
【００１０】
以下、上記目的を達成するための手段及びその作用効果について説明する。
請求項１に記載の発明は、複数の受信素子により電波又は音波を受信して受信波の波源の位置を推定する波源位置推定装置において、各前記受信素子が受信した受信波の入力ベクトルを算出する入力ベクトル算出手段と、代数的手法を用いた波源の位置算出の際に近傍波源の位置特定を可能とする変換行列を算出する変換行列算出手段と、前記入力ベクトルと前記変換行列とを基に前記代数的手法によって波源の位置を算出する波源位置算出手段とを備えることをその要旨としている。
【００１１】
同構成によれば、複数の受信素子により電波又は音波を受信して波源の位置を算出するにあたって、各受信素子が受信した受信波の入力ベクトルを入力ベクトル算出手段により算出し、代数的手法で近傍波源の位置特定を可能とする変換行列を変換行列算出手段によって算出し、波源位置算出手段により入力ベクトルと変換行列とを基に代数的手法としてＥＳＰＲＩＴ法における代数的処理によって波源の位置を算出する。従来、近傍波源の位置を推定する場合にはＭＵＳＩＣ法等のサーチ法が用いられるが、サーチに要する計算負荷が大きいという問題点がある。また、代数的手法により近傍波源の位置を推定する場合には、受信信号を近似する必要があり、近似に起因する誤差が発生する。また、受信素子配置（アレー）形状も等間隔リニアアレーかつ素子間隔が１／４波長以下である必要がある。そこで、本装置では、近傍波源の受信波の場合に、変換行列を使用することで代数的手法の適用を可能とするので、任意のアレー形状において計算負荷を低減できる。また、演算に際して入力値を近似することもないので、高い位置推定精度を確保することができる。
【００１２】
請求項２に記載の発明は、複数の受信素子により電波又は音波を受信して受信波の波源の方向を推定する波源方向推定装置において、各前記受信素子が受信した受信波の入力ベクトルを算出する入力ベクトル算出手段と、代数的手法を用いた波源の方向算出の際に波源の方向特定を可能とする変換行列を算出する変換行列算出手段と、前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の方向を算出する波源方向算出手段とを備えることをその要旨としている。
【００１３】
同構成によれば、複数の受信素子により電波又は音波を受信して波源の方向を算出するにあたって、入力ベクトル算出手段により各受信素子が受信した受信波の入力ベクトルを算出し、代数的手法で近傍波源の方向特定を可能とする変換行列を変換行列算出手段によって算出し、波源方向算出手段により入力ベクトルと変換行列とを基に代数的手法としてＥＳＰＲＩＴ法における代数的処理によって波源の方向を算出する。従来、到来波の方向を代数的手法により推定する場合には、受信素子配置（アレー）形状が、例えば等間隔リニアアレー等に制約される。しかしながら、本装置では、任意形状のアレーで受信した場合においても、変換行列を用いることにより代数的手法の適用を可能とするため、計算負荷を低減できる。
【００１４】
請求項３に記載の発明は、請求項１に記載の波源位置推定装置において、前記変換行列算出手段は、前記代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値と、波源の存在し得る位置とが、一対一対応となる前記変換行列を算出することをその要旨としている。
【００１５】
同構成によれば、波源位置を固有値から算出するため、簡単な代数計算を行えばよく、計算負荷を低減できる。
請求項４に記載の発明は、請求項３に記載の波源位置推定装置において、前記固有値は、前記波源の位置を表す複素数であることをその要旨としている。
【００１６】
同構成によれば、固有値の実部が波源位置のｘ座標、虚部がｙ座標を表すため、固有値から波源の位置を算出する際に、さらに簡単に計算することができる。
請求項５に記載の発明は、請求項１又は３に記載の波源位置推定装置において、前記変換行列は、前記入力ベクトルのうち、一部のみを変換することが可能な行列であることをその要旨としている。
【００１７】
同構成によれば、変換行列は入力ベクトルのうち、一部のみを変換することが可能である。したがって、変換行列の要素数が少なく、変換行列の算出に掛かる計算負荷を低減することができる。
【００１８】
請求項６に記載の発明は、請求項１に記載の波源位置推定装置において、前記変換行列算出手段は、前記代数的手法で用いられる回転不変式の対角行列の対角成分として算出される固有値を考慮することなく前記変換行列を算出することをその要旨としている。
【００１９】
同構成によれば、変換行列の計算に用いる固有値の値を、入力ベクトルから算出する必要がないため、この部分の計算負荷が軽減される。
請求項７に記載の発明は、請求項１、及び３〜５のいずれか一項に記載の波源位置推定装置において、前記波源位置算出手段は、前記代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値によって前記波源の位置を算出することをその要旨としている。
【００２０】
同構成によれば、代数的手法で用いられる回転不変式、すなわち受信素子全体の平行移動によって生じる各波の位相回転、及び強度の変化が成立する式の対角行列の成分として算出された固有値から波源位置を算出する際に、簡単な代数演算を行えばよいため、計算負荷が軽減される。
【００２１】
請求項８に記載の発明は、請求項１、及び３〜６のいずれか一項に記載の波源位置推定装置において、前記波源位置算出手段は、前記代数的手法で得られるモードベクトルによって前記波源の位置を算出することをその要旨としている。
【００２２】
同構成によれば、代数的手法で算出されたモードベクトルから波源位置を算出するため、変換行列を計算する際に固有値について考慮する必要がなく、変換行列算出のための計算負荷が軽減される。
【００２３】
請求項９に記載の発明は、請求項１、及び３〜８のいずれか一項に記載の波源位置推定装置において、前記波源位置算出手段は、前記入力ベクトルの自己相関行列を固有値展開して得られる信号部分空間に属する固有ベクトルを前記変換行列により変換して波源位置を算出することをその要旨としている。
【００２４】
同構成によれば、波源位置算出手段が行う固有値展開の演算量を削減することができる。
請求項１０に記載の発明は、請求項１、及び３〜８のいずれか一項に記載の波源位置推定装置において、前記波源位置算出手段は、前記入力ベクトルを前記変換行列により変換して波源位置を算出することをその要旨としている。
【００２５】
同構成によれば、波源位置算出手段で行われる入力ベクトルの自己相関行列に対する固有値展開において、一般固有値展開を適用することができる。したがって、変換行列がユニタリ行列でない場合においても、雑音による推定精度の劣化を最小限に留めることができる。
【００２６】
請求項１１に記載の発明は、請求項１、及び３〜１０のいずれか一項に記載の波源位置推定装置において、前記変換行列算出手段は、窓関数を用いて前記変換行列を算出することをその要旨としている。
【００２７】
同構成によれば、窓関数を用いて変換行列を算出するため、推定領域全体を変換対象として変換行列を算出した場合と比較し、変換の際に発生する変換誤差を抑制することができる。
【００２８】
請求項１２に記載の発明は、請求項１１に記載の波源位置推定装置において、前記窓関数の位置を初期推定に基づいて算出することをその要旨としている。
同構成によれば、窓関数の位置を初期推定に基づいて算出するため、初期推定を行わずに窓関数の位置を決定した場合と比較し、発生する変換誤差をより抑制することができる。
【００２９】
請求項１３に記載の発明は、請求項１、及び３〜１２のいずれか一項に記載の波源位置推定装置において、前記波源位置推定を繰り返し行うことをその要旨としている。
同構成によれば、波源位置推定を繰り返し行うため、位置の推定精度を更に向上させることができる。
【００３０】
請求項１４に記載の発明は、請求項１３に記載の波源位置推定装置において、前記波源位置推定を繰り返す際に、推定範囲を狭くすることをその要旨としている。
同構成によれば、波源位置推定を繰り返す際に推定範囲を狭くするため、繰り返す度に位置の推定精度を更に向上させることができる。
【００３１】
請求項１５に記載の発明は、請求項１、及び３〜１４のいずれか一項に記載の波源位置推定装置において、前記受信素子は、円形に配置されることをその要旨としている。
同構成によれば、窓関数内部における変換誤差の不均一性を抑制することができるため、位置の推定精度を更に向上させることができる。
【００３２】
請求項１６に記載の発明は、複数の受信素子により電波又は音波を受信して受信波の波源の位置を推定する波源位置推定方法において、各前記受信素子が受信した受信波の入力ベクトルを算出し、異なる前記受信素子又は異なる波源における受信波に対応するよう変換することが可能な変換行列を算出し、前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の位置を算出することをその要旨としている。
【００３３】
請求項１７に記載の発明は、複数の受信素子により電波又は音波を受信して受信波の波源の方向を推定する波源方向推定方法において、各前記受信素子が受信した受信波の入力ベクトルを算出し、異なる前記受信素子又は異なる波源における受信波に対応するよう変換することが可能な変換行列を算出し、前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の方向を算出することをその要旨としている。
【発明の効果】
【００３４】
本発明によれば、近傍波源の受信波の場合に、計算負荷を低減しつつ、推定精度の劣化が抑制された波源位置推定装置及び波源方向推定装置、並びに受信波の波源位置推定方法及び波源方向推定方法を提供することができる。
【図面の簡単な説明】
【００３５】
【図１】電子キーシステムの構成を示すブロック図。
【図２】車両の通信エリアを示す上面図。
【図３】アレーアンテナと波源との位置関係を示す図。
【図４】初期推定位置と窓関数の位置関係。
【図５】シミュレーション条件を示す図。
【図６】本実施形態の信号対雑音比と位置推定誤差との関係を示す図。
【図７】本実施形態の繰り返し回数と位置推定誤差との関係を示す図。
【図８】本実施形態の繰り返し回数と位置推定誤差との関係を示す図。
【発明を実施するための形態】
【００３６】
（第１の実施形態）
以下、本発明にかかる波源位置推定装置及び波源位置推定方法を車両の電子キーシステムに適用した第１の実施形態について図１〜図８を参照して説明する。
【００３７】
図１に示されるように、車両２には、例えば運転者が実際に車両キーを操作しなくてもドアロックの施解錠やエンジンの始動及び停止等の車両動作を行うことが可能な電子キーシステム、いわゆるキー操作フリーシステムが搭載されている。キー操作フリーシステムは、キー固有のＩＤコードを無線通信で発信可能な電子キー１が車両キーとして使用され、電子キー１と車両２との狭域無線通信によりＩＤ照合が実行される。
【００３８】
キー操作フリーシステムには、車両２に近づいたり離れたりした際に自動でＩＤ照合が行われてドアロックの施解錠が許可されるエントリーシステムがある。これを以下に説明すると、車両２には、電子キー１との間で狭域無線通信（以下、スマート通信という）を行う際にＩＤ照合を行う照合ＥＣＵ（Ｅｌｅｃｔｒｏｎｉｃ Ｃｏｎｔｒｏｌ Ｕｎｉｔ）２１と、車両２の状態を管理するメインボディＥＣＵ３１とが設けられている。照合ＥＣＵ２１には、車両２の室内天井に埋設されて（図２参照）車外及び車内にＵＨＦ（Ｕｌｔｒａ Ｈｉｇｈ Ｆｒｅｑｕｅｎｃｙ）帯の信号を発信及び受信可能な通信機２４が接続されている。通信機２４には、複数のアンテナ素子からなるアレーアンテナ２４ａが搭載されている。照合ＥＣＵ２１には、例えばドアロック施解錠等を管理するメインボディＥＣＵ３１が車内ＬＡＮ（Ｌｏｃａｌ Ａｒｅａ Ｎｅｔｗｏｒｋ）３０を介して接続されている。
【００３９】
一方、電子キー１には、車両２との間で電子キーシステムに準じた無線通信を行う際のコントロールユニットとして通信制御部１１が設けられている。通信制御部１１は、固有のキーコードとしてＩＤコードが記憶されたメモリ１１ａを備えている。通信制御部１１には、外部で発信されたＵＨＦ帯の信号を受信可能な受信部１２と、通信制御部１１の指令に従いＵＨＦ帯の信号を発信可能な発信部１３とが接続されている。
【００４０】
照合ＥＣＵ２１は、通信機２４からＵＨＦ帯のリクエスト信号Ｓｒｑを間欠的に発信させ、車両２を含む周辺に通信エリアＡを形成する（図２参照）。電子キー１がこの通信エリアＡに入り込んでリクエスト信号Ｓｒｑを受信部１２で受信すると、電子キー１はリクエスト信号Ｓｒｑに応答する形で、自身のメモリ１１ａに登録されたＩＤコードを含ませたＵＨＦ帯のＩＤコード信号Ｓｉｄを発信部１３から返信する。照合ＥＣＵ２１は、通信機２４でＩＤコード信号Ｓｉｄを受信すると、ＩＤコード信号Ｓｉｄを発信した電子キー１の位置を位置推定部２１ｂにより推定する。位置推定部２１ｂが受信波の波源位置推定装置に相当する。照合ＥＣＵ２１は、自身のメモリ２１ａに登録されたＩＤコードと電子キー１のＩＤコードとを照らし合わせてＩＤ照合を行う。照合ＥＣＵ２１は、電子キー１が車外（Ａ−Ａｉ）に位置すると推定して照合が成立する（車外照合成立）と、メモリ２１ａに車外照合フラグを一定時間立てて、メインボディＥＣＵ３１を介してドアロック装置３８によるドアロックの施解錠を許可又は実行する。
【００４１】
また、電子キーシステムには、エンジン始動停止操作の際に実際の車両キー操作を必要とせずに単なるスイッチ操作のみでエンジン始動停止操作を行うことが可能な機能としてワンプッシュエンジンスタートシステムがある。このワンプッシュエンジンスタートシステムを以下に説明すると、車両２には、照合ＥＣＵ２１のＩＤ照合成立結果を基に、エンジンの点火制御及び燃料噴射制御を行うエンジンＥＣＵ３２が設けられている。エンジンＥＣＵ３２は、車内ＬＡＮ３０を通じて照合ＥＣＵ２１等の各種ＥＣＵに接続されている。車両２の運転席には、車両２の電源状態（電源ポジション）を切り換える際に操作されるエンジンスイッチ３５が設けられている。エンジンスイッチ３５は、押し操作される度に車両２の電源状態をＡＣＣ（Ａｃｃｅｓｓｏｒｙ）オン→ＩＧ（Ｉｇｎｉｔｉｏｎ）オン→電源オフの純に繰り返し遷移させ、エンジン停止時にブレーキペダルが踏み込まれた状態で操作されると、エンジンを始動に切り換える。
【００４２】
照合ＥＣＵ２１は、車外照合が成立してドアロックが解錠された後、ドアが開けられて運転者が乗車したことを例えばカーテシスイッチ３７で認識すると、通信機２４からリクエスト信号Ｓｒｑを発信して通信エリアＡを形成する。照合ＥＣＵ２１は、電子キー１がこの通信エリアＡに入り込んで返信してきたＩＤコード信号Ｓｉｄを通信機２４で受信すると、ＩＤコード信号Ｓｉｄを発信した電子キー１の位置を位置推定部２１ｂにより推定する。照合ＥＣＵ２１は、自身に登録されたＩＤコードと電子キー１のＩＤコードとを照らし合わせてＩＤ照合を行う。照合ＥＣＵ２１は、電子キー１が車内（Ａｉ）に位置すると推定して照合の成立（車内照合成立）を確認すると、エンジンスイッチ３５による電源状態の切替を許可する。
【００４３】
次に、照合ＥＣＵ２１の位置推定部２１ｂによる位置推定方法について図３〜図８を参照して説明する。本発明では、回転不変式が成立する仮想的なアレーへの変換を行う変換行列を用いることで、代数的手法であるＥＳＰＲＩＴ法の適用を可能とする。よって、定式化、ＥＳＰＲＩＴ法、近傍波源に対するＥＳＰＲＩＴ法の適用、初期推定誤差と窓関数の順番で説明する。なお、位置推定部２１ｂが、入力ベクトル算出手段、変換行列算出手段、及び波源位置算出手段として機能する。
【００４４】
＜定式化＞
波源の個数をＬとする。図３は（ｘ，y）＝（ｘｌ，ｙｌ）＝βｌに位置する第ｌ波（ｌ＝１，２，…，Ｌ）の狭帯域波源と、ｘ軸上に配置された素子間隔ｄのＫ素子等間隔リニアアレーとの位置関係を示した図である。Ｌ波の到来波がアレーアンテナに入射した場合、ｋ番目（ｋ＝１，２，…，Ｋ）のアンテナ素子の受信信号ｘｋは次式のように表すことができる。
【００４５】
【数１】
ここに、Ｆｌはｌ番目の波源の複素振幅であり、ｄｋ（βｌ）はｌ番目の波源とｋ番目のアレー素子との距離、ηｋ（βｌ）はｌ番目の波源の位置βｌからアレー素子位置までの伝搬による位相遅れである。また、ｗｋはｋ番目のアレー素子における熱雑音、λｌはｌ番目の波源の波長を表す。式（１）を用い、入力ベクトルＸを次のように定義する。
【００４６】
【数２】
このとき、振幅ベクトルＦ、アレーのモード行列Ａ、熱雑音ベクトルＷを次のようにそれぞれ定義する。
【００４７】
【数３】
なお、Ｔは転置を表す。このとき、アレーアンテナへの入力ベクトルは次のように表すことができる。
【００４８】
【数４】
このようにして定義された入力ベクトルＸの相関行列Ｒｘｘは、次式のように表すことができる。
【００４９】
【数５】
ここで、Ｈは複素共役転置を、Ｅ［・］は期待値を表す。各アンテナの熱雑音は互いに無相関であり、入射波とも無相関であるので、スナップショット数が無限である場合、熱雑音電力をσ２とすると、Ｗ＝σ２Ｉ（Ｉ：単位行列）となる。
【００５０】
＜ＥＳＰＲＩＴ法＞
ＥＳＰＲＩＴ法は、波源が十分遠方であり、到来波が平面波とみなせることを前提として成立する。今、十分遠方に存在するｌ番目（ｌ＝１，２，・・・，Ｌ）の波源が、ｙ軸となす角θl（反時計回り方向を正とする）から到来する場合を考える。すると、波源が近傍である場合に式（８）で表されたモードベクトルａ（βｌ）の成分ａｋ（βｌ）は、次式のように変更される。
【００５１】
【数６】
ここで、第１，２，…，Ｋ−１素子からなるサブアレーを第１サブアレー、第２，３，…，Ｋ素子からなるサブアレーを第２サブアレーと呼ぶこととする。このとき、Ｊ１，Ｊ２を次のように定義する。
【００５２】
【数７】
第１サブアレーのモード行列はＪ１Ａ、第２サブアレーのモード行列はＪ２Ａと表すことができる。すると、これらサブアレーのモード行列の間には次の関係式が成り立つ。
【００５３】
【数８】
ここに、diag［・］は対角行列を表す。式（１７）は回転不変式（ｒｏｔａｔｉｏｎａｌ ｉｎｖａｒｉａｎｃｅ）と呼ばれ、ＥＳＰＲＩＴ法の基本原理を表す式である。これは第１サブアレーの第ｌ波のモードベクトルＪ１ａ（βｌ）をexp｛ｊφｌ｝倍すると、第２サブアレーのモードベクトルＪ２ａ（βｌ）となることを意味し、サブアレー同士が同形かつ、波源が遠方であるからこそ成り立つものである。
【００５４】
さて、ここで式（１１）により得られた相関行列Ｒｘｘに対する固有値問題を解く。
【００５５】
【数９】
上記を解いて得られるＫ個の固有値を、値の大きいものから順にλ１，λ２，…，λＫとする。言い換えれば、式（２１）のΛの対角成分λｋが固有値である。このうち、大きいものからＬ個の固有値に対応する固有ベクトルｖ１，ｖ２，…，ｖＬの張るＬ次元部分空間は、行列Ａを構成するＬ個の列ベクトルの張る部分空間に一致する。したがって、以下を満たすＬ次の正則な行列Ｔが唯一存在する。
【００５６】
【数１０】
したがって、式（２３）をＡ＝ＥＳＴ−１と変形し、式（１７）へ代入すると次式が導かれる。ここで、式（２３）におけるＥＳが信号部分空間に属する固有ベクトルである。
【００５７】
【数１１】
ここで、ＥＸ，ＥＹ，Ψはそれぞれ以下のように定義した。
【００５８】
【数１２】
ＥＸ，ＥＹはそれぞれ式（２６），（２７）より計算可能であるため、式（２５）よりΨを算出可能である。よって、式（２８）より、Ψを固有値展開することでＴ−１，Φを求めることができ、Φが求まれば式（１８），（１９）より到来角度θｌを算出することができる。
【００５９】
＜近傍波源に対するＥＳＰＲＩＴ法の適用＞
上記では、波源が遠方である場合にはＥＳＰＲＩＴ法の適用が可能であることを確認した。しかし、波源が近傍である場合には、等間隔アレーであるにも関わらず、隣接素子間の受信信号の位相差が一定ではない不等間隔アレーで受信したような状態となる。したがって、ＥＳＰＲＩＴ法の原理式である回転不変式（１７）を満たすような、同形のサブアレーのペアを取り出すことができない。
【００６０】
そこで、本方法では、次式を満たす変換行列Ｂを使用する。
【００６１】
【数１３】
式（２９）は、「位置βｌ＝(ｘｌ，ｙｌ)を表すモードベクトルａ（βｌ）を、変換行列Ｂにより変換すると、ａ（βｌ）の（ｘｌ＋ｊｙｌ）倍となる」ことを意味する。そこで、このような変換行列Ｂが得られれば，以下の方法により波源位置を推定することが可能である。
【００６２】
まず、モード行列Ａを用いて式（２９）を次のように表す。
【００６３】
【数１４】
ここに、変換行列Ｂはモード行列ＡをＡΩという行列に変換することがわかる。この式（３０）を回転不変式と呼ぶ。また、ω（βｌ）が波源位置βｌと１対１対応する。次に、式（２３）をＡ＝ＥＳＴ−１と変形し、式（３０）へ代入すると，次式が導かれる。
【００６４】
【数１５】
ここで、ＥＸ，ＥＹ，Ψはそれぞれ以下のように定義した。
【００６５】
【数１６】
以上より、式（３３）は式（２５）と等しく、ＥＳＰＲＩＴ法が適用可能であることがわかる。また、式（３３）より得られるΨは、Ωが対角行列であることから、式（３６）のようにΩの対角成分を固有値、Ｔ−１の列ベクトルを固有ベクトルとする行列である。したがって、Ψを固有値展開することによりΩが得られ、更に式（３１），（３２）を用いることで波源位置βｌを求めることができる。すなわち、式（３０）を満たす変換行列Ｂを得ることができれば、近傍波源に対してＥＳＰＲＩＴ法を適用することができ、波源位置が推定可能となる。
【００６６】
式（２９）における変換は非線形変換であるため、全ての位置に対して式（２９）を満たす行列Ｂは存在しない。そこで、本方法では次式を用いて変換行列Ｂを算出する。
【００６７】
【数１７】
ここで、leftは左逆行列を表す。すなわち、ＡleftＡ＝Ｉである。モード行列ＡはＫ×Ｌ行列であり、素子数Ｋ＞波源数Ｌより縦長の行列である。また、Ａの列ベクトル、すなわち各波源の位置βｌに対するモードベクトルａ（βｌ）はそれぞれ一次独立であるため、左逆行列Ａleftが存在する。したがって、式（３７）により変換行列Ｂを算出することが可能である。モード行列Ａは波源位置βｌより算出されるものであり、推定前の段階では未知の行列だということである。したがって、まず初期推定を行い、おおよその波源位置からモード行列を算出し、これを用いて変換行列の計算を行うこととなる。初期推定の方法としては、ビームフォーマ法等が挙げられる。
【００６８】
＜初期推定誤差と窓関数＞
上記では、近傍波源に対するＥＳＰＲＩＴ法を適用する方法について説明した。続いて、変換行列Ｂを求めるために行う初期推定について説明する。また、初期推定結果には一般に誤差が含まれるため、この初期推定に含まれる誤差の影響と、その対策について説明する。
【００６９】
初期推定により得られた誤差を含む推定位置を図４のようにβｌ′＝（ｘｌ′、ｙｌ′）とする。このβｌ′から式（６），（７），（８）を用いて算出されるモード行列をＡ′とし、Ａ′から式（３７）を用いて算出される変換行列をＢ′とする。変換行列Ｂ′は、理想的な変換行列Ｂに対して誤差を有する。この誤差をＰとおくと誤差を含む変換行列Ｂ′は次式のように表される。
【００７０】
【数１８】
そして、Ｂ′により真の位置を表すモード行列Ａを変換する。
【００７１】
【数１９】
上記のように変換誤差ＰＡが発生し、モード行列を正しく変換することができない。変換誤差が存在する場合には式（３０）が成立しないため、次のような窓関数を用いてこれを対策する。なお、窓関数については、「アレー補完処理に窓関数を導入した空間スムージング型超解像側角法」信学論（Ｂ）、vol.J82-B,no.6,pp.1185-1192,June 1999を参考にする。
【００７２】
まず、窓関数Γ(βｌ′)を次のように定義する。
【００７３】
【数２０】
ここで、μを窓関数幅、Ｎを窓関数幅分割数、γｌp,q(βｌ′)を格子点と呼ぶこととする。格子点γｌp,q(βｌ′)は、位置βｌ′に対してｘ方向、ｙ方向にそれぞれ距離pμ／（Ｎ−１），qμ／（Ｎ−１）離れた位置を表す。
【００７４】
例えば、格子点γｌ-（N-1）/2,（N-1）/2(βｌ′)は、位置βｌ′に対して−ｘ方向にμ／２、＋ｙ方向にμ／２離れた位置を表し、格子点γｌ0,0(βｌ′)は位置βｌ′を表すことになる。つまり、窓関数Γ(βｌ′)は、図４のように位置βｌ′を中心とした幅μの正方形格子状の点の集合となる。このような窓関数Γ(βｌ′)を用いて、モード行列Ａ（Ａの上に・が付く）を次式のように定義する。
【００７５】
【数２１】
そして、モード行列Ａ（Ａの上に・が付く）を用いて変換行列Ｂ（Ｂの上に・が付く）を求める。ただし、モード行列Ａ（Ａの上に・が付く）はＫ×ＮＬ行列であるため、Ｋ＜ＮＬの場合には左逆行列Ａleft（Ａの上に・が付く）が存在せず、式（３７）により変換行列Ｂ（Ｂの上に・が付く）を計算することはできない。この場合には、次式を最小化するＢ（Ｂの上に・が付く）を変換行列とする。
【００７６】
【数２２】
ここで、||・||ＦはＦｒｏｂｅｎｉｕｓ ｎｏｒｍを表し、Ω（Ωの上に・が付く）は式（３１）においてβｌの代わりにΓ(βｌ′)を代入し算出されるＮＬ×ＮＬの対角行列である。このようにして求められた変換行列Ｂ（Ｂの上に・が付く）は、モード行列Ａ（Ａの上に・が付く）を行列ＡΩ（ＡとΩとの上に・がそれぞれ付く）に変換する際に誤差が生じてしまう。しかし、窓関数幅μが十分小さく、また窓関数幅分割数Ｎが十分大きければ、窓関数内では誤差は小さく、誤差は発生しないとみなすことができる。したがって、窓関数Γ(βｌ′)の内部に真の波源位置βｌが存在すれば、次式が成立する。
【００７７】
【数２３】
すなわち、初期推定位置βｌ′に誤差が含まれていても、窓関数Γ(βｌ′)を用いて求めた変換行列Ｂ（Ｂの上に・が付く）により、モード行列ＡをＡΩに誤差なく変換することができる。よって、上記で説明したように、ＥＳＰＲＩＴ法により波源位置を推定することが可能となる。
【００７８】
実際には、真の波源位置βlが窓関数内に含まれるように、窓関数幅μは初期推定誤差よりも大きくとる必要がある。この場合、式（３９）のように、窓関数幅μの大きさに応じて変換誤差ＰＡが発生する。しかし、変換誤差ＰＡが十分小さい場合には、変換誤差ＰＡに起因する推定誤差も小さく、ＥＳＰＲＩＴ法による推定位置βｌ′′は、初期推定位置βｌ′よりも真の波源位置βｌに近づくことが期待される。また、変換誤差ＰＡは窓関数幅μの縮小に応じて小さくなる。このため、推定位置βｌ′′を初期推定位置とし、窓関数幅μを縮小して繰り返しＥＳＰＲＩＴ法を適用することにより、さらに推定精度を改善することも可能である。
【００７９】
＜シミュレーション＞
計算機シミュレーションにより上記推定方法の評価を行う。本評価では波源数は既知とし、真の波源位置から最も近い推定位置までの距離を推定誤差として、ＲＭＳＥによる定量評価を実施する。ここでは、それぞれの真の波源位置から最も近い推定位置が重複する場合を推定失敗として定義し、位置推定に失敗した事象はＲＭＳＥの計算に含めないこととする。また、ＳＮＲを計算する際の信号強度については、アレー中心（＝原点）における第１波源の信号強度を用いている。
【００８０】
ここでは、ＥＳＰＲＩＴ法による推定位置を初期推定位置とし、窓関数幅μを縮小して再びＥＳＰＲＩＴ法を適用する「繰り返し推定」を行うこととし、窓関数幅μの初期値として２λを与える。波源数は２個とし、初期推定誤差としては各波源に対して＋ｘ方向、＋ｙ方向にそれぞれ√（１／２）λの誤差を与える。つまり、ビームフォーマ法等を用いた初期推定は行わないこととする。距離にすると１λの誤差を各波源に与えることになるが、このとき真の波源位置βｌは窓関数Γ(βｌ′)の内部に含まれることとなる。
【００８１】
＜シミュレーション結果＞
それでは、具体的にシミュレーション結果を確認する。図６は、ＳＮＲと推定誤差の関係をプロットしたものである。繰り返し推定の回数はパラメータとし、４回〜３２回まで変化させた。参考のためＭＵＳＩＣ法による推定結果も併せて示す。図６より、ＭＵＳＩＣ法よりは若干劣るが、十分な繰り返しを行えば本方法は高い推定精度が得られることがわかる。また、本シミュレーション条件では、繰り返し回数は１２回程度で十分であることもわかる。しかし、図６は、シミュレーション前に初期窓関数幅、窓関数縮小率の最適化を行った上で得られた結果である。波源位置、アレー形状、初期推定誤差等が変化すれば、最適な初期窓関数幅、窓関数縮小率も変化することが予想される。
【００８２】
そこで、初期窓関数幅により推定結果がどの様に変化するかを確認する。図７は初期窓関数幅を２λから徐々に小さくした場合の、ある一回の独立試行における繰り返し回数と推定誤差の関係をプロットしたものである。なお、熱雑音による推定誤差が発生しないよう、熱雑音はないものとしてシミュレーションを行った。図７より、初期窓関数幅μ＝２λ，１．６λの場合には、繰り返し回数に応じて推定誤差が減少しており、波源位置を正しく推定できることがわかる。しかし、初期窓関数幅μ＝１．２λ，０．８λの場合には、繰り返しの途中で推定誤差が急激に増大しており、波源位置を推定できていないことが確認される。これは、初期窓関数幅μ＝１．２λ，０．８λの場合には、初期推定誤差に対して初期窓関数幅が小さく、初期窓関数内に真の波源位置が含まれない状態となっているためと考えられる。このような状態では、変換行列Ｂ（Ｂの上に・が付く）によるモード行列Ａの変換が正しく行われない、すなわち式（４５）が成立しなくなるため、正しく推定することができない。したがって初期窓関数幅は、初期推定誤差の大きさに応じてある程度広くとる必要がある。図８は、図７とは逆に、初期窓関数幅を１λから徐々に大きくした場合の、ある一回の独立試行における繰り返し回数と推定誤差の関係をプロットしたものである。ここでも、先程と同様に熱雑音はないものとした。また、初期推定誤差についても存在しない状態、つまり初期推定位置として真の波源位置を与えた状態でシミュレーションを行った。これは、変換行列Ｂ（Ｂの上に・が付く）によりモード行列Ａを変換する際に発生する誤差だけが、推定誤差の要因となる状態でシミュレーションを行うためである。したがって、繰り返し回数が１回目のときには推定誤差が０であるため、グラフ上に値がプロットされていない。繰り返し回数が２回目のときには推定誤差が発生しているが、これは変換誤差に起因するものである。図８より、初期窓関数幅μ＝１λ，２λの場合には、波源位置を正しく推定できることがわかる。
【００８３】
この結果は、初期窓関数幅μの大きさに上限があることを意味している。先程も述べたとおり、真の波源位置が窓関数内に含まれるためには、初期推定に含まれる誤差は、初期窓関数幅の上限より小さくなければならない。
【００８４】
さて、近傍波源に対する位置推定方法として、回転不変式の成立する仮想的アレーへの変換を行う変換行列を用いて、ＥＳＰＲＩＴ法を適用可能とする方法を説明した。シミュレーションでは７素子等間隔リニアアレーを用い、近傍に位置する２つの波源に対して本方法を適用し、ＭＵＳＩＣ法より若干劣る程度の高い推定精度が得られることを明らかにした。近傍波源をサーチすることなく推定でき、計算負荷を低減することができる。しかもアレー形状に対する制約が無いということは、非常に有用である。
【００８５】
以上、説明した実施形態によれば、以下の作用効果を奏することができる。
（１）アレーアンテナにより電波を受信して波源の位置を算出するにあたって、各アレーが受信した受信波の入力ベクトルＸを算出し、代数的手法で近傍波源の位置特定を可能とする変換行列Ｂを算出し、入力ベクトルＸと変換行列Ｂとを基に代数的手法としてＥＳＰＲＩＴ法における代数的処理によって波源の位置を算出する。なお、入力ベクトルＸは、自身の相関行列Ｒｘｘの固有ベクトルＶの中で信号部分空間に属する固有ベクトルＥＳとして用いる。従来、近傍波源の位置を推定する場合にはＭＵＳＩＣ法等のサーチ法が用いられるが、サーチに要する計算負荷が大きいという問題点がある。また、代数的手法により近傍波源の位置を推定する場合には、受信信号を近似する必要があり、近似に起因する誤差が発生する。また、アレー形状も等間隔リニアアレーかつ素子間隔が１／４波長以下である必要がある。そこで、本装置では、近傍波源の受信波の場合に、変換行列Ｂを使用することで代数的手法の適用を可能とするので、計算負荷を低減できる。また、演算に際して入力値を近似することもないので、高い位置推定精度を確保することができる。
【００８６】
（２）固有値は、式（３１）における対角行列Ωの対角成分ω（βｌ）であって、ω（βｌ）は式（３２）に示されるように、波源位置と１対１対応する。波源位置を固有値から算出するため、簡単な代数計算を行えばよく、計算負荷を低減できる。
【００８７】
（３）固有値ω（βｌ）の実部が波源位置βｌのｘ座標、虚部がｙ座標を表すため、固有値から波源の位置を算出する際に、さらに簡単に計算することができる。
（４）回転不変式（３０）の中の対角行列Ωの成分ω（βｌ）が固有値に相当する。この対角行列Ωは、式（３６）のΨを固有値展開することにより得られる。すなわち、代数的手法で用いられる回転不変式、すなわちアレー全体の平行移動によって生じる各波の位相回転、及び強度の変化が成立する式の対角行列の成分として算出された固有値から波源位置を算出する際に、簡単な代数演算を行えばよいため、計算負荷が軽減される。
【００８８】
（５）入力ベクトルＸの自己相関行列Ｒｘｘを固有値展開して得られる信号部分空間に属する固有ベクトルＥＳを変換行列Ｂにより変換して波源位置を算出する。このため、固有値展開の演算量を削減することができる。
【００８９】
（６）窓関数を用いて変換行列を算出するため、推定領域全体を変換対象として変換行列を算出した場合と比較し、変換の際に発生する変換誤差を抑制することができる。
（７）窓関数の位置を初期推定に基づいて算出するため、初期推定を行わずに窓関数の位置を決定した場合と比較し、発生する変換誤差をより抑制することができる。
【００９０】
（８）波源位置推定を繰り返し行うため、位置の推定精度を更に向上させることができる。
（９）波源位置推定を繰り返す際に推定範囲を狭くするため、繰り返す度に位置の推定精度を更に向上させることができる。
【００９１】
（第２の実施形態）
以下、本発明にかかる波源位置推定装置及び波源位置推定方法を車両の電子キーシステムに適用した第２の実施形態について説明する。この実施形態では、代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値を考慮することなく変換行列Ｂを算出する点、変換行列Ｂを入力ベクトルＸに掛けて変換する点、及びモードベクトルによって波源の位置を算出する点が上記第１の実施形態と異なっている。以下、第１の実施形態との相違点を中心に説明する。
【００９２】
まず、＜ＥＳＰＲＩＴ法＞まで第１の実施形態と同様である。
次に、近傍波源に対するＥＳＰＲＩＴ法の適用について説明する。
本方法では、次式を満たす変換行列Ｂを使用する。
【００９３】
【数２４】
式（４６）は、「位置βｌ＝(ｘｌ，ｙｌ)を表すモードベクトルａ（βｌ）を、変換行列Ｂにより変換すると、ａ（βｌ）のω（βｌ）倍となる」ことを意味する。そこで、このような変換行列Ｂが得られれば，以下の方法により波源位置を推定することが可能である。
【００９４】
まず、モード行列Ａを用いて式（４６）を次のように表す。
【００９５】
【数２５】
ここに、変換行列Ｂはモード行列ＡをＡΩという行列に変換することがわかる。この式（４８）を回転不変式と呼ぶ。
【００９６】
次に、この変換行列Ｂを用いて、次式で表される入力ベクトルＸＢを定義する。
【００９７】
【数２６】
このように定義された入力ベクトルＸＢの相関行列ＲｘｘＢは次式のように表される。
【００９８】
【数２７】
この相関行列ＲｘｘＢを、Ｂ２を用いて次のように一般固有値展開する。
【００９９】
【数２８】
上記を解いて得られるＫ個の固有値を、値の大きいものから順にλ１，λ２，…，λＫとする。言い換えれば、式（５４）のΛの対角成分λＫが固有値である。このうち、大きいものからＬ個の固有値に対応する固有ベクトルｖ１，ｖ２，…，ｖＬからなる行列をＥＳと定義する。このとき、ＥＳとモード行列Ａの間には、次の関係式が成り立つ。
【０１００】
【数２９】
次に、式（５６）をＡ＝ＢＨＥＳＴ−１と変形し、式（４８）に代入すると、次式が導かれる。
【０１０１】
【数３０】
ここで、ＥＸ，ＥＹ，Ψはそれぞれ以下のように定義した。
【０１０２】
【数３１】
以上より、式（５７）は式（２５）と等しく、ＥＳＰＲＩＴ法が適用可能であることがわかる。また、式（５７）より得られるΨは、Ωが対角行列であることから、式（６０）のようにΩの対角成分を固有値、Ｔ−１の列ベクトルを固有ベクトルとする行列である。したがって、Ψを固有値展開することによりＴが得られる。このＴと、信号部分空間に属する固有ベクトルＥＳを式（２３）に代入することにより、モード行列Ａを求めることができる。よって、モード行列Ａの各列ベクトルであるモードベクトルａ（βｌ）から、式（６）、（７）、（８）を用いて波源位置βｌを算出することが可能である。すなわち、式（４８）を満たす変換行列Ｂを得ることができれば、近傍波源に対してＥＳＰＲＩＴ法を適用することができ、波源位置が推定可能となる。以下、第１の実施形態と同様であるため割愛する。
【０１０３】
さて、本実施形態では、代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値を考慮することなく変換行列Ｂを算出し、変換行列Ｂを入力ベクトルＸに掛けて変換し、モードベクトルから波源の位置を算出した。このように第１の実施形態と異なる方法にしても近傍の波源位置を算出することができ、変換行列を使用することで代数的手法の適用を可能とするので、任意のアレー形状において計算負荷を低減できる。また、演算に際して入力値を近似することもないので、高い位置推定精度を確保することができる。
【０１０４】
以上、説明した実施形態によれば、第１の実施形態の（１）及び（６）〜（９）の作用効果に加え、以下の作用効果を奏することができる。
（１０）代数的手法で用いられる回転不変式の対角行列の対角成分として算出される固有値を考慮することなく変換行列を算出するため、変換行列Ｂの計算に用いる固有値の値を、入力ベクトルから算出する必要がなく、この部分の計算負荷が軽減される。
【０１０５】
（１１）代数的手法で算出されたモードベクトルから波源位置を算出するため、変換行列Ｂを計算する際に固有値について考慮する必要がなく、変換行列算出のための計算負荷が軽減される。
【０１０６】
（１２）入力ベクトルＸを変換行列Ｂにより変換して波源位置を算出する。このため、入力ベクトルＸの自己相関行列Ｒｘｘに対する固有値展開において、一般固有値展開を適用することができる。したがって、変換行列Ｂがユニタリ行列でない場合においても、雑音による推定精度の劣化を最小限に留めることができる。
【０１０７】
（第３の実施形態）
以下、本発明にかかる波源位置推定装置及び波源位置推定方法を車両の電子キーシステムに適用した第３の実施形態について説明する。この実施形態では、変換行列を入力ベクトルのうち一部のみを変換することが可能な行列とした点が上記第１の実施形態と異なっている。以下、第１の実施形態との相違点を中心に説明する。
【０１０８】
まず、＜ＥＳＰＲＩＴ法＞まで第１の実施形態と同様である。
次に、近傍波源に対するＥＳＰＲＩＴ法の適用について説明する。
本方法では、次式を満たす変換行列Ｂを使用する。
【０１０９】
【数３２】
式（６１）は、「位置βｌ＝（ｘｌ，ｙｌ）を表す第２サブアレーのモードベクトルＪ２ａ（βｌ）を変換行列Ｂにより変換すると、第１サブアレーのモードベクトルＪ１ａ（βｌ）と同一形状となり、その際に第２サブアレーの基準素子位置は変化させない」ことを意味する。そこで、このような変換行列Ｂが得られれば，以下の方法により波源位置を推定することが可能である。
【０１１０】
まず、モード行列Ａを用いて式（６１）を次のように表す。
【０１１１】
【数３３】
ここに、変換行列Ｂは第２サブアレーのモード行列Ｊ２Ａを第１サブアレーのモード行列Ｊ１ＡにΩを乗じた行列に変換することがわかる。この式（６２）を回転不変式と呼ぶ。本実施形態においては第１の実施形態と異なり、変換行列Ｂはモード行列Ａのうち、第２、３、・・・、Ｋ素子に対応するモード行列Ｊ２Ａの変換を行う。すなわち、モード行列Ａの一部のみを変換する。このような変換行列Ｂを用いると、ＢＪ２Ｘ＝ＢＪ２ＡＦ＝Ｊ１ＡΩＦのように入力ベクトルの一部のみを変換することが可能である。次に、式（２３）をＡ＝ＥＳＴ−１と変形し、式（６２）へ代入すると，次式が導かれる。
【０１１２】
【数３４】
ここで、ＥＸ，ＥＹ，Ψはそれぞれ以下のように定義した。
【０１１３】
【数３５】
以上より、式（６５）は式（２５）と等しく、ＥＳＰＲＩＴ法が適用可能であることがわかる。また、式（３３）より得られるΨは、Ωが対角行列であることから、式（６８）のようにΩの対角成分を固有値、Ｔ−１の列ベクトルを固有ベクトルとする行列である。したがって、Ψを固有値展開することによりΩが得られ、更に式（６３），（６４）を用いることで波源位置βｌを求めることができる。すなわち、式（６２）を満たす変換行列Ｂを得ることができれば、任意形状のアレーに対してＥＳＰＲＩＴ法を適用することができ、波源位置が推定可能となる。
【０１１４】
本方法では次式を用いて変換行列Ｂを算出する。
【０１１５】
【数３６】
ここで、leftは左逆行列を表す。すなわち、｛Ｊ２Ａ｝left｛Ｊ２Ａ｝＝Ｉである。第２サブアレーのモード行列｛Ｊ２Ａ｝は（Ｋ−１）×Ｌ行列であり、（Ｋ−１）≧Ｌより縦長の行列、あるいは正方行列である。また、Ｊ２Ａの列ベクトル、すなわち各波源の位置βｌに対する第２サブアレーのモードベクトルＪ２ａ（βｌ）はそれぞれ一次独立であるため、左逆行列｛Ｊ２Ａ｝leftが存在する。したがって、式（６９）により変換行列Ｂを算出することが可能である。モード行列Ａは波源位置βｌより算出されるものであり、推定前の段階では未知の行列だということである。したがって、まず初期推定を行い、おおよその波源位置からモード行列を算出し、これを用いて変換行列の計算を行うこととなる。初期推定の方法としては、ビームフォーマ法等が挙げられる。
【０１１６】
以下、第１の実施形態と同様であるため割愛する。
さて、本実施形態では、変換行列を入力ベクトルのうち一部のみを変換することが可能な行列とした。このように、変換行列は入力ベクトルのうち、一部のみを変換することが可能である。したがって、変換行列の要素数が少なく、変換行列の算出に掛かる計算負荷を低減することができる。
【０１１７】
以上、説明した実施形態によれば、第１の実施形態の（１）、（２）及び（４）〜（９）の作用効果に加え、以下の作用効果を奏することができる。
（１３）変換行列は、入力ベクトルのうち一部のみを変換することが可能な行列とした。したがって、変換行列の要素数が少なく、変換行列の算出に掛かる計算負荷を低減することができる。
【０１１８】
（第４の実施形態）
以下、本発明にかかる波源方向推定装置及び波源方向推定方法を車両の電子キーシステムに適用した第４の実施形態について説明する。この実施形態では、位置推定ではなく、方向推定である点が上記第１の実施形態と異なっている。以下、第１の実施形態との相違点を中心に説明する。なお、位置推定部２１ｂが、入力ベクトル算出手段、変換行列算出手段、及び波源方向算出手段として機能する。
【０１１９】
まず、＜ＥＳＰＲＩＴ法＞まで第１の実施形態と同様である。
次に、任意形状アレーに対するＥＳＰＲＩＴ法の適用について説明する。
＜ＥＳＰＲＩＴ法＞では、アレーアンテナの形状が等間隔リニアアレーである場合にＥＳＰＲＩＴ法の適用が可能であることを確認した。しかし、アレーアンテナの形状が任意である場合には、隣接素子間の受信信号の位相差が一定でないため、ＥＳＰＲＩＴ法の原理式である回転不変式（１７）を満たすような、同形のサブアレーのペアを取り出すことができない。
【０１２０】
そこで、本方法では、次式を満たす変換行列Ｂを使用する。
【０１２１】
【数３７】
式（７１）は、「到来波の方向θｌを表すモードベクトルａ（βｌ）を変換行列Ｂにより変換すると、ａ（βｌ）の（ｅｘp｛ｊθｌ｝）倍となる」ことを意味する。そこで、このような変換行列Ｂが得られれば，以下の方法により波源方向を推定することが可能である。
【０１２２】
まず、モード行列Ａを用いて式（６１）を次のように表す。
【０１２３】
【数３８】
ここに、変換行列Ｂはモード行列ＡをＡΩという行列に変換することがわかる。この式（７２）を回転不変式と呼ぶ。また、ω（βｌ）が波源方向と１対１対応する。次に、式（２３）をＡ＝ＥＳＴ−１と変形し、式（７２）へ代入すると，次式が導かれる。
【０１２４】
【数３９】
ここで、ＥＸ，ＥＹ，Ψはそれぞれ以下のように定義した。
【０１２５】
【数４０】
以上より、式（７５）は式（２５）と等しく、ＥＳＰＲＩＴ法が適用可能であることがわかる。また、式（３３）より得られるΨは、Ωが対角行列であることから、式（７８）のようにΩの対角成分を固有値、Ｔ−１の列ベクトルを固有ベクトルとする行列である。したがって、Ψを固有値展開することによりΩが得られ、更に式（７３），（７４）を用いることで波源方向θｌを求めることができる。すなわち、式（７２）を満たす変換行列Ｂを得ることができれば、任意形状のアレーに対してＥＳＰＲＩＴ法を適用することができ、波源方向が推定可能となる。以下、第１の実施形態の波源位置βｌを波源方向θｌと置き換えれば同様であるため割愛する。
【０１２６】
さて、本実施形態では、等間隔リニアアレーではなく任意形状のアレーにおいて、近傍波源の方向を算出した。このように、任意形状のアレーにおいても波源の方向を算出することができ、変換行列を使用することで代数的手法の適用を可能とするので、計算負荷を低減することができる。また、演算に際して入力値を近似することもないので、高い方向推定精度を確保することができる。
【０１２７】
以上、説明した実施形態によれば、以下の作用効果を奏することができる。
（１）アレーアンテナにより電波を受信して波源の方向を算出するにあたって、各アレーが受信した受信波の入力ベクトルＸを算出し、代数的手法で近傍波源の方向特定を可能とする変換行列Ｂを算出し、入力ベクトルＸと変換行列Ｂとを基に代数的手法としてＥＳＰＲＩＴ法における代数的処理によって波源の方向を算出する。なお、入力ベクトルＸは、自身の相関行列Ｒｘｘの固有ベクトルＶの中で信号部分空間に属する固有ベクトルＥＳとして用いる。従来、到来波の方向を代数的手法により推定する場合には、アレー形状が、例えば等間隔リニアアレー等に制約される。しかしながら、本装置では、任意形状のアレーで受信した場合においても、変換行列を用いることにより代数的手法の適用を可能とするため、計算負荷を低減できる。
【０１２８】
（２）固有値は、式（７２）における対角行列Ωの対角成分ω（βｌ）であって、ω（βｌ）は式（７４）に示されるように、波源方向と１対１対応する。波源方向を固有値から算出するため、簡単な代数計算を行えばよく、計算負荷を低減できる。
【０１２９】
（３）固有値ω（βｌ）の位相成分が波源方向θｌを表すため、固有値から波源の方向を算出する際に、さらに簡単に計算することができる。
（４）回転不変式（７２）の中の対角行列Ωの成分ω（βｌ）が固有値に相当する。この対角行列Ωは、式（７８）のΨを固有値展開することにより得られる。すなわち、代数的手法で用いられる回転不変式、すなわちアレー全体の平行移動によって生じる各波の位相回転、及び強度の変化が成立する式の対角行列の成分として算出された固有値から波源方向を算出する際に、簡単な代数演算を行えばよいため、計算負荷が軽減される。
【０１３０】
（５）入力ベクトルＸの自己相関行列Ｒｘｘを固有値展開して得られる信号部分空間に属する固有ベクトルＥｓを変換行列Ｂにより変換して波源方向を算出する。このため、固有値展開の演算量を削減することができる。
【０１３１】
（６）窓関数を用いて変換行列を算出するため、推定領域全体を変換対象として変換行列を算出した場合と比較し、変換の際に発生する変換誤差を抑制することができる。
（７）窓関数の位置を初期推定に基づいて算出するため、初期推定を行わずに窓関数の位置を決定した場合と比較し、発生する変換誤差をより抑制することができる。
【０１３２】
（８）波源方向推定を繰り返し行うため、方向の推定精度を更に向上させることができる。
（９）波源方向推定を繰り返す際に推定範囲を狭くするため、繰り返す度に方向の推定精度を更に向上させることができる。
【０１３３】
なお、上記実施形態は、これを適宜変更した以下の形態にて実施することができる。
・上記第１、第３、及び第４の実施形態では、変換行列を、入力ベクトルすべてを変換可能な行列としたが、入力ベクトルのうち一部のみを変換することが可能な行列としてもよい。入力ベクトルのうち一部のみを変換することが可能な変換行列は、行列の要素数が少ない。したがって、変換行列の算出に掛かる計算負荷を低減することができる。
【０１３４】
・上記第２の実施形態では、変換行列を入力ベクトルのうち一部のみを変換可能な行列としたが、入力ベクトルのすべてを変換することが可能な行列としてもよい。
・上記第１、第３、及び第４の実施形態では、回転不変式の対角行列の成分として算出される固有値と、波源の存在し得る位置とが、一対一対応となる変換行列を算出したが、回転不変式の対角行列の対角成分として算出される固有値を考慮することなく変換行列を算出してもよい。
【０１３５】
・上記第１、第３、及び第４の実施形態では、位置の算出において固有値を使用したが、モードベクトルを使用してもよい。すなわち、式（３３）により得られるΨは、式（３６）となる。そして、このΨを固有値展開し、固有値Ω、固有ベクトルＴ−１を得る。次に、入力ベクトルＸの自己相関行列Ｒｘｘの固有ベクトルＶの中で信号部分空間に属する固有ベクトルＥＳと、Ψの固有ベクトルＴ−１とから式（２３）を用いてモード行列を算出する。そして、モード行列Ａの各列成分であるモードベクトルａ（βｌ）から、波源位置βｌを例えば非線形計画法により算出する。
【０１３６】
・上記第２の実施形態では、回転不変式の対角行列の対角成分として算出される固有値を考慮することなく変換行列を算出したが、回転不変式の対角行列の成分として算出される固有値と、波源の存在し得る位置とが、一対一対応となる変換行列を算出してもよい。
【０１３７】
・上記第１の実施形態では、変換行列は式（２９）の変換を行う行列としたが、式（２９）以外の変換を行う行列を変換行列としてもよい。
・上記第２の実施形態では、変換行列は式（４６）の変換を行う行列としたが、式（４６）以外の変換を行う行列を変換行列としてもよい。
【０１３８】
・上記第３の実施形態では、変換行列は式（６１）の変換を行う行列としたが、式（６１）以外の変換を行う行列を変換行列としてもよい。
・上記第１の実施形態では、変換行列は式（７１）の変換を行う行列としたが、式（７１）以外の変換を行う行列を変換行列としてもよい。
【０１３９】
・上記第１、第３、及び第４の実施形態では、変換行列は式（３７）により変換行列を算出したが、式（３７）以外の式を用いて変換行列を算出してもよい。
・上記第３の実施形態では、変換行列は式（６９）により変換行列を算出したが、式（６９）以外の式を用いて変換行列を算出してもよい。
【０１４０】
・上記第１、第３、及び第４の実施形態では、変換行列Ｂと固有ベクトルＥｓとを掛けたが、変換行列Ｂと入力ベクトルＸとを掛けて算出するようにしてもよい。
・上記第２の実施形態では、変換行列Ｂと入力ベクトルＸとを掛けたが、変換行列Ｂと固有ベクトルＥｓとを掛けて算出するようにしてもよい。
【０１４１】
・上記構成において、アレー配置を円形としてもよい。このようにすれば、窓関数内部における変換誤差の不均一性を抑制することができるため、位置又は方向の推定精度を更に向上させることができる。
【０１４２】
・上記第１〜３の実施形態では、アレー配置は等間隔リニアアレーとしたが、これに限らず任意のアレー配置としてよい。
・上記実施形態では、変換行列算出時に窓関数を用いたが、所望の推定精度を得ることができれば、窓関数を用いなくともよい。
【０１４３】
・上記実施形態では、所望の推定精度を得ることができれば、繰り返し時に推定範囲を狭くしなくてもよい。
・上記実施形態では、代数的手法としてＥＳＰＲＩＴ法を用いたが、他の代数的手法、たとえばＤＯＡ−Ｍａｔｒｉｘ法等を使用してもよい。
【０１４４】
・上記実施形態では、波源位置推定の算出処理を１２回程度繰り返し行うのがよいとしたが、所望の推定精度に合わせて適宜変更可能である。また、１回で所望の推定精度を得ることができれば、繰り返さなくてもよい。
【０１４５】
・上記実施形態では、波源位置推定装置を電子キーシステムの電子キー１の位置推定に適用したが、電子キーシステムに限らず、種々の通信機器に適用してもよい。
・また、波源位置推定手段を備えた波源位置推定装置に限らず、波源方向推定手段を備えた波源方向推定装置として種々の通信機器に適用してもよい。
【０１４６】
・波源位置推定装置及び波源方向推定装置等の装置に限らず、波源位置推定方法及び波源方向推定方法を、アレーアンテナを備えた種々の通信機器の演算に適用してもよい。
・上記構成において、図１に示されるように、電子キー１に施錠スイッチ１９ａ及び解錠スイッチ１９ｂを設け、この施錠スイッチ１９ａ及び解錠スイッチ１９ｂを操作することで、車両２へ施錠信号Ｓｌ又は解錠信号Ｓｕｌを発信して、車両２のドアロックの解錠又は施錠を可能としてもよい。
【０１４７】
・電子キーシステムで使用する各通信の周波数はどの周波数の電波を使用してもよい。
・上記実施形態において、波源は電波を発信する装置であったが、これに限らず、例えばスピーカーのような音波を発信する装置であっても良い。
【０１４８】
・上記実施形態において、受信素子は電波を受信する装置であったが、これに限らず、例えばマイクロホンのような音波を受信する装置であっても良い。
・上記実施形態では、電波に適用したが、音波に適用してもよい。
【０１４９】
次に、上記実施形態から把握できる技術的思想をその効果と共に記載する。
（イ）請求項２に記載の波源方向推定装置において、前記変換行列算出手段は、前記代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値と、波源の存在し得る方向とが、一対一対応となる前記変換行列を算出することを特徴とする波源方向推定装置。
【０１５０】
同構成によれば、波源位置を固有値から算出するため、簡単な代数計算を行えばよく、計算負荷を低減できる。
（ロ）（イ）に記載の波源方向推定装置において、前記固有値は、前記波源の方向を表す複素数であることを特徴とする波源方向推定装置。
【０１５１】
同構成によれば、固有値の位相成分が波源方向を表すため、固有値から波源の方向を算出する際に、さらに簡単に計算することができる。
（ハ）請求項２、又は（イ）に記載の波源方向推定装置において、前記変換行列は、前記入力ベクトルのうち、一部のみを変換することが可能な行列であることを特徴とする波源方向推定装置。
【０１５２】
同構成によれば、変換行列は入力ベクトルのうち、一部のみを変換することが可能である。したがって、変換行列の要素数が少なく、変換行列の算出に掛かる計算負荷を低減することができる。
【０１５３】
（ニ）請求項２に記載の波源方向推定装置において、前記変換行列算出手段は、前記代数的手法で用いられる回転不変式の対角行列の対角成分として算出される固有値を考慮することなく前記変換行列を算出することを特徴とする波源方向推定装置。
【０１５４】
同構成によれば、変換行列の計算に用いる固有値の値を、入力ベクトルから算出する必要がないため、この部分の計算負荷が軽減される。
（ホ）請求項２、及び（イ）〜（ハ）のいずれかに記載の波源方向推定装置において、前記波源方向算出手段は、前記代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値によって前記波源の方向を算出することを特徴とする波源方向推定装置。
【０１５５】
同構成によれば、代数的手法で用いられる回転不変式、すなわち受信素子全体の平行移動によって生じる各波の位相回転、及び強度の変化が成立する式の対角行列の成分として算出された固有値から波源方向を算出する際に、簡単な代数演算を行えばよいため、計算負荷が軽減される。
【０１５６】
（へ）請求項２、及び（イ）〜（ニ）のいずれかに記載の波源方向推定装置において、前記波源方向算出手段は、前記代数的手法で得られるモードベクトルによって前記波源の方向を算出することを特徴とする波源方向推定装置。
【０１５７】
同構成によれば、代数的手法で算出されたモードベクトルから波源方向を算出するため、変換行列を計算する際に固有値について考慮する必要がなく、変換行列算出のための計算負荷が軽減される。
【０１５８】
（ト）請求項２、及び（イ）〜（へ）のいずれかに記載の波源方向推定装置において、前記波源方向算出手段は、前記入力ベクトルの自己相関行列を固有値展開して得られる信号部分空間に属する固有ベクトルを前記変換行列により変換して波源方向を算出することを特徴とする波源方向推定装置。
【０１５９】
同構成によれば、波源方向算出手段が行う固有値展開の演算量を削減することができる。
（チ）請求項２、及び（イ）〜（へ）のいずれかに記載の波源方向推定装置において、前記波源方向算出手段は、前記入力ベクトルを前記変換行列により変換して波源方向を算出することを特徴とする波源方向推定装置。
【０１６０】
同構成によれば、波源方向算出手段で行われる入力ベクトルの自己相関行列に対する固有値展開において、一般固有値展開を適用することができる。したがって、変換行列がユニタリ行列でない場合においても、雑音による推定精度の劣化を最小限に留めることができる。
【０１６１】
（リ）請求項２、及び（イ）〜（チ）のいずれかに記載の波源方向推定装置において、前記変換行列算出手段は、窓関数を用いて前記変換行列を算出することを特徴とする波源方向推定装置。
【０１６２】
同構成によれば、窓関数を用いて変換行列を算出するため、推定領域全体を変換対象として変換行列を算出した場合と比較し、変換の際に発生する変換誤差を抑制することができる。
【０１６３】
（ヌ）（リ）に記載の波源方向推定装置において、前記窓関数の位置を初期推定に基づいて算出することを特徴とする波源方向推定装置。
同構成によれば、窓関数の位置を初期推定に基づいて算出するため、初期推定を行わずに窓関数の位置を決定した場合と比較し、発生する変換誤差をより抑制することができる。
【０１６４】
（ル）請求項２、及び（イ）〜（ヌ）のいずれかに記載の波源方向推定装置において、前記波源方向推定を繰り返し行うことを特徴とする波源方向推定装置。
同構成によれば、波源方向推定を繰り返し行うため、方向の推定精度を更に向上させることができる。
【０１６５】
（ヲ）（ル）に記載の波源方向推定装置において、前記波源方向推定を繰り返す際に、推定範囲を狭くすることを特徴とする波源方向推定装置。
同構成によれば、波源方向推定を繰り返す際に推定範囲を狭くするため、繰り返す度に方向の推定精度を更に向上させることができる。
【０１６６】
（ワ）請求項２、及び（イ）〜（ヲ）のいずれかに記載の波源方向推定装置において、前記受信素子は、円形に配置されることを特徴とする波源方向推定装置。
同構成によれば、窓関数内部における変換誤差の不均一性を抑制することができるため、方向の推定精度を更に向上させることができる。
【符号の説明】
【０１６７】
１…電子キー、２…車両、１１…通信制御部、１１ａ…メモリ、１２…受信部、１３…発信部、２１…照合ＥＣＵ、２１ａ…メモリ、２１ｂ…位置推定部、２４…通信機、２４ａ…アレーアンテナ、３５…エンジンスイッチ、３８…ドアロック装置、Ａ…通信エリア、Ａｉ…車内エリア、Ｓｉｄ…ＩＤコード信号、Ｓｒｑ…リクエスト信号。
【技術分野】
【０００１】
この発明は、近傍波源である受信波の波源位置推定装置及び波源方向推定装置、並びに受信波の波源位置推定方法及び波源方向推定方法に関する。
【背景技術】
【０００２】
近年、地上デジタル放送機器、携帯電話、電子キー等の様々な無線機器の使用が日常化しており、希望の相手の位置を把握するため、また希望の相手からの電波を受信するために伝播特性測定の必要性が高まっている。
【０００３】
例えば、車両においては、無線通信により車両に固有のＩＤコードを発信可能な電子キーを車両のキーとして使用する電子キーシステムが搭載されている。この電子キーシステムでは、電子キーから発信されたＩＤコードを車両が受信すると、車両が電子キーのＩＤコードと自身に予め登録されているＩＤコードとを照合するＩＤ照合を行い、これらＩＤコードが一致してＩＤ照合が成立すれば、車載機器の制御、例えば車両のドアロックの施解錠やエンジンの始動が許可又は実行される。
【０００４】
この電子キーシステムにおいて、電子キーの位置が車内と車外とのいずれにあるかで車載機器の制御を切り分けている場合には、車両からのＩＤコードの返信要求を車内と車外とで異なるタイミングで発信して対応している。そこで、車両が電子キーから発信された電波の波源の位置を推定できれば、ＩＤコードの返信要求を異なるタイミングで発信する動作を車両に実行させる必要がなくなるので、この利点を享受するために、波源位置推定によりキー位置を特定する技術が求められている。
【０００５】
波源位置推定を行う方法としては、ＭＵＳＩＣ（ＭＵｌｔｉｐｌｅ ＳＩｇｎａｌ Ｃｌａｓｓｉｆｉｃａｔｉｏｎ）法やＥＳＰＲＩＴ（Ｅｓｔｉｍａｔｉｏｎ ｏｆ Ｓｉｇｎａｌ Ｐａｒａｍｅｔｅｒｓ ｖｉｓａ Ｒｏｔａｔｉｏｎａｌ Ｉｎｖａｒｉａｎｃｅ Ｔｅｃｈｎｉｑｕｅｓ）法等が挙げられるが、これらの方法では波源が十分遠方であり、到来波が平面波として近似できることを前提としている（例えば、非特許文献１参照）。
【０００６】
そこで、モードベクトルを近傍波源に対応したものに変更することで近傍波源位置推定を可能とした、球面波モードベクトルＭＵＳＩＣ法が提案されている。しかし、ＭＵＳＩＣ法を近傍波源に対して適用する場合、２次元のピークサーチが必要となるため、その計算負荷が懸念される。
【先行技術文献】
【非特許文献】
【０００７】
【非特許文献１】菊間信良著「アレーアンテナによる適応信号処理」科学技術出版
【発明の概要】
【発明が解決しようとする課題】
【０００８】
一方、ＥＳＰＲＩＴ法に基づいた近傍波源位置推定法も提案されている。これは、ＥＳＰＲＩＴ法に類似したＤＯＡ（Ｄｉｒｅｃｔｉｏｎ Ｏｆ Ａｒｒｉｖａｌ）−Ｍａｔｒｉｘ法を近傍波源に対して適用可能としたものであり、代数的に位置を算出することができる。したがって、ピークサーチが不要であり、球面波モードベクトルＭＵＳＩＣ法で懸念される計算負荷の低減に有効である。しかし、この方法ではアレー形状が等間隔リニアアレーであることが要求され、また受信波の位相項を２次の級数近似として扱うため、アレー素子間隔が１／４波長以下であることが要求される。更に、アレー開口の増大に伴い近似誤差も増大するという課題がある。そこで、アレー形状に制約が無く、波源が近傍である受信波の位置や方向の推定において計算負荷が低減されつつ、推定精度の劣化が抑制された推定装置が求められていた。
【０００９】
この発明は、こうした実情に鑑みてなされたものであり、その目的は、近傍波源の受信波の場合に、計算負荷を低減しつつ、推定精度の劣化が抑制された波源位置推定装置及び波源方向推定装置、並びに受信波の波源位置推定方法及び波源方向推定方法を提供することにある。
【課題を解決するための手段】
【００１０】
以下、上記目的を達成するための手段及びその作用効果について説明する。
請求項１に記載の発明は、複数の受信素子により電波又は音波を受信して受信波の波源の位置を推定する波源位置推定装置において、各前記受信素子が受信した受信波の入力ベクトルを算出する入力ベクトル算出手段と、代数的手法を用いた波源の位置算出の際に近傍波源の位置特定を可能とする変換行列を算出する変換行列算出手段と、前記入力ベクトルと前記変換行列とを基に前記代数的手法によって波源の位置を算出する波源位置算出手段とを備えることをその要旨としている。
【００１１】
同構成によれば、複数の受信素子により電波又は音波を受信して波源の位置を算出するにあたって、各受信素子が受信した受信波の入力ベクトルを入力ベクトル算出手段により算出し、代数的手法で近傍波源の位置特定を可能とする変換行列を変換行列算出手段によって算出し、波源位置算出手段により入力ベクトルと変換行列とを基に代数的手法としてＥＳＰＲＩＴ法における代数的処理によって波源の位置を算出する。従来、近傍波源の位置を推定する場合にはＭＵＳＩＣ法等のサーチ法が用いられるが、サーチに要する計算負荷が大きいという問題点がある。また、代数的手法により近傍波源の位置を推定する場合には、受信信号を近似する必要があり、近似に起因する誤差が発生する。また、受信素子配置（アレー）形状も等間隔リニアアレーかつ素子間隔が１／４波長以下である必要がある。そこで、本装置では、近傍波源の受信波の場合に、変換行列を使用することで代数的手法の適用を可能とするので、任意のアレー形状において計算負荷を低減できる。また、演算に際して入力値を近似することもないので、高い位置推定精度を確保することができる。
【００１２】
請求項２に記載の発明は、複数の受信素子により電波又は音波を受信して受信波の波源の方向を推定する波源方向推定装置において、各前記受信素子が受信した受信波の入力ベクトルを算出する入力ベクトル算出手段と、代数的手法を用いた波源の方向算出の際に波源の方向特定を可能とする変換行列を算出する変換行列算出手段と、前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の方向を算出する波源方向算出手段とを備えることをその要旨としている。
【００１３】
同構成によれば、複数の受信素子により電波又は音波を受信して波源の方向を算出するにあたって、入力ベクトル算出手段により各受信素子が受信した受信波の入力ベクトルを算出し、代数的手法で近傍波源の方向特定を可能とする変換行列を変換行列算出手段によって算出し、波源方向算出手段により入力ベクトルと変換行列とを基に代数的手法としてＥＳＰＲＩＴ法における代数的処理によって波源の方向を算出する。従来、到来波の方向を代数的手法により推定する場合には、受信素子配置（アレー）形状が、例えば等間隔リニアアレー等に制約される。しかしながら、本装置では、任意形状のアレーで受信した場合においても、変換行列を用いることにより代数的手法の適用を可能とするため、計算負荷を低減できる。
【００１４】
請求項３に記載の発明は、請求項１に記載の波源位置推定装置において、前記変換行列算出手段は、前記代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値と、波源の存在し得る位置とが、一対一対応となる前記変換行列を算出することをその要旨としている。
【００１５】
同構成によれば、波源位置を固有値から算出するため、簡単な代数計算を行えばよく、計算負荷を低減できる。
請求項４に記載の発明は、請求項３に記載の波源位置推定装置において、前記固有値は、前記波源の位置を表す複素数であることをその要旨としている。
【００１６】
同構成によれば、固有値の実部が波源位置のｘ座標、虚部がｙ座標を表すため、固有値から波源の位置を算出する際に、さらに簡単に計算することができる。
請求項５に記載の発明は、請求項１又は３に記載の波源位置推定装置において、前記変換行列は、前記入力ベクトルのうち、一部のみを変換することが可能な行列であることをその要旨としている。
【００１７】
同構成によれば、変換行列は入力ベクトルのうち、一部のみを変換することが可能である。したがって、変換行列の要素数が少なく、変換行列の算出に掛かる計算負荷を低減することができる。
【００１８】
請求項６に記載の発明は、請求項１に記載の波源位置推定装置において、前記変換行列算出手段は、前記代数的手法で用いられる回転不変式の対角行列の対角成分として算出される固有値を考慮することなく前記変換行列を算出することをその要旨としている。
【００１９】
同構成によれば、変換行列の計算に用いる固有値の値を、入力ベクトルから算出する必要がないため、この部分の計算負荷が軽減される。
請求項７に記載の発明は、請求項１、及び３〜５のいずれか一項に記載の波源位置推定装置において、前記波源位置算出手段は、前記代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値によって前記波源の位置を算出することをその要旨としている。
【００２０】
同構成によれば、代数的手法で用いられる回転不変式、すなわち受信素子全体の平行移動によって生じる各波の位相回転、及び強度の変化が成立する式の対角行列の成分として算出された固有値から波源位置を算出する際に、簡単な代数演算を行えばよいため、計算負荷が軽減される。
【００２１】
請求項８に記載の発明は、請求項１、及び３〜６のいずれか一項に記載の波源位置推定装置において、前記波源位置算出手段は、前記代数的手法で得られるモードベクトルによって前記波源の位置を算出することをその要旨としている。
【００２２】
同構成によれば、代数的手法で算出されたモードベクトルから波源位置を算出するため、変換行列を計算する際に固有値について考慮する必要がなく、変換行列算出のための計算負荷が軽減される。
【００２３】
請求項９に記載の発明は、請求項１、及び３〜８のいずれか一項に記載の波源位置推定装置において、前記波源位置算出手段は、前記入力ベクトルの自己相関行列を固有値展開して得られる信号部分空間に属する固有ベクトルを前記変換行列により変換して波源位置を算出することをその要旨としている。
【００２４】
同構成によれば、波源位置算出手段が行う固有値展開の演算量を削減することができる。
請求項１０に記載の発明は、請求項１、及び３〜８のいずれか一項に記載の波源位置推定装置において、前記波源位置算出手段は、前記入力ベクトルを前記変換行列により変換して波源位置を算出することをその要旨としている。
【００２５】
同構成によれば、波源位置算出手段で行われる入力ベクトルの自己相関行列に対する固有値展開において、一般固有値展開を適用することができる。したがって、変換行列がユニタリ行列でない場合においても、雑音による推定精度の劣化を最小限に留めることができる。
【００２６】
請求項１１に記載の発明は、請求項１、及び３〜１０のいずれか一項に記載の波源位置推定装置において、前記変換行列算出手段は、窓関数を用いて前記変換行列を算出することをその要旨としている。
【００２７】
同構成によれば、窓関数を用いて変換行列を算出するため、推定領域全体を変換対象として変換行列を算出した場合と比較し、変換の際に発生する変換誤差を抑制することができる。
【００２８】
請求項１２に記載の発明は、請求項１１に記載の波源位置推定装置において、前記窓関数の位置を初期推定に基づいて算出することをその要旨としている。
同構成によれば、窓関数の位置を初期推定に基づいて算出するため、初期推定を行わずに窓関数の位置を決定した場合と比較し、発生する変換誤差をより抑制することができる。
【００２９】
請求項１３に記載の発明は、請求項１、及び３〜１２のいずれか一項に記載の波源位置推定装置において、前記波源位置推定を繰り返し行うことをその要旨としている。
同構成によれば、波源位置推定を繰り返し行うため、位置の推定精度を更に向上させることができる。
【００３０】
請求項１４に記載の発明は、請求項１３に記載の波源位置推定装置において、前記波源位置推定を繰り返す際に、推定範囲を狭くすることをその要旨としている。
同構成によれば、波源位置推定を繰り返す際に推定範囲を狭くするため、繰り返す度に位置の推定精度を更に向上させることができる。
【００３１】
請求項１５に記載の発明は、請求項１、及び３〜１４のいずれか一項に記載の波源位置推定装置において、前記受信素子は、円形に配置されることをその要旨としている。
同構成によれば、窓関数内部における変換誤差の不均一性を抑制することができるため、位置の推定精度を更に向上させることができる。
【００３２】
請求項１６に記載の発明は、複数の受信素子により電波又は音波を受信して受信波の波源の位置を推定する波源位置推定方法において、各前記受信素子が受信した受信波の入力ベクトルを算出し、異なる前記受信素子又は異なる波源における受信波に対応するよう変換することが可能な変換行列を算出し、前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の位置を算出することをその要旨としている。
【００３３】
請求項１７に記載の発明は、複数の受信素子により電波又は音波を受信して受信波の波源の方向を推定する波源方向推定方法において、各前記受信素子が受信した受信波の入力ベクトルを算出し、異なる前記受信素子又は異なる波源における受信波に対応するよう変換することが可能な変換行列を算出し、前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の方向を算出することをその要旨としている。
【発明の効果】
【００３４】
本発明によれば、近傍波源の受信波の場合に、計算負荷を低減しつつ、推定精度の劣化が抑制された波源位置推定装置及び波源方向推定装置、並びに受信波の波源位置推定方法及び波源方向推定方法を提供することができる。
【図面の簡単な説明】
【００３５】
【図１】電子キーシステムの構成を示すブロック図。
【図２】車両の通信エリアを示す上面図。
【図３】アレーアンテナと波源との位置関係を示す図。
【図４】初期推定位置と窓関数の位置関係。
【図５】シミュレーション条件を示す図。
【図６】本実施形態の信号対雑音比と位置推定誤差との関係を示す図。
【図７】本実施形態の繰り返し回数と位置推定誤差との関係を示す図。
【図８】本実施形態の繰り返し回数と位置推定誤差との関係を示す図。
【発明を実施するための形態】
【００３６】
（第１の実施形態）
以下、本発明にかかる波源位置推定装置及び波源位置推定方法を車両の電子キーシステムに適用した第１の実施形態について図１〜図８を参照して説明する。
【００３７】
図１に示されるように、車両２には、例えば運転者が実際に車両キーを操作しなくてもドアロックの施解錠やエンジンの始動及び停止等の車両動作を行うことが可能な電子キーシステム、いわゆるキー操作フリーシステムが搭載されている。キー操作フリーシステムは、キー固有のＩＤコードを無線通信で発信可能な電子キー１が車両キーとして使用され、電子キー１と車両２との狭域無線通信によりＩＤ照合が実行される。
【００３８】
キー操作フリーシステムには、車両２に近づいたり離れたりした際に自動でＩＤ照合が行われてドアロックの施解錠が許可されるエントリーシステムがある。これを以下に説明すると、車両２には、電子キー１との間で狭域無線通信（以下、スマート通信という）を行う際にＩＤ照合を行う照合ＥＣＵ（Ｅｌｅｃｔｒｏｎｉｃ Ｃｏｎｔｒｏｌ Ｕｎｉｔ）２１と、車両２の状態を管理するメインボディＥＣＵ３１とが設けられている。照合ＥＣＵ２１には、車両２の室内天井に埋設されて（図２参照）車外及び車内にＵＨＦ（Ｕｌｔｒａ Ｈｉｇｈ Ｆｒｅｑｕｅｎｃｙ）帯の信号を発信及び受信可能な通信機２４が接続されている。通信機２４には、複数のアンテナ素子からなるアレーアンテナ２４ａが搭載されている。照合ＥＣＵ２１には、例えばドアロック施解錠等を管理するメインボディＥＣＵ３１が車内ＬＡＮ（Ｌｏｃａｌ Ａｒｅａ Ｎｅｔｗｏｒｋ）３０を介して接続されている。
【００３９】
一方、電子キー１には、車両２との間で電子キーシステムに準じた無線通信を行う際のコントロールユニットとして通信制御部１１が設けられている。通信制御部１１は、固有のキーコードとしてＩＤコードが記憶されたメモリ１１ａを備えている。通信制御部１１には、外部で発信されたＵＨＦ帯の信号を受信可能な受信部１２と、通信制御部１１の指令に従いＵＨＦ帯の信号を発信可能な発信部１３とが接続されている。
【００４０】
照合ＥＣＵ２１は、通信機２４からＵＨＦ帯のリクエスト信号Ｓｒｑを間欠的に発信させ、車両２を含む周辺に通信エリアＡを形成する（図２参照）。電子キー１がこの通信エリアＡに入り込んでリクエスト信号Ｓｒｑを受信部１２で受信すると、電子キー１はリクエスト信号Ｓｒｑに応答する形で、自身のメモリ１１ａに登録されたＩＤコードを含ませたＵＨＦ帯のＩＤコード信号Ｓｉｄを発信部１３から返信する。照合ＥＣＵ２１は、通信機２４でＩＤコード信号Ｓｉｄを受信すると、ＩＤコード信号Ｓｉｄを発信した電子キー１の位置を位置推定部２１ｂにより推定する。位置推定部２１ｂが受信波の波源位置推定装置に相当する。照合ＥＣＵ２１は、自身のメモリ２１ａに登録されたＩＤコードと電子キー１のＩＤコードとを照らし合わせてＩＤ照合を行う。照合ＥＣＵ２１は、電子キー１が車外（Ａ−Ａｉ）に位置すると推定して照合が成立する（車外照合成立）と、メモリ２１ａに車外照合フラグを一定時間立てて、メインボディＥＣＵ３１を介してドアロック装置３８によるドアロックの施解錠を許可又は実行する。
【００４１】
また、電子キーシステムには、エンジン始動停止操作の際に実際の車両キー操作を必要とせずに単なるスイッチ操作のみでエンジン始動停止操作を行うことが可能な機能としてワンプッシュエンジンスタートシステムがある。このワンプッシュエンジンスタートシステムを以下に説明すると、車両２には、照合ＥＣＵ２１のＩＤ照合成立結果を基に、エンジンの点火制御及び燃料噴射制御を行うエンジンＥＣＵ３２が設けられている。エンジンＥＣＵ３２は、車内ＬＡＮ３０を通じて照合ＥＣＵ２１等の各種ＥＣＵに接続されている。車両２の運転席には、車両２の電源状態（電源ポジション）を切り換える際に操作されるエンジンスイッチ３５が設けられている。エンジンスイッチ３５は、押し操作される度に車両２の電源状態をＡＣＣ（Ａｃｃｅｓｓｏｒｙ）オン→ＩＧ（Ｉｇｎｉｔｉｏｎ）オン→電源オフの純に繰り返し遷移させ、エンジン停止時にブレーキペダルが踏み込まれた状態で操作されると、エンジンを始動に切り換える。
【００４２】
照合ＥＣＵ２１は、車外照合が成立してドアロックが解錠された後、ドアが開けられて運転者が乗車したことを例えばカーテシスイッチ３７で認識すると、通信機２４からリクエスト信号Ｓｒｑを発信して通信エリアＡを形成する。照合ＥＣＵ２１は、電子キー１がこの通信エリアＡに入り込んで返信してきたＩＤコード信号Ｓｉｄを通信機２４で受信すると、ＩＤコード信号Ｓｉｄを発信した電子キー１の位置を位置推定部２１ｂにより推定する。照合ＥＣＵ２１は、自身に登録されたＩＤコードと電子キー１のＩＤコードとを照らし合わせてＩＤ照合を行う。照合ＥＣＵ２１は、電子キー１が車内（Ａｉ）に位置すると推定して照合の成立（車内照合成立）を確認すると、エンジンスイッチ３５による電源状態の切替を許可する。
【００４３】
次に、照合ＥＣＵ２１の位置推定部２１ｂによる位置推定方法について図３〜図８を参照して説明する。本発明では、回転不変式が成立する仮想的なアレーへの変換を行う変換行列を用いることで、代数的手法であるＥＳＰＲＩＴ法の適用を可能とする。よって、定式化、ＥＳＰＲＩＴ法、近傍波源に対するＥＳＰＲＩＴ法の適用、初期推定誤差と窓関数の順番で説明する。なお、位置推定部２１ｂが、入力ベクトル算出手段、変換行列算出手段、及び波源位置算出手段として機能する。
【００４４】
＜定式化＞
波源の個数をＬとする。図３は（ｘ，y）＝（ｘｌ，ｙｌ）＝βｌに位置する第ｌ波（ｌ＝１，２，…，Ｌ）の狭帯域波源と、ｘ軸上に配置された素子間隔ｄのＫ素子等間隔リニアアレーとの位置関係を示した図である。Ｌ波の到来波がアレーアンテナに入射した場合、ｋ番目（ｋ＝１，２，…，Ｋ）のアンテナ素子の受信信号ｘｋは次式のように表すことができる。
【００４５】
【数１】
ここに、Ｆｌはｌ番目の波源の複素振幅であり、ｄｋ（βｌ）はｌ番目の波源とｋ番目のアレー素子との距離、ηｋ（βｌ）はｌ番目の波源の位置βｌからアレー素子位置までの伝搬による位相遅れである。また、ｗｋはｋ番目のアレー素子における熱雑音、λｌはｌ番目の波源の波長を表す。式（１）を用い、入力ベクトルＸを次のように定義する。
【００４６】
【数２】
このとき、振幅ベクトルＦ、アレーのモード行列Ａ、熱雑音ベクトルＷを次のようにそれぞれ定義する。
【００４７】
【数３】
なお、Ｔは転置を表す。このとき、アレーアンテナへの入力ベクトルは次のように表すことができる。
【００４８】
【数４】
このようにして定義された入力ベクトルＸの相関行列Ｒｘｘは、次式のように表すことができる。
【００４９】
【数５】
ここで、Ｈは複素共役転置を、Ｅ［・］は期待値を表す。各アンテナの熱雑音は互いに無相関であり、入射波とも無相関であるので、スナップショット数が無限である場合、熱雑音電力をσ２とすると、Ｗ＝σ２Ｉ（Ｉ：単位行列）となる。
【００５０】
＜ＥＳＰＲＩＴ法＞
ＥＳＰＲＩＴ法は、波源が十分遠方であり、到来波が平面波とみなせることを前提として成立する。今、十分遠方に存在するｌ番目（ｌ＝１，２，・・・，Ｌ）の波源が、ｙ軸となす角θl（反時計回り方向を正とする）から到来する場合を考える。すると、波源が近傍である場合に式（８）で表されたモードベクトルａ（βｌ）の成分ａｋ（βｌ）は、次式のように変更される。
【００５１】
【数６】
ここで、第１，２，…，Ｋ−１素子からなるサブアレーを第１サブアレー、第２，３，…，Ｋ素子からなるサブアレーを第２サブアレーと呼ぶこととする。このとき、Ｊ１，Ｊ２を次のように定義する。
【００５２】
【数７】
第１サブアレーのモード行列はＪ１Ａ、第２サブアレーのモード行列はＪ２Ａと表すことができる。すると、これらサブアレーのモード行列の間には次の関係式が成り立つ。
【００５３】
【数８】
ここに、diag［・］は対角行列を表す。式（１７）は回転不変式（ｒｏｔａｔｉｏｎａｌ ｉｎｖａｒｉａｎｃｅ）と呼ばれ、ＥＳＰＲＩＴ法の基本原理を表す式である。これは第１サブアレーの第ｌ波のモードベクトルＪ１ａ（βｌ）をexp｛ｊφｌ｝倍すると、第２サブアレーのモードベクトルＪ２ａ（βｌ）となることを意味し、サブアレー同士が同形かつ、波源が遠方であるからこそ成り立つものである。
【００５４】
さて、ここで式（１１）により得られた相関行列Ｒｘｘに対する固有値問題を解く。
【００５５】
【数９】
上記を解いて得られるＫ個の固有値を、値の大きいものから順にλ１，λ２，…，λＫとする。言い換えれば、式（２１）のΛの対角成分λｋが固有値である。このうち、大きいものからＬ個の固有値に対応する固有ベクトルｖ１，ｖ２，…，ｖＬの張るＬ次元部分空間は、行列Ａを構成するＬ個の列ベクトルの張る部分空間に一致する。したがって、以下を満たすＬ次の正則な行列Ｔが唯一存在する。
【００５６】
【数１０】
したがって、式（２３）をＡ＝ＥＳＴ−１と変形し、式（１７）へ代入すると次式が導かれる。ここで、式（２３）におけるＥＳが信号部分空間に属する固有ベクトルである。
【００５７】
【数１１】
ここで、ＥＸ，ＥＹ，Ψはそれぞれ以下のように定義した。
【００５８】
【数１２】
ＥＸ，ＥＹはそれぞれ式（２６），（２７）より計算可能であるため、式（２５）よりΨを算出可能である。よって、式（２８）より、Ψを固有値展開することでＴ−１，Φを求めることができ、Φが求まれば式（１８），（１９）より到来角度θｌを算出することができる。
【００５９】
＜近傍波源に対するＥＳＰＲＩＴ法の適用＞
上記では、波源が遠方である場合にはＥＳＰＲＩＴ法の適用が可能であることを確認した。しかし、波源が近傍である場合には、等間隔アレーであるにも関わらず、隣接素子間の受信信号の位相差が一定ではない不等間隔アレーで受信したような状態となる。したがって、ＥＳＰＲＩＴ法の原理式である回転不変式（１７）を満たすような、同形のサブアレーのペアを取り出すことができない。
【００６０】
そこで、本方法では、次式を満たす変換行列Ｂを使用する。
【００６１】
【数１３】
式（２９）は、「位置βｌ＝(ｘｌ，ｙｌ)を表すモードベクトルａ（βｌ）を、変換行列Ｂにより変換すると、ａ（βｌ）の（ｘｌ＋ｊｙｌ）倍となる」ことを意味する。そこで、このような変換行列Ｂが得られれば，以下の方法により波源位置を推定することが可能である。
【００６２】
まず、モード行列Ａを用いて式（２９）を次のように表す。
【００６３】
【数１４】
ここに、変換行列Ｂはモード行列ＡをＡΩという行列に変換することがわかる。この式（３０）を回転不変式と呼ぶ。また、ω（βｌ）が波源位置βｌと１対１対応する。次に、式（２３）をＡ＝ＥＳＴ−１と変形し、式（３０）へ代入すると，次式が導かれる。
【００６４】
【数１５】
ここで、ＥＸ，ＥＹ，Ψはそれぞれ以下のように定義した。
【００６５】
【数１６】
以上より、式（３３）は式（２５）と等しく、ＥＳＰＲＩＴ法が適用可能であることがわかる。また、式（３３）より得られるΨは、Ωが対角行列であることから、式（３６）のようにΩの対角成分を固有値、Ｔ−１の列ベクトルを固有ベクトルとする行列である。したがって、Ψを固有値展開することによりΩが得られ、更に式（３１），（３２）を用いることで波源位置βｌを求めることができる。すなわち、式（３０）を満たす変換行列Ｂを得ることができれば、近傍波源に対してＥＳＰＲＩＴ法を適用することができ、波源位置が推定可能となる。
【００６６】
式（２９）における変換は非線形変換であるため、全ての位置に対して式（２９）を満たす行列Ｂは存在しない。そこで、本方法では次式を用いて変換行列Ｂを算出する。
【００６７】
【数１７】
ここで、leftは左逆行列を表す。すなわち、ＡleftＡ＝Ｉである。モード行列ＡはＫ×Ｌ行列であり、素子数Ｋ＞波源数Ｌより縦長の行列である。また、Ａの列ベクトル、すなわち各波源の位置βｌに対するモードベクトルａ（βｌ）はそれぞれ一次独立であるため、左逆行列Ａleftが存在する。したがって、式（３７）により変換行列Ｂを算出することが可能である。モード行列Ａは波源位置βｌより算出されるものであり、推定前の段階では未知の行列だということである。したがって、まず初期推定を行い、おおよその波源位置からモード行列を算出し、これを用いて変換行列の計算を行うこととなる。初期推定の方法としては、ビームフォーマ法等が挙げられる。
【００６８】
＜初期推定誤差と窓関数＞
上記では、近傍波源に対するＥＳＰＲＩＴ法を適用する方法について説明した。続いて、変換行列Ｂを求めるために行う初期推定について説明する。また、初期推定結果には一般に誤差が含まれるため、この初期推定に含まれる誤差の影響と、その対策について説明する。
【００６９】
初期推定により得られた誤差を含む推定位置を図４のようにβｌ′＝（ｘｌ′、ｙｌ′）とする。このβｌ′から式（６），（７），（８）を用いて算出されるモード行列をＡ′とし、Ａ′から式（３７）を用いて算出される変換行列をＢ′とする。変換行列Ｂ′は、理想的な変換行列Ｂに対して誤差を有する。この誤差をＰとおくと誤差を含む変換行列Ｂ′は次式のように表される。
【００７０】
【数１８】
そして、Ｂ′により真の位置を表すモード行列Ａを変換する。
【００７１】
【数１９】
上記のように変換誤差ＰＡが発生し、モード行列を正しく変換することができない。変換誤差が存在する場合には式（３０）が成立しないため、次のような窓関数を用いてこれを対策する。なお、窓関数については、「アレー補完処理に窓関数を導入した空間スムージング型超解像側角法」信学論（Ｂ）、vol.J82-B,no.6,pp.1185-1192,June 1999を参考にする。
【００７２】
まず、窓関数Γ(βｌ′)を次のように定義する。
【００７３】
【数２０】
ここで、μを窓関数幅、Ｎを窓関数幅分割数、γｌp,q(βｌ′)を格子点と呼ぶこととする。格子点γｌp,q(βｌ′)は、位置βｌ′に対してｘ方向、ｙ方向にそれぞれ距離pμ／（Ｎ−１），qμ／（Ｎ−１）離れた位置を表す。
【００７４】
例えば、格子点γｌ-（N-1）/2,（N-1）/2(βｌ′)は、位置βｌ′に対して−ｘ方向にμ／２、＋ｙ方向にμ／２離れた位置を表し、格子点γｌ0,0(βｌ′)は位置βｌ′を表すことになる。つまり、窓関数Γ(βｌ′)は、図４のように位置βｌ′を中心とした幅μの正方形格子状の点の集合となる。このような窓関数Γ(βｌ′)を用いて、モード行列Ａ（Ａの上に・が付く）を次式のように定義する。
【００７５】
【数２１】
そして、モード行列Ａ（Ａの上に・が付く）を用いて変換行列Ｂ（Ｂの上に・が付く）を求める。ただし、モード行列Ａ（Ａの上に・が付く）はＫ×ＮＬ行列であるため、Ｋ＜ＮＬの場合には左逆行列Ａleft（Ａの上に・が付く）が存在せず、式（３７）により変換行列Ｂ（Ｂの上に・が付く）を計算することはできない。この場合には、次式を最小化するＢ（Ｂの上に・が付く）を変換行列とする。
【００７６】
【数２２】
ここで、||・||ＦはＦｒｏｂｅｎｉｕｓ ｎｏｒｍを表し、Ω（Ωの上に・が付く）は式（３１）においてβｌの代わりにΓ(βｌ′)を代入し算出されるＮＬ×ＮＬの対角行列である。このようにして求められた変換行列Ｂ（Ｂの上に・が付く）は、モード行列Ａ（Ａの上に・が付く）を行列ＡΩ（ＡとΩとの上に・がそれぞれ付く）に変換する際に誤差が生じてしまう。しかし、窓関数幅μが十分小さく、また窓関数幅分割数Ｎが十分大きければ、窓関数内では誤差は小さく、誤差は発生しないとみなすことができる。したがって、窓関数Γ(βｌ′)の内部に真の波源位置βｌが存在すれば、次式が成立する。
【００７７】
【数２３】
すなわち、初期推定位置βｌ′に誤差が含まれていても、窓関数Γ(βｌ′)を用いて求めた変換行列Ｂ（Ｂの上に・が付く）により、モード行列ＡをＡΩに誤差なく変換することができる。よって、上記で説明したように、ＥＳＰＲＩＴ法により波源位置を推定することが可能となる。
【００７８】
実際には、真の波源位置βlが窓関数内に含まれるように、窓関数幅μは初期推定誤差よりも大きくとる必要がある。この場合、式（３９）のように、窓関数幅μの大きさに応じて変換誤差ＰＡが発生する。しかし、変換誤差ＰＡが十分小さい場合には、変換誤差ＰＡに起因する推定誤差も小さく、ＥＳＰＲＩＴ法による推定位置βｌ′′は、初期推定位置βｌ′よりも真の波源位置βｌに近づくことが期待される。また、変換誤差ＰＡは窓関数幅μの縮小に応じて小さくなる。このため、推定位置βｌ′′を初期推定位置とし、窓関数幅μを縮小して繰り返しＥＳＰＲＩＴ法を適用することにより、さらに推定精度を改善することも可能である。
【００７９】
＜シミュレーション＞
計算機シミュレーションにより上記推定方法の評価を行う。本評価では波源数は既知とし、真の波源位置から最も近い推定位置までの距離を推定誤差として、ＲＭＳＥによる定量評価を実施する。ここでは、それぞれの真の波源位置から最も近い推定位置が重複する場合を推定失敗として定義し、位置推定に失敗した事象はＲＭＳＥの計算に含めないこととする。また、ＳＮＲを計算する際の信号強度については、アレー中心（＝原点）における第１波源の信号強度を用いている。
【００８０】
ここでは、ＥＳＰＲＩＴ法による推定位置を初期推定位置とし、窓関数幅μを縮小して再びＥＳＰＲＩＴ法を適用する「繰り返し推定」を行うこととし、窓関数幅μの初期値として２λを与える。波源数は２個とし、初期推定誤差としては各波源に対して＋ｘ方向、＋ｙ方向にそれぞれ√（１／２）λの誤差を与える。つまり、ビームフォーマ法等を用いた初期推定は行わないこととする。距離にすると１λの誤差を各波源に与えることになるが、このとき真の波源位置βｌは窓関数Γ(βｌ′)の内部に含まれることとなる。
【００８１】
＜シミュレーション結果＞
それでは、具体的にシミュレーション結果を確認する。図６は、ＳＮＲと推定誤差の関係をプロットしたものである。繰り返し推定の回数はパラメータとし、４回〜３２回まで変化させた。参考のためＭＵＳＩＣ法による推定結果も併せて示す。図６より、ＭＵＳＩＣ法よりは若干劣るが、十分な繰り返しを行えば本方法は高い推定精度が得られることがわかる。また、本シミュレーション条件では、繰り返し回数は１２回程度で十分であることもわかる。しかし、図６は、シミュレーション前に初期窓関数幅、窓関数縮小率の最適化を行った上で得られた結果である。波源位置、アレー形状、初期推定誤差等が変化すれば、最適な初期窓関数幅、窓関数縮小率も変化することが予想される。
【００８２】
そこで、初期窓関数幅により推定結果がどの様に変化するかを確認する。図７は初期窓関数幅を２λから徐々に小さくした場合の、ある一回の独立試行における繰り返し回数と推定誤差の関係をプロットしたものである。なお、熱雑音による推定誤差が発生しないよう、熱雑音はないものとしてシミュレーションを行った。図７より、初期窓関数幅μ＝２λ，１．６λの場合には、繰り返し回数に応じて推定誤差が減少しており、波源位置を正しく推定できることがわかる。しかし、初期窓関数幅μ＝１．２λ，０．８λの場合には、繰り返しの途中で推定誤差が急激に増大しており、波源位置を推定できていないことが確認される。これは、初期窓関数幅μ＝１．２λ，０．８λの場合には、初期推定誤差に対して初期窓関数幅が小さく、初期窓関数内に真の波源位置が含まれない状態となっているためと考えられる。このような状態では、変換行列Ｂ（Ｂの上に・が付く）によるモード行列Ａの変換が正しく行われない、すなわち式（４５）が成立しなくなるため、正しく推定することができない。したがって初期窓関数幅は、初期推定誤差の大きさに応じてある程度広くとる必要がある。図８は、図７とは逆に、初期窓関数幅を１λから徐々に大きくした場合の、ある一回の独立試行における繰り返し回数と推定誤差の関係をプロットしたものである。ここでも、先程と同様に熱雑音はないものとした。また、初期推定誤差についても存在しない状態、つまり初期推定位置として真の波源位置を与えた状態でシミュレーションを行った。これは、変換行列Ｂ（Ｂの上に・が付く）によりモード行列Ａを変換する際に発生する誤差だけが、推定誤差の要因となる状態でシミュレーションを行うためである。したがって、繰り返し回数が１回目のときには推定誤差が０であるため、グラフ上に値がプロットされていない。繰り返し回数が２回目のときには推定誤差が発生しているが、これは変換誤差に起因するものである。図８より、初期窓関数幅μ＝１λ，２λの場合には、波源位置を正しく推定できることがわかる。
【００８３】
この結果は、初期窓関数幅μの大きさに上限があることを意味している。先程も述べたとおり、真の波源位置が窓関数内に含まれるためには、初期推定に含まれる誤差は、初期窓関数幅の上限より小さくなければならない。
【００８４】
さて、近傍波源に対する位置推定方法として、回転不変式の成立する仮想的アレーへの変換を行う変換行列を用いて、ＥＳＰＲＩＴ法を適用可能とする方法を説明した。シミュレーションでは７素子等間隔リニアアレーを用い、近傍に位置する２つの波源に対して本方法を適用し、ＭＵＳＩＣ法より若干劣る程度の高い推定精度が得られることを明らかにした。近傍波源をサーチすることなく推定でき、計算負荷を低減することができる。しかもアレー形状に対する制約が無いということは、非常に有用である。
【００８５】
以上、説明した実施形態によれば、以下の作用効果を奏することができる。
（１）アレーアンテナにより電波を受信して波源の位置を算出するにあたって、各アレーが受信した受信波の入力ベクトルＸを算出し、代数的手法で近傍波源の位置特定を可能とする変換行列Ｂを算出し、入力ベクトルＸと変換行列Ｂとを基に代数的手法としてＥＳＰＲＩＴ法における代数的処理によって波源の位置を算出する。なお、入力ベクトルＸは、自身の相関行列Ｒｘｘの固有ベクトルＶの中で信号部分空間に属する固有ベクトルＥＳとして用いる。従来、近傍波源の位置を推定する場合にはＭＵＳＩＣ法等のサーチ法が用いられるが、サーチに要する計算負荷が大きいという問題点がある。また、代数的手法により近傍波源の位置を推定する場合には、受信信号を近似する必要があり、近似に起因する誤差が発生する。また、アレー形状も等間隔リニアアレーかつ素子間隔が１／４波長以下である必要がある。そこで、本装置では、近傍波源の受信波の場合に、変換行列Ｂを使用することで代数的手法の適用を可能とするので、計算負荷を低減できる。また、演算に際して入力値を近似することもないので、高い位置推定精度を確保することができる。
【００８６】
（２）固有値は、式（３１）における対角行列Ωの対角成分ω（βｌ）であって、ω（βｌ）は式（３２）に示されるように、波源位置と１対１対応する。波源位置を固有値から算出するため、簡単な代数計算を行えばよく、計算負荷を低減できる。
【００８７】
（３）固有値ω（βｌ）の実部が波源位置βｌのｘ座標、虚部がｙ座標を表すため、固有値から波源の位置を算出する際に、さらに簡単に計算することができる。
（４）回転不変式（３０）の中の対角行列Ωの成分ω（βｌ）が固有値に相当する。この対角行列Ωは、式（３６）のΨを固有値展開することにより得られる。すなわち、代数的手法で用いられる回転不変式、すなわちアレー全体の平行移動によって生じる各波の位相回転、及び強度の変化が成立する式の対角行列の成分として算出された固有値から波源位置を算出する際に、簡単な代数演算を行えばよいため、計算負荷が軽減される。
【００８８】
（５）入力ベクトルＸの自己相関行列Ｒｘｘを固有値展開して得られる信号部分空間に属する固有ベクトルＥＳを変換行列Ｂにより変換して波源位置を算出する。このため、固有値展開の演算量を削減することができる。
【００８９】
（６）窓関数を用いて変換行列を算出するため、推定領域全体を変換対象として変換行列を算出した場合と比較し、変換の際に発生する変換誤差を抑制することができる。
（７）窓関数の位置を初期推定に基づいて算出するため、初期推定を行わずに窓関数の位置を決定した場合と比較し、発生する変換誤差をより抑制することができる。
【００９０】
（８）波源位置推定を繰り返し行うため、位置の推定精度を更に向上させることができる。
（９）波源位置推定を繰り返す際に推定範囲を狭くするため、繰り返す度に位置の推定精度を更に向上させることができる。
【００９１】
（第２の実施形態）
以下、本発明にかかる波源位置推定装置及び波源位置推定方法を車両の電子キーシステムに適用した第２の実施形態について説明する。この実施形態では、代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値を考慮することなく変換行列Ｂを算出する点、変換行列Ｂを入力ベクトルＸに掛けて変換する点、及びモードベクトルによって波源の位置を算出する点が上記第１の実施形態と異なっている。以下、第１の実施形態との相違点を中心に説明する。
【００９２】
まず、＜ＥＳＰＲＩＴ法＞まで第１の実施形態と同様である。
次に、近傍波源に対するＥＳＰＲＩＴ法の適用について説明する。
本方法では、次式を満たす変換行列Ｂを使用する。
【００９３】
【数２４】
式（４６）は、「位置βｌ＝(ｘｌ，ｙｌ)を表すモードベクトルａ（βｌ）を、変換行列Ｂにより変換すると、ａ（βｌ）のω（βｌ）倍となる」ことを意味する。そこで、このような変換行列Ｂが得られれば，以下の方法により波源位置を推定することが可能である。
【００９４】
まず、モード行列Ａを用いて式（４６）を次のように表す。
【００９５】
【数２５】
ここに、変換行列Ｂはモード行列ＡをＡΩという行列に変換することがわかる。この式（４８）を回転不変式と呼ぶ。
【００９６】
次に、この変換行列Ｂを用いて、次式で表される入力ベクトルＸＢを定義する。
【００９７】
【数２６】
このように定義された入力ベクトルＸＢの相関行列ＲｘｘＢは次式のように表される。
【００９８】
【数２７】
この相関行列ＲｘｘＢを、Ｂ２を用いて次のように一般固有値展開する。
【００９９】
【数２８】
上記を解いて得られるＫ個の固有値を、値の大きいものから順にλ１，λ２，…，λＫとする。言い換えれば、式（５４）のΛの対角成分λＫが固有値である。このうち、大きいものからＬ個の固有値に対応する固有ベクトルｖ１，ｖ２，…，ｖＬからなる行列をＥＳと定義する。このとき、ＥＳとモード行列Ａの間には、次の関係式が成り立つ。
【０１００】
【数２９】
次に、式（５６）をＡ＝ＢＨＥＳＴ−１と変形し、式（４８）に代入すると、次式が導かれる。
【０１０１】
【数３０】
ここで、ＥＸ，ＥＹ，Ψはそれぞれ以下のように定義した。
【０１０２】
【数３１】
以上より、式（５７）は式（２５）と等しく、ＥＳＰＲＩＴ法が適用可能であることがわかる。また、式（５７）より得られるΨは、Ωが対角行列であることから、式（６０）のようにΩの対角成分を固有値、Ｔ−１の列ベクトルを固有ベクトルとする行列である。したがって、Ψを固有値展開することによりＴが得られる。このＴと、信号部分空間に属する固有ベクトルＥＳを式（２３）に代入することにより、モード行列Ａを求めることができる。よって、モード行列Ａの各列ベクトルであるモードベクトルａ（βｌ）から、式（６）、（７）、（８）を用いて波源位置βｌを算出することが可能である。すなわち、式（４８）を満たす変換行列Ｂを得ることができれば、近傍波源に対してＥＳＰＲＩＴ法を適用することができ、波源位置が推定可能となる。以下、第１の実施形態と同様であるため割愛する。
【０１０３】
さて、本実施形態では、代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値を考慮することなく変換行列Ｂを算出し、変換行列Ｂを入力ベクトルＸに掛けて変換し、モードベクトルから波源の位置を算出した。このように第１の実施形態と異なる方法にしても近傍の波源位置を算出することができ、変換行列を使用することで代数的手法の適用を可能とするので、任意のアレー形状において計算負荷を低減できる。また、演算に際して入力値を近似することもないので、高い位置推定精度を確保することができる。
【０１０４】
以上、説明した実施形態によれば、第１の実施形態の（１）及び（６）〜（９）の作用効果に加え、以下の作用効果を奏することができる。
（１０）代数的手法で用いられる回転不変式の対角行列の対角成分として算出される固有値を考慮することなく変換行列を算出するため、変換行列Ｂの計算に用いる固有値の値を、入力ベクトルから算出する必要がなく、この部分の計算負荷が軽減される。
【０１０５】
（１１）代数的手法で算出されたモードベクトルから波源位置を算出するため、変換行列Ｂを計算する際に固有値について考慮する必要がなく、変換行列算出のための計算負荷が軽減される。
【０１０６】
（１２）入力ベクトルＸを変換行列Ｂにより変換して波源位置を算出する。このため、入力ベクトルＸの自己相関行列Ｒｘｘに対する固有値展開において、一般固有値展開を適用することができる。したがって、変換行列Ｂがユニタリ行列でない場合においても、雑音による推定精度の劣化を最小限に留めることができる。
【０１０７】
（第３の実施形態）
以下、本発明にかかる波源位置推定装置及び波源位置推定方法を車両の電子キーシステムに適用した第３の実施形態について説明する。この実施形態では、変換行列を入力ベクトルのうち一部のみを変換することが可能な行列とした点が上記第１の実施形態と異なっている。以下、第１の実施形態との相違点を中心に説明する。
【０１０８】
まず、＜ＥＳＰＲＩＴ法＞まで第１の実施形態と同様である。
次に、近傍波源に対するＥＳＰＲＩＴ法の適用について説明する。
本方法では、次式を満たす変換行列Ｂを使用する。
【０１０９】
【数３２】
式（６１）は、「位置βｌ＝（ｘｌ，ｙｌ）を表す第２サブアレーのモードベクトルＪ２ａ（βｌ）を変換行列Ｂにより変換すると、第１サブアレーのモードベクトルＪ１ａ（βｌ）と同一形状となり、その際に第２サブアレーの基準素子位置は変化させない」ことを意味する。そこで、このような変換行列Ｂが得られれば，以下の方法により波源位置を推定することが可能である。
【０１１０】
まず、モード行列Ａを用いて式（６１）を次のように表す。
【０１１１】
【数３３】
ここに、変換行列Ｂは第２サブアレーのモード行列Ｊ２Ａを第１サブアレーのモード行列Ｊ１ＡにΩを乗じた行列に変換することがわかる。この式（６２）を回転不変式と呼ぶ。本実施形態においては第１の実施形態と異なり、変換行列Ｂはモード行列Ａのうち、第２、３、・・・、Ｋ素子に対応するモード行列Ｊ２Ａの変換を行う。すなわち、モード行列Ａの一部のみを変換する。このような変換行列Ｂを用いると、ＢＪ２Ｘ＝ＢＪ２ＡＦ＝Ｊ１ＡΩＦのように入力ベクトルの一部のみを変換することが可能である。次に、式（２３）をＡ＝ＥＳＴ−１と変形し、式（６２）へ代入すると，次式が導かれる。
【０１１２】
【数３４】
ここで、ＥＸ，ＥＹ，Ψはそれぞれ以下のように定義した。
【０１１３】
【数３５】
以上より、式（６５）は式（２５）と等しく、ＥＳＰＲＩＴ法が適用可能であることがわかる。また、式（３３）より得られるΨは、Ωが対角行列であることから、式（６８）のようにΩの対角成分を固有値、Ｔ−１の列ベクトルを固有ベクトルとする行列である。したがって、Ψを固有値展開することによりΩが得られ、更に式（６３），（６４）を用いることで波源位置βｌを求めることができる。すなわち、式（６２）を満たす変換行列Ｂを得ることができれば、任意形状のアレーに対してＥＳＰＲＩＴ法を適用することができ、波源位置が推定可能となる。
【０１１４】
本方法では次式を用いて変換行列Ｂを算出する。
【０１１５】
【数３６】
ここで、leftは左逆行列を表す。すなわち、｛Ｊ２Ａ｝left｛Ｊ２Ａ｝＝Ｉである。第２サブアレーのモード行列｛Ｊ２Ａ｝は（Ｋ−１）×Ｌ行列であり、（Ｋ−１）≧Ｌより縦長の行列、あるいは正方行列である。また、Ｊ２Ａの列ベクトル、すなわち各波源の位置βｌに対する第２サブアレーのモードベクトルＪ２ａ（βｌ）はそれぞれ一次独立であるため、左逆行列｛Ｊ２Ａ｝leftが存在する。したがって、式（６９）により変換行列Ｂを算出することが可能である。モード行列Ａは波源位置βｌより算出されるものであり、推定前の段階では未知の行列だということである。したがって、まず初期推定を行い、おおよその波源位置からモード行列を算出し、これを用いて変換行列の計算を行うこととなる。初期推定の方法としては、ビームフォーマ法等が挙げられる。
【０１１６】
以下、第１の実施形態と同様であるため割愛する。
さて、本実施形態では、変換行列を入力ベクトルのうち一部のみを変換することが可能な行列とした。このように、変換行列は入力ベクトルのうち、一部のみを変換することが可能である。したがって、変換行列の要素数が少なく、変換行列の算出に掛かる計算負荷を低減することができる。
【０１１７】
以上、説明した実施形態によれば、第１の実施形態の（１）、（２）及び（４）〜（９）の作用効果に加え、以下の作用効果を奏することができる。
（１３）変換行列は、入力ベクトルのうち一部のみを変換することが可能な行列とした。したがって、変換行列の要素数が少なく、変換行列の算出に掛かる計算負荷を低減することができる。
【０１１８】
（第４の実施形態）
以下、本発明にかかる波源方向推定装置及び波源方向推定方法を車両の電子キーシステムに適用した第４の実施形態について説明する。この実施形態では、位置推定ではなく、方向推定である点が上記第１の実施形態と異なっている。以下、第１の実施形態との相違点を中心に説明する。なお、位置推定部２１ｂが、入力ベクトル算出手段、変換行列算出手段、及び波源方向算出手段として機能する。
【０１１９】
まず、＜ＥＳＰＲＩＴ法＞まで第１の実施形態と同様である。
次に、任意形状アレーに対するＥＳＰＲＩＴ法の適用について説明する。
＜ＥＳＰＲＩＴ法＞では、アレーアンテナの形状が等間隔リニアアレーである場合にＥＳＰＲＩＴ法の適用が可能であることを確認した。しかし、アレーアンテナの形状が任意である場合には、隣接素子間の受信信号の位相差が一定でないため、ＥＳＰＲＩＴ法の原理式である回転不変式（１７）を満たすような、同形のサブアレーのペアを取り出すことができない。
【０１２０】
そこで、本方法では、次式を満たす変換行列Ｂを使用する。
【０１２１】
【数３７】
式（７１）は、「到来波の方向θｌを表すモードベクトルａ（βｌ）を変換行列Ｂにより変換すると、ａ（βｌ）の（ｅｘp｛ｊθｌ｝）倍となる」ことを意味する。そこで、このような変換行列Ｂが得られれば，以下の方法により波源方向を推定することが可能である。
【０１２２】
まず、モード行列Ａを用いて式（６１）を次のように表す。
【０１２３】
【数３８】
ここに、変換行列Ｂはモード行列ＡをＡΩという行列に変換することがわかる。この式（７２）を回転不変式と呼ぶ。また、ω（βｌ）が波源方向と１対１対応する。次に、式（２３）をＡ＝ＥＳＴ−１と変形し、式（７２）へ代入すると，次式が導かれる。
【０１２４】
【数３９】
ここで、ＥＸ，ＥＹ，Ψはそれぞれ以下のように定義した。
【０１２５】
【数４０】
以上より、式（７５）は式（２５）と等しく、ＥＳＰＲＩＴ法が適用可能であることがわかる。また、式（３３）より得られるΨは、Ωが対角行列であることから、式（７８）のようにΩの対角成分を固有値、Ｔ−１の列ベクトルを固有ベクトルとする行列である。したがって、Ψを固有値展開することによりΩが得られ、更に式（７３），（７４）を用いることで波源方向θｌを求めることができる。すなわち、式（７２）を満たす変換行列Ｂを得ることができれば、任意形状のアレーに対してＥＳＰＲＩＴ法を適用することができ、波源方向が推定可能となる。以下、第１の実施形態の波源位置βｌを波源方向θｌと置き換えれば同様であるため割愛する。
【０１２６】
さて、本実施形態では、等間隔リニアアレーではなく任意形状のアレーにおいて、近傍波源の方向を算出した。このように、任意形状のアレーにおいても波源の方向を算出することができ、変換行列を使用することで代数的手法の適用を可能とするので、計算負荷を低減することができる。また、演算に際して入力値を近似することもないので、高い方向推定精度を確保することができる。
【０１２７】
以上、説明した実施形態によれば、以下の作用効果を奏することができる。
（１）アレーアンテナにより電波を受信して波源の方向を算出するにあたって、各アレーが受信した受信波の入力ベクトルＸを算出し、代数的手法で近傍波源の方向特定を可能とする変換行列Ｂを算出し、入力ベクトルＸと変換行列Ｂとを基に代数的手法としてＥＳＰＲＩＴ法における代数的処理によって波源の方向を算出する。なお、入力ベクトルＸは、自身の相関行列Ｒｘｘの固有ベクトルＶの中で信号部分空間に属する固有ベクトルＥＳとして用いる。従来、到来波の方向を代数的手法により推定する場合には、アレー形状が、例えば等間隔リニアアレー等に制約される。しかしながら、本装置では、任意形状のアレーで受信した場合においても、変換行列を用いることにより代数的手法の適用を可能とするため、計算負荷を低減できる。
【０１２８】
（２）固有値は、式（７２）における対角行列Ωの対角成分ω（βｌ）であって、ω（βｌ）は式（７４）に示されるように、波源方向と１対１対応する。波源方向を固有値から算出するため、簡単な代数計算を行えばよく、計算負荷を低減できる。
【０１２９】
（３）固有値ω（βｌ）の位相成分が波源方向θｌを表すため、固有値から波源の方向を算出する際に、さらに簡単に計算することができる。
（４）回転不変式（７２）の中の対角行列Ωの成分ω（βｌ）が固有値に相当する。この対角行列Ωは、式（７８）のΨを固有値展開することにより得られる。すなわち、代数的手法で用いられる回転不変式、すなわちアレー全体の平行移動によって生じる各波の位相回転、及び強度の変化が成立する式の対角行列の成分として算出された固有値から波源方向を算出する際に、簡単な代数演算を行えばよいため、計算負荷が軽減される。
【０１３０】
（５）入力ベクトルＸの自己相関行列Ｒｘｘを固有値展開して得られる信号部分空間に属する固有ベクトルＥｓを変換行列Ｂにより変換して波源方向を算出する。このため、固有値展開の演算量を削減することができる。
【０１３１】
（６）窓関数を用いて変換行列を算出するため、推定領域全体を変換対象として変換行列を算出した場合と比較し、変換の際に発生する変換誤差を抑制することができる。
（７）窓関数の位置を初期推定に基づいて算出するため、初期推定を行わずに窓関数の位置を決定した場合と比較し、発生する変換誤差をより抑制することができる。
【０１３２】
（８）波源方向推定を繰り返し行うため、方向の推定精度を更に向上させることができる。
（９）波源方向推定を繰り返す際に推定範囲を狭くするため、繰り返す度に方向の推定精度を更に向上させることができる。
【０１３３】
なお、上記実施形態は、これを適宜変更した以下の形態にて実施することができる。
・上記第１、第３、及び第４の実施形態では、変換行列を、入力ベクトルすべてを変換可能な行列としたが、入力ベクトルのうち一部のみを変換することが可能な行列としてもよい。入力ベクトルのうち一部のみを変換することが可能な変換行列は、行列の要素数が少ない。したがって、変換行列の算出に掛かる計算負荷を低減することができる。
【０１３４】
・上記第２の実施形態では、変換行列を入力ベクトルのうち一部のみを変換可能な行列としたが、入力ベクトルのすべてを変換することが可能な行列としてもよい。
・上記第１、第３、及び第４の実施形態では、回転不変式の対角行列の成分として算出される固有値と、波源の存在し得る位置とが、一対一対応となる変換行列を算出したが、回転不変式の対角行列の対角成分として算出される固有値を考慮することなく変換行列を算出してもよい。
【０１３５】
・上記第１、第３、及び第４の実施形態では、位置の算出において固有値を使用したが、モードベクトルを使用してもよい。すなわち、式（３３）により得られるΨは、式（３６）となる。そして、このΨを固有値展開し、固有値Ω、固有ベクトルＴ−１を得る。次に、入力ベクトルＸの自己相関行列Ｒｘｘの固有ベクトルＶの中で信号部分空間に属する固有ベクトルＥＳと、Ψの固有ベクトルＴ−１とから式（２３）を用いてモード行列を算出する。そして、モード行列Ａの各列成分であるモードベクトルａ（βｌ）から、波源位置βｌを例えば非線形計画法により算出する。
【０１３６】
・上記第２の実施形態では、回転不変式の対角行列の対角成分として算出される固有値を考慮することなく変換行列を算出したが、回転不変式の対角行列の成分として算出される固有値と、波源の存在し得る位置とが、一対一対応となる変換行列を算出してもよい。
【０１３７】
・上記第１の実施形態では、変換行列は式（２９）の変換を行う行列としたが、式（２９）以外の変換を行う行列を変換行列としてもよい。
・上記第２の実施形態では、変換行列は式（４６）の変換を行う行列としたが、式（４６）以外の変換を行う行列を変換行列としてもよい。
【０１３８】
・上記第３の実施形態では、変換行列は式（６１）の変換を行う行列としたが、式（６１）以外の変換を行う行列を変換行列としてもよい。
・上記第１の実施形態では、変換行列は式（７１）の変換を行う行列としたが、式（７１）以外の変換を行う行列を変換行列としてもよい。
【０１３９】
・上記第１、第３、及び第４の実施形態では、変換行列は式（３７）により変換行列を算出したが、式（３７）以外の式を用いて変換行列を算出してもよい。
・上記第３の実施形態では、変換行列は式（６９）により変換行列を算出したが、式（６９）以外の式を用いて変換行列を算出してもよい。
【０１４０】
・上記第１、第３、及び第４の実施形態では、変換行列Ｂと固有ベクトルＥｓとを掛けたが、変換行列Ｂと入力ベクトルＸとを掛けて算出するようにしてもよい。
・上記第２の実施形態では、変換行列Ｂと入力ベクトルＸとを掛けたが、変換行列Ｂと固有ベクトルＥｓとを掛けて算出するようにしてもよい。
【０１４１】
・上記構成において、アレー配置を円形としてもよい。このようにすれば、窓関数内部における変換誤差の不均一性を抑制することができるため、位置又は方向の推定精度を更に向上させることができる。
【０１４２】
・上記第１〜３の実施形態では、アレー配置は等間隔リニアアレーとしたが、これに限らず任意のアレー配置としてよい。
・上記実施形態では、変換行列算出時に窓関数を用いたが、所望の推定精度を得ることができれば、窓関数を用いなくともよい。
【０１４３】
・上記実施形態では、所望の推定精度を得ることができれば、繰り返し時に推定範囲を狭くしなくてもよい。
・上記実施形態では、代数的手法としてＥＳＰＲＩＴ法を用いたが、他の代数的手法、たとえばＤＯＡ−Ｍａｔｒｉｘ法等を使用してもよい。
【０１４４】
・上記実施形態では、波源位置推定の算出処理を１２回程度繰り返し行うのがよいとしたが、所望の推定精度に合わせて適宜変更可能である。また、１回で所望の推定精度を得ることができれば、繰り返さなくてもよい。
【０１４５】
・上記実施形態では、波源位置推定装置を電子キーシステムの電子キー１の位置推定に適用したが、電子キーシステムに限らず、種々の通信機器に適用してもよい。
・また、波源位置推定手段を備えた波源位置推定装置に限らず、波源方向推定手段を備えた波源方向推定装置として種々の通信機器に適用してもよい。
【０１４６】
・波源位置推定装置及び波源方向推定装置等の装置に限らず、波源位置推定方法及び波源方向推定方法を、アレーアンテナを備えた種々の通信機器の演算に適用してもよい。
・上記構成において、図１に示されるように、電子キー１に施錠スイッチ１９ａ及び解錠スイッチ１９ｂを設け、この施錠スイッチ１９ａ及び解錠スイッチ１９ｂを操作することで、車両２へ施錠信号Ｓｌ又は解錠信号Ｓｕｌを発信して、車両２のドアロックの解錠又は施錠を可能としてもよい。
【０１４７】
・電子キーシステムで使用する各通信の周波数はどの周波数の電波を使用してもよい。
・上記実施形態において、波源は電波を発信する装置であったが、これに限らず、例えばスピーカーのような音波を発信する装置であっても良い。
【０１４８】
・上記実施形態において、受信素子は電波を受信する装置であったが、これに限らず、例えばマイクロホンのような音波を受信する装置であっても良い。
・上記実施形態では、電波に適用したが、音波に適用してもよい。
【０１４９】
次に、上記実施形態から把握できる技術的思想をその効果と共に記載する。
（イ）請求項２に記載の波源方向推定装置において、前記変換行列算出手段は、前記代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値と、波源の存在し得る方向とが、一対一対応となる前記変換行列を算出することを特徴とする波源方向推定装置。
【０１５０】
同構成によれば、波源位置を固有値から算出するため、簡単な代数計算を行えばよく、計算負荷を低減できる。
（ロ）（イ）に記載の波源方向推定装置において、前記固有値は、前記波源の方向を表す複素数であることを特徴とする波源方向推定装置。
【０１５１】
同構成によれば、固有値の位相成分が波源方向を表すため、固有値から波源の方向を算出する際に、さらに簡単に計算することができる。
（ハ）請求項２、又は（イ）に記載の波源方向推定装置において、前記変換行列は、前記入力ベクトルのうち、一部のみを変換することが可能な行列であることを特徴とする波源方向推定装置。
【０１５２】
同構成によれば、変換行列は入力ベクトルのうち、一部のみを変換することが可能である。したがって、変換行列の要素数が少なく、変換行列の算出に掛かる計算負荷を低減することができる。
【０１５３】
（ニ）請求項２に記載の波源方向推定装置において、前記変換行列算出手段は、前記代数的手法で用いられる回転不変式の対角行列の対角成分として算出される固有値を考慮することなく前記変換行列を算出することを特徴とする波源方向推定装置。
【０１５４】
同構成によれば、変換行列の計算に用いる固有値の値を、入力ベクトルから算出する必要がないため、この部分の計算負荷が軽減される。
（ホ）請求項２、及び（イ）〜（ハ）のいずれかに記載の波源方向推定装置において、前記波源方向算出手段は、前記代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値によって前記波源の方向を算出することを特徴とする波源方向推定装置。
【０１５５】
同構成によれば、代数的手法で用いられる回転不変式、すなわち受信素子全体の平行移動によって生じる各波の位相回転、及び強度の変化が成立する式の対角行列の成分として算出された固有値から波源方向を算出する際に、簡単な代数演算を行えばよいため、計算負荷が軽減される。
【０１５６】
（へ）請求項２、及び（イ）〜（ニ）のいずれかに記載の波源方向推定装置において、前記波源方向算出手段は、前記代数的手法で得られるモードベクトルによって前記波源の方向を算出することを特徴とする波源方向推定装置。
【０１５７】
同構成によれば、代数的手法で算出されたモードベクトルから波源方向を算出するため、変換行列を計算する際に固有値について考慮する必要がなく、変換行列算出のための計算負荷が軽減される。
【０１５８】
（ト）請求項２、及び（イ）〜（へ）のいずれかに記載の波源方向推定装置において、前記波源方向算出手段は、前記入力ベクトルの自己相関行列を固有値展開して得られる信号部分空間に属する固有ベクトルを前記変換行列により変換して波源方向を算出することを特徴とする波源方向推定装置。
【０１５９】
同構成によれば、波源方向算出手段が行う固有値展開の演算量を削減することができる。
（チ）請求項２、及び（イ）〜（へ）のいずれかに記載の波源方向推定装置において、前記波源方向算出手段は、前記入力ベクトルを前記変換行列により変換して波源方向を算出することを特徴とする波源方向推定装置。
【０１６０】
同構成によれば、波源方向算出手段で行われる入力ベクトルの自己相関行列に対する固有値展開において、一般固有値展開を適用することができる。したがって、変換行列がユニタリ行列でない場合においても、雑音による推定精度の劣化を最小限に留めることができる。
【０１６１】
（リ）請求項２、及び（イ）〜（チ）のいずれかに記載の波源方向推定装置において、前記変換行列算出手段は、窓関数を用いて前記変換行列を算出することを特徴とする波源方向推定装置。
【０１６２】
同構成によれば、窓関数を用いて変換行列を算出するため、推定領域全体を変換対象として変換行列を算出した場合と比較し、変換の際に発生する変換誤差を抑制することができる。
【０１６３】
（ヌ）（リ）に記載の波源方向推定装置において、前記窓関数の位置を初期推定に基づいて算出することを特徴とする波源方向推定装置。
同構成によれば、窓関数の位置を初期推定に基づいて算出するため、初期推定を行わずに窓関数の位置を決定した場合と比較し、発生する変換誤差をより抑制することができる。
【０１６４】
（ル）請求項２、及び（イ）〜（ヌ）のいずれかに記載の波源方向推定装置において、前記波源方向推定を繰り返し行うことを特徴とする波源方向推定装置。
同構成によれば、波源方向推定を繰り返し行うため、方向の推定精度を更に向上させることができる。
【０１６５】
（ヲ）（ル）に記載の波源方向推定装置において、前記波源方向推定を繰り返す際に、推定範囲を狭くすることを特徴とする波源方向推定装置。
同構成によれば、波源方向推定を繰り返す際に推定範囲を狭くするため、繰り返す度に方向の推定精度を更に向上させることができる。
【０１６６】
（ワ）請求項２、及び（イ）〜（ヲ）のいずれかに記載の波源方向推定装置において、前記受信素子は、円形に配置されることを特徴とする波源方向推定装置。
同構成によれば、窓関数内部における変換誤差の不均一性を抑制することができるため、方向の推定精度を更に向上させることができる。
【符号の説明】
【０１６７】
１…電子キー、２…車両、１１…通信制御部、１１ａ…メモリ、１２…受信部、１３…発信部、２１…照合ＥＣＵ、２１ａ…メモリ、２１ｂ…位置推定部、２４…通信機、２４ａ…アレーアンテナ、３５…エンジンスイッチ、３８…ドアロック装置、Ａ…通信エリア、Ａｉ…車内エリア、Ｓｉｄ…ＩＤコード信号、Ｓｒｑ…リクエスト信号。
【特許請求の範囲】
【請求項１】
複数の受信素子により電波又は音波を受信して受信波の波源の位置を推定する波源位置推定装置において、
各前記受信素子が受信した受信波の入力ベクトルを算出する入力ベクトル算出手段と、
代数的手法を用いた波源の位置算出の際に近傍波源の位置特定を可能とする変換行列を算出する変換行列算出手段と、
前記入力ベクトルと前記変換行列とを基に前記代数的手法によって波源の位置を算出する波源位置算出手段とを備える
ことを特徴とする波源位置推定装置。
【請求項２】
複数の受信素子により電波又は音波を受信して受信波の波源の方向を推定する波源方向推定装置において、
各前記受信素子が受信した受信波の入力ベクトルを算出する入力ベクトル算出手段と、
代数的手法を用いた波源の方向算出の際に波源の方向特定を可能とする変換行列を算出する変換行列算出手段と、
前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の方向を算出する波源方向算出手段とを備える
ことを特徴とする波源方向推定装置。
【請求項３】
請求項１に記載の波源位置推定装置において、
前記変換行列算出手段は、前記代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値と、波源の存在し得る位置とが、一対一対応となる前記変換行列を算出する
ことを特徴とする波源位置推定装置。
【請求項４】
請求項３に記載の波源位置推定装置において、
前記固有値は、前記波源の位置を表す複素数である
ことを特徴とする波源位置推定装置。
【請求項５】
請求項１又は３に記載の波源位置推定装置において、
前記変換行列は、前記入力ベクトルのうち、一部のみを変換することが可能な行列である
ことを特徴とする波源位置推定装置。
【請求項６】
請求項１に記載の波源位置推定装置において、
前記変換行列算出手段は、前記代数的手法で用いられる回転不変式の対角行列の対角成分として算出される固有値を考慮することなく前記変換行列を算出する
ことを特徴とする波源位置推定装置。
【請求項７】
請求項１、及び３〜５のいずれか一項に記載の波源位置推定装置において、
前記波源位置算出手段は、前記代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値によって前記波源の位置を算出する
ことを特徴とする波源位置推定装置。
【請求項８】
請求項１、及び３〜６のいずれか一項に記載の波源位置推定装置において、
前記波源位置算出手段は、前記代数的手法で得られるモードベクトルによって前記波源の位置を算出する
ことを特徴とする波源位置推定装置。
【請求項９】
請求項１、及び３〜８のいずれか一項に記載の波源位置推定装置において、
前記波源位置算出手段は、前記入力ベクトルの自己相関行列を固有値展開して得られる信号部分空間に属する固有ベクトルを前記変換行列により変換して波源位置を算出する
ことを特徴とする波源位置推定装置。
【請求項１０】
請求項１、及び３〜８のいずれか一項に記載の波源位置推定装置において、
前記波源位置算出手段は、前記入力ベクトルを前記変換行列により変換して波源位置を算出する
ことを特徴とする波源位置推定装置。
【請求項１１】
請求項１、及び３〜１０のいずれか一項に記載の波源位置推定装置において、
前記変換行列算出手段は、窓関数を用いて前記変換行列を算出する
ことを特徴とする波源位置推定装置。
【請求項１２】
請求項１１に記載の波源位置推定装置において、
前記窓関数の位置を初期推定に基づいて算出する
ことを特徴とする波源位置推定装置。
【請求項１３】
請求項１、及び３〜１２のいずれか一項に記載の波源位置推定装置において、
前記波源位置推定を繰り返し行う
ことを特徴とする波源位置推定装置。
【請求項１４】
請求項１３に記載の波源位置推定装置において、
前記波源位置推定を繰り返す際に、推定範囲を狭くする
ことを特徴とする波源位置推定装置。
【請求項１５】
請求項１、及び３〜１４のいずれか一項に記載の波源位置推定装置において、
前記受信素子は、円形に配置される
ことを特徴とする波源位置推定装置。
【請求項１６】
複数の受信素子により電波又は音波を受信して受信波の波源の位置を推定する波源位置推定方法において、
各前記受信素子が受信した受信波の入力ベクトルを算出し、
異なる前記受信素子又は異なる波源における受信波に対応するよう変換することが可能な変換行列を算出し、
前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の位置を算出する
ことを特徴とする波源位置推定方法。
【請求項１７】
複数の受信素子により電波又は音波を受信して受信波の波源の方向を推定する波源方向推定方法において、
各前記受信素子が受信した受信波の入力ベクトルを算出し、
異なる前記受信素子又は異なる波源における受信波に対応するよう変換することが可能な変換行列を算出し、
前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の方向を算出する
ことを特徴とする波源方向推定方法。
【請求項１】
複数の受信素子により電波又は音波を受信して受信波の波源の位置を推定する波源位置推定装置において、
各前記受信素子が受信した受信波の入力ベクトルを算出する入力ベクトル算出手段と、
代数的手法を用いた波源の位置算出の際に近傍波源の位置特定を可能とする変換行列を算出する変換行列算出手段と、
前記入力ベクトルと前記変換行列とを基に前記代数的手法によって波源の位置を算出する波源位置算出手段とを備える
ことを特徴とする波源位置推定装置。
【請求項２】
複数の受信素子により電波又は音波を受信して受信波の波源の方向を推定する波源方向推定装置において、
各前記受信素子が受信した受信波の入力ベクトルを算出する入力ベクトル算出手段と、
代数的手法を用いた波源の方向算出の際に波源の方向特定を可能とする変換行列を算出する変換行列算出手段と、
前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の方向を算出する波源方向算出手段とを備える
ことを特徴とする波源方向推定装置。
【請求項３】
請求項１に記載の波源位置推定装置において、
前記変換行列算出手段は、前記代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値と、波源の存在し得る位置とが、一対一対応となる前記変換行列を算出する
ことを特徴とする波源位置推定装置。
【請求項４】
請求項３に記載の波源位置推定装置において、
前記固有値は、前記波源の位置を表す複素数である
ことを特徴とする波源位置推定装置。
【請求項５】
請求項１又は３に記載の波源位置推定装置において、
前記変換行列は、前記入力ベクトルのうち、一部のみを変換することが可能な行列である
ことを特徴とする波源位置推定装置。
【請求項６】
請求項１に記載の波源位置推定装置において、
前記変換行列算出手段は、前記代数的手法で用いられる回転不変式の対角行列の対角成分として算出される固有値を考慮することなく前記変換行列を算出する
ことを特徴とする波源位置推定装置。
【請求項７】
請求項１、及び３〜５のいずれか一項に記載の波源位置推定装置において、
前記波源位置算出手段は、前記代数的手法で用いられる回転不変式の対角行列の成分として算出される固有値によって前記波源の位置を算出する
ことを特徴とする波源位置推定装置。
【請求項８】
請求項１、及び３〜６のいずれか一項に記載の波源位置推定装置において、
前記波源位置算出手段は、前記代数的手法で得られるモードベクトルによって前記波源の位置を算出する
ことを特徴とする波源位置推定装置。
【請求項９】
請求項１、及び３〜８のいずれか一項に記載の波源位置推定装置において、
前記波源位置算出手段は、前記入力ベクトルの自己相関行列を固有値展開して得られる信号部分空間に属する固有ベクトルを前記変換行列により変換して波源位置を算出する
ことを特徴とする波源位置推定装置。
【請求項１０】
請求項１、及び３〜８のいずれか一項に記載の波源位置推定装置において、
前記波源位置算出手段は、前記入力ベクトルを前記変換行列により変換して波源位置を算出する
ことを特徴とする波源位置推定装置。
【請求項１１】
請求項１、及び３〜１０のいずれか一項に記載の波源位置推定装置において、
前記変換行列算出手段は、窓関数を用いて前記変換行列を算出する
ことを特徴とする波源位置推定装置。
【請求項１２】
請求項１１に記載の波源位置推定装置において、
前記窓関数の位置を初期推定に基づいて算出する
ことを特徴とする波源位置推定装置。
【請求項１３】
請求項１、及び３〜１２のいずれか一項に記載の波源位置推定装置において、
前記波源位置推定を繰り返し行う
ことを特徴とする波源位置推定装置。
【請求項１４】
請求項１３に記載の波源位置推定装置において、
前記波源位置推定を繰り返す際に、推定範囲を狭くする
ことを特徴とする波源位置推定装置。
【請求項１５】
請求項１、及び３〜１４のいずれか一項に記載の波源位置推定装置において、
前記受信素子は、円形に配置される
ことを特徴とする波源位置推定装置。
【請求項１６】
複数の受信素子により電波又は音波を受信して受信波の波源の位置を推定する波源位置推定方法において、
各前記受信素子が受信した受信波の入力ベクトルを算出し、
異なる前記受信素子又は異なる波源における受信波に対応するよう変換することが可能な変換行列を算出し、
前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の位置を算出する
ことを特徴とする波源位置推定方法。
【請求項１７】
複数の受信素子により電波又は音波を受信して受信波の波源の方向を推定する波源方向推定方法において、
各前記受信素子が受信した受信波の入力ベクトルを算出し、
異なる前記受信素子又は異なる波源における受信波に対応するよう変換することが可能な変換行列を算出し、
前記入力ベクトルと前記変換行列とを基に代数的手法によって波源の方向を算出する
ことを特徴とする波源方向推定方法。
【図１】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【公開番号】特開２０１１−７５４４８（Ｐ２０１１−７５４４８Ａ）
【公開日】平成２３年４月１４日（２０１１．４．１４）
【国際特許分類】
【出願番号】特願２００９−２２８５１２（Ｐ２００９−２２８５１２）
【出願日】平成２１年９月３０日（２００９．９．３０）
【出願人】（０００００３５５１）株式会社東海理化電機製作所 (3,198)
【出願人】（３０４０２１２７７）国立大学法人 名古屋工業大学 (784)
【Ｆターム（参考）】
【公開日】平成２３年４月１４日（２０１１．４．１４）
【国際特許分類】
【出願日】平成２１年９月３０日（２００９．９．３０）
【出願人】（０００００３５５１）株式会社東海理化電機製作所 (3,198)
【出願人】（３０４０２１２７７）国立大学法人 名古屋工業大学 (784)
【Ｆターム（参考）】
[ Back to top ]