特許 6239664 周辺環境推定装置及び周辺環境推定方法

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】特許公報(B2)

(11)【特許番号】6239664

(24)【登録日】2017年11月10日

(45)【発行日】2017年11月29日

(54)【発明の名称】周辺環境推定装置及び周辺環境推定方法

(51)【国際特許分類】

   G01S 13/89 20060101AFI20171120BHJP

   G08G 1/16 20060101ALI20171120BHJP

【ＦＩ】

   G01S13/89

   G08G1/16 C

【請求項の数】14

【全頁数】26

(21)【出願番号】特願2016-52221(P2016-52221)

(22)【出願日】2016年3月16日

(65)【公開番号】特開2017-166966(P2017-166966A)

(43)【公開日】2017年9月21日

【審査請求日】2016年12月21日

(73)【特許権者】

【識別番号】502324066

【氏名又は名称】株式会社デンソーアイティーラボラトリ

(74)【代理人】

【識別番号】100093067

【弁理士】

【氏名又は名称】二瓶 正敬

(72)【発明者】

【氏名】安達 仁吾

(72)【発明者】

【氏名】山野 千晴

【審査官】 大▲瀬▼ 裕久

(56)【参考文献】

【文献】 特開２０１２−４８６４２（ＪＰ，Ａ）

【文献】 米国特許出願公開第２０１４／７１２４０（ＵＳ，Ａ１）

【文献】 KIMURA Nobutaka et al.，”Real-Time Updating of 2D Map for Autonomous Robot Locomotion Based on Distinction Between Static and Semi-Static Objects”，ＡＤＶＡＮＣＥＤ ＲＯＢＯＴＩＣＳ １１，Robotics Society of Japan，２０１２年，第26巻，p.1343-1368

【文献】 北島 健太 他，”局所地図の時系列統合による大域地図の生成 Global Map Generation by Temporal Integration of Local Maps”，ロボティクス・メカトロニクス講演会 ’０９ 講演論文集 Ｐｒｏｃｅｅｄｉｎｇｓ ｏｆ ｔｈｅ ２００９ ＪＳＭＥ Ｃｏｎｆｅｒｅｎｃｅ ｏｎ Ｒｏｂｏｔｉｃｓ ａｎｄ Ｍｅｃｈａｔｒｏｎｉｃｓ，社団法人日本機械学会，２００９年，2P1-F16 p.1-4

【文献】 MERALI S. Rehman et al.，”Optimizing Online Occupancy Grid Mapping to Capture the Residual Uncertainty”，IEEE International Conference on Robotics & Automaion(ICRA)，２０１４年，p.6070-6076

(58)【調査した分野】（Int.Cl.，ＤＢ名）

Ｇ０１Ｓ １／００−１９／５５

Ｇ０８Ｇ １／１６

ＪＳＴＰｌｕｓ／ＪＭＥＤＰｌｕｓ／ＪＳＴ７５８０（ＪＤｒｅａｍＩＩＩ）

(57)【特許請求の範囲】

【請求項1】

センシング信号を送信する送信機と、前記センシング信号が反射して戻ってくる信号を受信する受信機とを有するセンシング機器を用いて得られる観測データから、周辺の環境の物体の存在を推定する周辺環境推定装置であって、
前記観測データに基づいて前記受信器で受信した信号の強度値を計算し、前記強度値から前記受信器で受信した信号が被検出物から反射してきたものであるかどうかを表す尤度を計算する強度値尤度計算部と、
前記観測データに基づいて前記被検出物の方位角を計算し、前記強度値から前記方位角を中心とした分散の程度を表す角度の分布関数を計算する分布関数計算部と、
前記周辺の空間を格子状に分割して各格子における物体の占有確率を記述する占有格子地図に関して、前記尤度を前記分布関数によって分散させた値を用いて、前記占有格子地図内の複数の格子における前記物体の占有確率を更新する計算を行う占有格子地図計算部とを、
有し、
前記複数の観測データについて、前記強度値尤度計算部による計算、前記分布関数計算部による計算、並びに、前記占有格子地図計算部による計算を行うことで、前記複数の観測データに基づく前記占有格子地図を作成するよう構成されている周辺環境推定装置。

【請求項2】

前記占有格子地図計算部は、前記占有格子地図内の格子のうち、前記観測データに含まれる前記被検出物との距離と等しい距離を持つ複数の格子の配列に沿って、前記尤度を前記分布関数によって分散させるよう構成されている請求項１に記載の周辺環境推定装置。

【請求項3】

前記強度値に対して前記被検出物が物体又はノイズのいずれであるかを示す確率が事前に設定されており、前記強度値尤度計算部は、前記設定された確率から得られる物体の占有確率の対数オッズを計算し、その計算結果を前記尤度とするよう構成されており、
前記占有格子地図計算部は、前記尤度と前記分布関数との積を計算することで、前記占有格子地図内の複数の格子における前記物体の占有確率を更新するよう構成されている請求項１又は２に記載の周辺環境推定装置。

【請求項4】

前記占有格子地図計算部は、前記複数の観測データのそれぞれから得られる前記尤度と前記分布関数との積の総和を計算することによって、前記占有格子地図を作成するよう構成されている請求項３に記載の周辺環境推定装置。

【請求項5】

前記占有格子地図計算部は、前記複数の観測データのそれぞれから得られる前記尤度と前記分布関数との積を計算して、前記計算された積のうちの最大値をとることによって、前記占有格子地図を作成するよう構成されている請求項３に記載の周辺環境推定装置。

【請求項6】

前記観測データに基づいて前記被検出物が静止物であるか移動体であるかを判別する静止物／移動体判別部を更に有し、
前記占有格子地図計算部は、前記尤度と前記分布関数との積を計算する際に、前記被検出物が静止物であると判断された観測データについては、前記尤度と前記分布関数との積に正の値を掛け、前記被検出物が移動体であると判断された観測データについては、前記尤度と前記分布関数との積に負の値を掛けるよう構成されている請求項４又は５に記載の周辺環境推定装置。

【請求項7】

前記センシング機器によって複数の異なる局所的な領域に関して作成された複数の前記占有格子地図から、前記複数の異なる局所的な領域を含む大局的な領域に関する占有格子地図を作成するよう構成されている請求項１から６のいずれか１つに記載の周辺環境推定装置。

【請求項8】

センシング信号を送信する送信機と、前記センシング信号が反射して戻ってくる信号を受信する受信機とを有するセンシング機器を用いて得られる観測データから、周辺の環境の物体の存在を推定する周辺環境推定方法であって、
前記観測データに基づいて前記受信機で受信した信号の強度値を計算し、前記強度値から前記受信機で受信した信号が被検出物から反射してきたものであるかどうかを表す尤度を計算する強度値尤度計算ステップと、
前記観測データに基づいて前記被検出物の方位角を計算し、前記強度値から前記方位角を中心とした分散の程度を表す角度の分布関数を計算する分布関数計算ステップと、
前記周辺の空間を格子状に分割して各格子における物体の占有確率を記述する占有格子地図に関して、前記尤度を前記分布関数によって分散させた値を用いて、前記占有格子地図内の複数の格子における前記物体の占有確率を更新する計算を行う占有格子地図計算ステップとを、
有し、
前記複数の観測データについて、前記強度値尤度計算ステップによる計算、前記分布関数計算ステップによる計算、並びに、前記占有格子地図計算ステップによる計算を行うことで、前記複数の観測データに基づく前記占有格子地図を作成する周辺環境推定方法。

【請求項9】

前記占有格子地図計算ステップにおいては、前記占有格子地図内の格子のうち、前記観測データに含まれる前記被検出物との距離と等しい距離を持つ複数の格子の配列に沿って、前記尤度を前記分布関数によって分散させる請求項８に記載の周辺環境推定方法。

【請求項10】

前記強度値に対して前記被検出物が物体又はノイズのいずれであるかを示す確率が事前に設定されており、前記強度値尤度計算ステップにおいて、前記設定された確率から得られる物体の占有確率の対数オッズを計算し、その計算結果を前記尤度とし、
前記占有格子地図計算ステップにおいて、前記尤度と前記分布関数との積を計算することで、前記占有格子地図内の複数の格子における前記物体の占有確率を更新する請求項８又は９に記載の周辺環境推定方法。

【請求項11】

前記占有格子地図計算ステップにおいて、前記複数の観測データのそれぞれから得られる前記尤度と前記分布関数との積の総和を計算することによって、前記占有格子地図を作成する請求項１０に記載の周辺環境推定方法。

【請求項12】

前記占有格子地図計算ステップにおいて、前記複数の観測データのそれぞれから得られる前記尤度と前記分布関数との積を計算して、前記計算された積のうちの最大値をとることによって、前記占有格子地図を作成する請求項１０に記載の周辺環境推定方法。

【請求項13】

前記観測データに基づいて前記被検出物が静止物であるか移動体であるかを判別する静止物／移動体判別ステップを更に有し、
前記占有格子地図計算ステップは、前記尤度と前記分布関数との積を計算する際に、前記被検出物が静止物であると判断された観測データについては、前記尤度と前記分布関数との積に正の値を掛け、前記被検出物が移動体であると判断された観測データについては、前記尤度と前記分布関数との積に負の値を掛ける請求項１１又は１２に記載の周辺環境推定方法。

【請求項14】

前記センシング機器によって複数の異なる局所的な領域に関して作成された複数の前記占有格子地図から、前記複数の異なる局所的な領域を含む大局的な領域に関する占有格子地図を作成する請求項８から１３のいずれか１つに記載の周辺環境推定方法。

【発明の詳細な説明】

【技術分野】

【0001】

本発明は、センシング機器（センサー）を用いて周辺をセンシングすることで得られた観測データから、周辺の環境を推定する技術に関する。本発明は、特に、車両に搭載されたセンシング機器を用いて得られた観測データから、車両周辺の環境を推定する技術に関する。

【背景技術】

【0002】

レーダーなどのセンシング機器を用いてセンシングを行うことによって得られた観測データは、センシング機器の周辺の物体の位置や形状、物体の移動の態様（移動体又は静止物）などの推定に用いられる。観測データは、観測誤差や外乱などの不確定な要素を多く含む。そのため、観測する物体の正確な情報を得るためには、何らかのフィルタ処理が必要となる。

【0003】

センシング機器を用いて得られた観測データから、周辺の物体を認識及び推定する従来の技術の多くは、確率的な手法に基づくものである。こうした確率的な手法の１つとして、下記の非特許文献１によって最初に提案された占有格子地図（Occupancy Grid Map）を用いた手法が存在する。以下、非特許文献１によって最初に提案された手法に基づいて、下記の非特許文献２で再提案された占有格子地図の手法について簡単に説明する。

【0004】

占有格子地図を用いた手法は、観測する周辺の空間を格子（グリッド）に分けた占有格子地図を作成し、それぞれの格子１つ１つに対して、占有状態（占有状態は「ｍ」で表される）か非占有状態か、すなわち、その格子内の空間が物理的に埋まっているかどうかを表す占有確率ｐ（ｍ）を計算するというものである。図１には、従来の技術における占有格子地図の一例が示されている。

【0005】

観測という行為がt時間フレーム（時刻ｔ）において行われた場合、センサーによって検知された空間に対応する格子（図１の格子状に区切られた各領域）に対し、占有か非占有かの確率を更新計算する。非特許文献２によれば、ある格子iに関して、占有確率ｐ（ｍ）の更新式は下記の式（１）で表され、この数式によって、ｔ−１時間フレーム（時刻ｔ−１）における状態から、ｔ時間フレーム（時刻ｔ）における状態を導出することができる。

【0006】

【数1】

【0007】

式（１）では、占有確率ｐ（ｍ）に関して、占有確率ｐ（ｍ）を非占有確率１−ｐ（ｍ）で割った値の対数（ＬｏｇＯｄｄｓ：対数オッズ）を単位として加算が行われる。

【0008】

式（１）の左辺は、格子ｉにおいて、ｔ時間フレーム（時刻ｔ）までに得られた占有確率ｐ（ｍ_ｉ｜ｚ_１，．．．，ｚ_ｔ）の対数オッズであり、占有確率ｐ（ｍ）に関する最新の推定結果を表している。式（１）の左辺（すなわち、占有確率ｐ（ｍ）に関する最新の推定結果）は、式（１）の右辺を計算することによって得られる。なお、式（１）の左辺は、時刻ｔまでの観測における対数尤度比と呼ばれる。

【0009】

式（１）の右辺第１項は、格子ｉにおけるｔ時間フレーム（時刻ｔ）の観測強度ｚ_ｔの占有確率ｐ（ｍ_ｉ｜ｚ_ｔ）の対数オッズである。センサーによって得られた格子ｉにおけるｔ時間フレーム（時刻ｔ）の観測結果は、この右辺第１項に代入される。なお、式（１）の右辺第１項は、時刻ｔの観測における対数尤度比と呼ばれる。

【0010】

また、式（１）の右辺第２項は初期値であり、初期値として設定された格子ｉにおける占有確率ｐ（ｍ_ｉ）の対数オッズである。なお、観測が全く行われていない初期状態における占有確率ｐ（ｍ_ｉ）は、例えばｐ（ｍ_ｉ）＝０．５とすることができる。

【0011】

また、式（１）の右辺第３項は、格子ｉにおいて、ｔ−１時間フレーム（時刻ｔ−１）までに得られた占有確率ｐ（ｍ_ｉ｜ｚ_１，．．．，ｚ_ｔ−１）の対数オッズであり、直前の計算によって得られた推定結果（直前の計算で得られた式（１）の左辺）の値が、この右辺第３項に代入される。なお、式（１）の右辺第３項は、時刻ｔ−１までの観測における対数尤度比と呼ばれる。

【0012】

上記のように、初期状態における占有確率ｐ（ｍ_ｉ）をｐ（ｍ_ｉ）＝０．５とした場合、非占有確率１−ｐ（ｍ_ｉ）も１−ｐ（ｍ_ｉ）＝０．５となるため、式（１）の右辺第２項はゼロとなる。すなわち、初期値としてｐ（ｍ_ｉ）＝０．５とした場合、単純に、時刻ｔ−１までの観測における対数尤度比の値に、時刻ｔの観測における対数尤度比の値を足すだけで、時刻ｔまでの観測における対数尤度比の値を得ることができる。この更新式を繰り返して計算することから、占有確率ｐ（ｍ_ｉ）の更新式を、下記の式（２）のような総和で表すことができる。

【0013】

【数2】

【0014】

上記の式（２）は時刻ｔごとに１回観測されるとして計算されていることから、観測回数をパラメータとして用いた場合も同様の式によって表すことができる。すなわち、観測回数のパラメータをｋとし、Ｋ回の観測が行われた場合も同様に式（１）をＫ回繰り返して計算すると考え、占有確率ｐ（ｍ_ｉ）の更新式を、下記の式（３）のような総和で表すことができる。

【0015】

【数3】

【0016】

なお、式（３）の左辺は、格子ｉにおいて、Ｋ回の観測回数までに得られた占有確率ｐ（ｍ_ｉ｜ｚ_１，．．．，ｚ_Ｋ）の対数オッズであり、すなわち、Ｋ回の観測回数までの観測における対数尤度比である。

【0017】

また、下記の特許文献１には、格子占有地図を用いて、周辺の物体を認識及び推定する技術が開示されている。特許文献１の開示技術では、格子占有地図上の静的死角領域と、現在の検出位置に基づく現在死角領域との差分に基づいて、移動している物体によって形成された動的死角領域の推定が行われる。

【先行技術文献】

【特許文献】

【0018】

【特許文献1】特開２０１１−１２３５５１号公報（要約書、段落［００４６］）

【非特許文献】

【0019】

【非特許文献1】“Using Occupancy Grids for Mobile Robot Perception and Navigation”, Alberto Elfes, Carnegie Mellon University, IEEE 0018-9162/89/0600-0046 1989

【非特許文献2】“Learning Occupancy Grid Maps with Forward Sensor Models”, Sebastian Thrun, Autonomous Robots - Springer Journals, DOI: 10.1023/A:1025584807625, Published: Sep 1, 2003

【非特許文献3】“Optimum Array Processing: Part IV of Detection, Estimation and Modulation Theory”, Wiley-Interscience, ISBN 0-471-09390-4, 2002

【発明の概要】

【発明が解決しようとする課題】

【0020】

従来の技術におけるセンサーモデルでは、センサーから伝搬されるレーダーの波を直進性の強いものと仮定し、センサーから照射され物体に反射して戻ってきた波を、あたかもボールを真っ直ぐ投げて何かにぶつかって返ってきたかのように考える。レーダーの波がぶつかって反射した地点には何らかの物体があると認識して、例えば図１に示されているように、その地点に対応する占有格子地図上の格子を占有状態（図１内の黒い格子）とし、また、センサーと何らかの物体が認識された地点とを結ぶ直線上に存在する格子は、レーダーの波が通った経路と見なして、無の空間、すなわち非占有状態（図１内の白い格子）とする。

【0021】

例えば、上述した特許文献１の開示技術は、ある物体が認識された場合、自車両からその物体までの直線上には他の物体は存在しないという前提に基づき、自車両から物体までの直線上に存在するブロック（格子）における物体の存在確率を減少させる処理が行われる。

【0022】

しかしながら、レーダーの波は、ボールを投げて戻ってくる場合と大きく異なり、散乱、回折、反射などの波の性質を有するため、様々な経路（いわゆるマルチパス）を通り得る。その結果、物体の裏側に回り込んだり、様々な物体において散乱や反射を繰り返したりする場合もある。したがって、センサーと何らかの物体が認識された地点とを結ぶ直線上に別の物体が存在する場合であっても、この別の物体を認識できないことがあり、センサーと認識された物体とを単純に直線で結んで、その直線上を非占有状態であると断定することはできない。

【0023】

また、従来の手法では、観測データ空間内のピクセルのうち、例えば、事前に決められたある閾値より強度の高いものだけが選別され、選別されたピクセルと１対１に対応付けられる物理的空間の格子だけについて、上述した式（１）などを用いた更新計算が行われている。しかしながら、実際の観測データは、上述したように観測環境の様々な要因によって不確実であり、こうした不確実性を包含した新たな手法が必要である。

【0024】

上記の問題に鑑み、本発明は、新たなセンサーモデルを提案し、センシング機器を用いて得られた観測データから、より高い精度で周辺の環境を推定できるようにする周辺環境推定装置及び周辺環境推定方法を提供することを目的とする。

【課題を解決するための手段】

【0025】

上記の目的を達成するため、本発明によれば、センシング信号を送信する送信機と、前記センシング信号が反射して戻ってくる信号を受信する受信機とを有するセンシング機器を用いて得られる観測データから、周辺の環境の物体の存在を推定する周辺環境推定装置であって、
前記観測データに基づいて前記受信器で受信した信号の強度値を計算し、前記強度値から前記受信器で受信した信号が被検出物から反射してきたものであるかどうかを表す尤度を計算する強度値尤度計算部と、
前記観測データに基づいて前記被検出物の方位角を計算し、前記強度値から前記方位角を中心とした分散の程度を表す角度の分布関数を計算する分布関数計算部と、
前記周辺の空間を格子状に分割して各格子における物体の占有確率を記述する占有格子地図に関して、前記尤度を前記分布関数によって分散させた値を用いて、前記占有格子地図内の複数の格子における前記物体の占有確率を更新する計算を行う占有格子地図計算部とを、
有し、
前記複数の観測データについて、前記強度値尤度計算部による計算、前記分布関数計算部による計算、並びに、前記占有格子地図計算部による計算を行うことで、前記複数の観測データに基づく前記占有格子地図を作成するよう構成されている周辺環境推定装置が提供される。

【0026】

また、上記の目的を達成するため、本発明によれば、センシング信号を送信する送信機と、前記センシング信号が反射して戻ってくる信号を受信する受信機とを有するセンシング機器を用いて得られる観測データから、周辺の環境の物体の存在を推定する周辺環境推定方法であって、
前記観測データに基づいて前記受信機で受信した信号の強度値を計算し、前記強度値から前記受信機で受信した信号が被検出物から反射してきたものであるかどうかを表す尤度を計算する強度値尤度計算ステップと、
前記観測データに基づいて前記被検出物の方位角を計算し、前記強度値から前記方位角を中心とした分散の程度を表す角度の分布関数を計算する分布関数計算ステップと、
前記周辺の空間を格子状に分割して各格子における物体の占有確率を記述する占有格子地図に関して、前記尤度を前記分布関数によって分散させた値を用いて、前記占有格子地図内の複数の格子における前記物体の占有確率を更新する計算を行う占有格子地図計算ステップとを、
有し、
前記複数の観測データについて、前記強度値尤度計算ステップによる計算、前記分布関数計算ステップによる計算、並びに、前記占有格子地図計算ステップによる計算を行うことで、前記複数の観測データに基づく前記占有格子地図を作成する周辺環境推定方法が提供される。

【発明の効果】

【0027】

本発明は、センシング機器を用いて得られた観測データから、より高い精度で周辺の環境を推定できるようになるという効果を有する。また、特に、センシング機器が車両に搭載されている場合には、より高い精度で車両周辺の環境を推定できるようになるという効果を有する。

【図面の簡単な説明】

【0028】

【図1】従来の技術における占有格子地図の一例を示す図である。

【図2】本発明の第１の実施の形態における周辺環境推定装置の構成の一例を示す図である。

【図3】本発明の第１の実施の形態において用いられる、センサーのアンテナ構成の一例を示す図である。

【図4】（ａ）は、センサーを用いたセンシングによって得られた観測データがマッピングされる観測データ空間を模式的に示す図であり、（ｂ）は、観測データに基づいて作成される占有格子地図が配置される物理的空間を模式的に示す図である。

【図5】（ａ）は、占有格子地図が配置される２次元空間を示す図であり、（ｂ）は、占有格子地図が配置される３次元空間を示す図である。

【図6】本発明の第１の実施の形態における処理の一例を示すフローチャートである。

【図7】（ａ）は、本発明の第１の実施の形態において用いられる、物体及びノイズに関する強度の確率分布グラフの一例であり、（ｂ）は、（ａ）に示すグラフをより詳細に説明するための部分拡大図である。

【図8】アレイセンサーにおいて観測された信号のＳＮＲと物体の推定角度の分散の変化を表す図である。

【図9】（ａ）は、図８を用いて説明した分散の様子に関して、角度方向を示す横軸を、角度方向に配置された格子の配列と重ね合わせて示す図であり、（ｂ）及び（ｃ）は、物体の存在する確率が広がりを持つことを考慮して、ピクセルｊが複数の格子ｉにマッピングされる様子を模式的に示す図である。

【図10】（ａ）は、センサーが搭載された車両から前方を見た風景を示す図であり、（ｂ）は、このときセンサーによってセンシングした観測データから、静止物のみを抽出して作成された占有格子地図の一例を示す図である。

【図11】（ａ）は、センサーが搭載された車両から前方を見た風景を示す図であり、（ｂ）は、このときセンサーによってセンシングした観測データから、移動体のみを抽出して作成された占有格子地図の一例を示す図である。

【図12】（ａ）及び（ｂ）は、本発明の第１の実施の形態における手法を用いて実際に占有格子地図を作成した例を示す図である。

【図13】本発明の第３の実施の形態における周辺環境推定装置の構成の一例を示す図である。

【図14】本発明の第３の実施の形態における特徴的な処理の一例を示すフローチャートである。

【図15】（ａ）及び（ｂ）は、１回の観測データを含んだローカル占有格子地図がマッピングされたグローバル占有格子地図の一例を示す図であり、（ａ）及び（ｂ）は、複数のローカル占有格子地図が組み合わされたグローバル占有格子地図の一例を示す図である。

【図16】本発明の第３の実施の形態における手法に基づいて実際に作成されたグローバル占有格子地図の一例を示す図である。

【発明を実施するための形態】

【0029】

以下、図面を参照しながら、本発明の実施の形態について説明する。

【0030】

＜第１の実施の形態＞
以下、本発明の第１の実施の形態について説明する。まず、図２を参照しながら、本発明の第１の実施の形態における周辺環境推定装置の構成について説明する。図２には、本発明の第１の実施の形態における周辺環境推定装置の構成の一例が図示されている。図２に示す周辺環境推定装置１０は、観測データ空間処理部１０１、強度値尤度計算部１０３、分布関数計算部１０５、占有格子地図計算部１０７を有する。

【0031】

観測データ空間処理部１０１は、センサー（センシング機器）５０から得られた観測データを、観測データ空間上のデータとして扱うことができるよう処理を行う機能を有する。すなわち、観測データ空間処理部１０１は、センサー５０から得られた観測データに基づき、観測データ空間を構築する。

【0032】

強度値尤度計算部１０３は、観測データ空間に存在する特定のピクセルについて、ピクセルに関連する強度値を取得し、その強度値に対応した尤度比を計算する機能を有する。強度値尤度計算部１０３では、例えば図６に示すステップＳ２０５及びＳ２０７における処理が行われるが、詳細については後述する。

【0033】

分布関数計算部１０５は、観測データ空間に存在する特定のピクセルについて、ピクセルに関連する位相差を取得し、その位相差の分布関数を計算する機能を有する。分布関数計算部１０５では、例えば図６に示すステップＳ２０９及びＳ２１１における処理が行われるが、詳細については後述する。

【0034】

占有格子地図計算部１０７は、強度値尤度計算部１０３によって計算された尤度比と、分布関数計算部１０５によって計算された位相差の分布関数とを用いて、センサー５０の周辺の物理的空間に対応する占有格子地図を計算する機能を有する。占有格子地図計算部１０７では、例えば図６に示すステップＳ２１３における処理が行われるが、詳細については後述する。占有格子地図計算部１０７によって計算された占有格子地図は、占有格子地図データとして出力され、任意の記録媒体に格納されてもよく、他の任意の装置に供給され利用されてもよい。

【0035】

本発明の第１の実施の形態では、観測データを得るためのセンサー５０は、例えば図３に示すように、１つの送信アンテナ１１と、複数の受信アンテナ１２とを有するアンテナ１０を備えたものが用いられる。送信アンテナ１１からは、センシング信号を含む波（電波）が照射され、センサー５０の周辺に存在する物体に反射して戻ってきた波は、複数の受信アンテナ１２でセンシングされる。このとき、複数の受信アンテナ１２が、アレイ型（受信位相アレイ）に配置されていることによって、物体の方位を検知することが可能となる。

【0036】

なお、本発明の第１の実施の形態で用いられるアンテナ１０は、図３に示すように受信位相アレイなどの構成を有するアンテナに限定されるものではなく、物体の方位を検知できる機能を有する他の構成のアンテナであってもよい。また、図３では、１つの送信アンテナ１１を有する構成のアンテナ１０が図示されているが、２つ以上の送信アンテナ１１を有していてもよい。

【0037】

ここで、図４（ａ）及び図４（ｂ）を参照して、センサー５０のセンシングによって得られた観測データが配置される観測データ空間、及び、観測データに基づいて作成される占有格子地図が配置される物理的空間（センサーモデルとなるセンサー５０の周辺の物理的空間に対応）について説明する。

【0038】

センシングによって得られた観測データは、観測データ空間として表現される。検知されたレーダーの波の位相から「距離−ドップラー」（すなわち、物体までの距離と、物体との間のドップラー速度）の２次元からなる複素数データが取得される。その様子を図４（ａ）に示す。

【0039】

図４（ａ）には、センサー５０のセンシングによって得られた観測データがマッピングされる観測データ空間が模式的に示されている。図４（ａ）に示す観測データ空間は、縦軸を距離（縦軸：例えば、範囲１２８）とし、横軸をドップラー速度（例えば、範囲１２８）とする２次元空間により構成される。この２次元空間は、特定の距離及び特定のドップラー速度に対応する『ピクセル』により構成され、各ピクセルには「距離−ドップラー」の２次元からなる複素数の値が入る。なお、複素数の絶対値の２乗が、強度（信号の強度）となる。

【0040】

また、図４（ａ）では、アンテナ１０の受信位相アレイに対応して、複数の２次元空間（図３に示すアンテナ１０の構成に合わせて、図４（ａ）では６つの２次元空間）が図示されている。受信位相アレイを構成する複数の受信アンテナ１２のそれぞれにおいて観測データがセンシングされることにより、それらの観測データの位相差から方位角φを推定することができる。

【0041】

一方、センサーモデルとなるセンサー５０の周辺の物理的空間を図４（ｂ）に示す。センサー５０の周辺の物理的空間は、例えば図４（ｂ）に示すように、センサー５０（図４（ｂ）では、車両に設置されていることを想定している）の前方に、距離ｒと方位角φからなる円座標で表現されている。この円座標からなる２次元空間は、特定の距離及び特定の方位角に対応する『格子』により構成される。したがって、各格子は「距離−方位角」の２次元の値を表している。

【0042】

なお、混乱を避けるため、本明細書では、図４（ａ）に示すような観測データ空間における微小領域（単位領域）を『ピクセル』と呼び、図４（ｂ）に示すような物理的空間における微小領域（単位領域）を『格子』又は『占有格子』と呼ぶことで、両者の区別を明確にする。すなわち、ピクセルは観測データ空間を区切った離散的な微小領域であり、格子は物理的空間を区切った離散的な微小領域である。また、本明細書では、ピクセルのインデックスには文字『ｊ』を用い、格子又は占有格子のインデックスには文字『ｉ』を用いる。

【0043】

また、占有格子地図は、図５（ａ）に示すように、例えば円座標（ｒ，φ）で表される２次元占有格子（２次元空間）によって構成されてもよく、図５（ｂ）に示すように、例えば球座標（ｒ，φ，θ）で表される３次元占有格子（３次元空間）によって構成されてもよい。以下、本明細書では、本発明をより良く理解できるようにするため、円座標（ｒ，φ）によって表される２次元占有格子（２次元空間）を持つ占有格子地図を一例に挙げて説明を行うが、これに限定されるものではない。

【0044】

さらに、本明細書では、説明を簡略化するため、作成する占有格子地図を構成する円座標の２次元空間の視野角の範囲を約９０度とする。また、円座標の２次元空間全体の形は、センサー５０によってセンシング可能な領域に合わせて、例えば、センシング距離の範囲を６０ｍまで、センシング方位角の範囲を−４０°から４０°とし、図４（ｂ）に示す物理的空間の格子も、占有格子地図を構成する円座標の２次元空間に合わせたサイズの扇形にしてフィットさせる。本発明の第１の実施の形態では、このような観測データ空間及び物理的空間の構成によって、図４（ａ）に示す観測データ空間内のある１ピクセルのデータが、物理的空間内ではどのような確率分布で複数の格子にまたがって存在しているのかを反映できるようにする。

【0045】

図４（ａ）に示す観測データ空間について、センサーモデルとなるセンサー５０の周辺の物理的空間との対応をまとめると、以下のようになる。

【0046】

図４（ａ）に示す観測データ空間の「距離−ドップラー」の２次元の離散データの各ピクセルを占める複素数に関して、その複素数の絶対値の２乗が強度となり、強度値の大きいものほど反射の強い反射点と考えられる。

【0047】

図４（ａ）に示す観測データ空間の「距離−ドップラー」の縦軸は、観測地点（センサー５０の位置）から反射点までの距離であり、例えば図４（ｂ）の扇形における、車両前方からの距離に相当し、０ｍから６０ｍの距離の範囲となる。

【0048】

図４（ａ）に示す観測データ空間の「距離−ドップラー」の横軸は、ドップラー効果によるセンサーから反射点までを結んだ方向に沿った速度に関する値となる。この速度は、観測地点の速度（例えば、センサー５０が搭載された車両の速度）と、反射点に存在する物体との相対速度に相当する。

【0049】

図４（ａ）に示す観測データ空間の「距離−ドップラー」の２次元空間は、受信位相アレイを構成する受信アンテナ１２の数の分だけ存在しており、各２次元空間において同一の位置に対応するピクセル同士の複素数の値の位相差から、そのピクセルに相当する反射点を観測地点から見た方位角φを推定できる。

【0050】

本発明の第１の実施の形態では、観測データ空間内のある１つのピクセルの観測データを物理的空間内のある格子に１対１に対応させて計算を行うのではなく、観測データ空間内のある１つのピクセルが複数の格子に影響を及ぼすように、複数の格子にまたがって分布させて計算することに特徴がある。このとき、観測データ空間内のある１つのピクセルの強度が弱いほど、物理的空間での対応する方位角φの推定値の分散が広がっていくように、つまり物理的空間内のある格子との対応があいまいとなっていくように、複数の格子にまたがって分布させる。

【0051】

従来の手法では、観測データ空間内のピクセルのうち、例えば、事前に決められたある閾値より強度の高いものだけを選別し、選別されたピクセルと１対１に対応付けられる物理的空間の格子だけについて、上述した式（１）を用いた更新計算を行っていた。

【0052】

それに対し、本発明では、観測データ空間内の各ピクセルの強度を考慮した選別は行わず、すべてのピクセルについて計算を行うようにする。本発明では、あるピクセルｊに対応する物理的空間の格子ｉを１つに絞るのではなく、物理的空間に存在する確率分布に含まれる領域のすべての格子群に対して、本発明によって提案される数式を用いて更新計算を行う。

【0053】

以下に、本発明の第１の実施の形態における具体的な処理の手順について説明する。図６は、本発明の第１の実施の形態における処理の一例を示すフローチャートであり、図２に示す周辺環境推定装置１０において処理が行われる場合を想定して説明する。以下では、まず、図６に基づいて、本発明の第１の実施の形態における処理の概略について説明し、その後、それぞれのステップの処理や計算に用いられる数式などについて詳細に説明する。

【0054】

図６において、観測データ空間処理部１０１は、センサー５０によってセンシングされた観測データを受け取り、観測データ空間を構築する（ステップＳ２０１）。具体的には、センサー５０のアンテナ１０に備えられた複数の受信アンテナ１２のそれぞれから受け取った観測データを、図４（ａ）に示すような観測データ空間上に配置し、各ピクセルに関連付けられた「距離−ドップラー」の複素数の値を取得する。

【0055】

次に、周辺環境推定装置１０は、観測データ空間からあるピクセルをサンプリングする（ステップＳ２０３）。サンプリングするピクセルは、観測データ空間全体（あるいは、観測データ空間のうちの関心のある空間全体）内の任意のピクセルであり、後述するステップＳ２１５の処理と合わせて、観測データ空間全体のピクセルを順次サンプリングする。すなわち、周辺環境推定装置１０は、観測データ空間全体の各ピクセルに対して、ステップＳ２０５とステップＳ２１５との間の処理（ステップＳ２０５〜Ｓ２１３の処理）を順次行う。

【0056】

ステップＳ２０３でサンプリングしたピクセルに含まれているデータ（複素数の値）は、強度値尤度計算部１０３及び分布関数計算部１０５のそれぞれにおいて処理される。ステップＳ２０５、Ｓ２０７は、強度値尤度計算部１０３で行われる処理であり、ステップＳ２０９、Ｓ２１１は、分布関数計算部１０５で行われる処理であり、ステップＳ２１３は、占有格子地図計算部１０７で行われる処理である。

【0057】

強度値尤度計算部１０３は、サンプリングしたピクセルの強度値を得る（ステップＳ２０５）。例えば、サンプリングしたピクセルに含まれる複素数の絶対値の２乗をとることで受信信号の強度値を得ることができ、ここではさらに、その強度値の対数をとりデシベル（ｄＢ）に変換する。そして、強度値尤度計算部１０３は、ステップＳ２０５で計算した強度値（ｄＢ）と、物体の反射強度及びノイズの強度の尤度分布表とに基づいて、強度値（ｄＢ）に対応した対数尤度比を計算する（ステップＳ２０７）。

【0058】

一方、分布関数計算部１０５は、サンプリングしたピクセルの位相差から方位角φを得る（ステップＳ２０９）。例えば、サンプリングしたピクセルには、アレイ型の複数の受信アンテナ１２の観測データが関連付けられており、これらの観測データ同士の位相差から、ピクセルｊに対応している物理的空間上の方位角φを推定することができる。そして、分布関数計算部１０５は、方位角φに対する分散の程度を表す分布関数を計算し、この分布関数から重みＷ（Φ_ｉ，ｊ）を計算する（ステップＳ２１１）。後述するように、分散の程度は、例えば、強度値（ｄＢ）とＣＲＢ（クラメール・ラオ限界：Cramer Rao Bound）から求めることができ、方位角φ（すなわち、格子ｉの左右方向の配列）をパラメータとする分布関数を規格化することで重みとすることができる。

【0059】

そして、占有格子地図計算部１０７は、強度値尤度計算部１０３による計算結果（対数尤度比）、分布関数計算部１０５による計算結果（重み）を利用して、複数の格子に対して数式を用いて更新計算を行う（ステップＳ２１３）。後述するように、このとき更新される対象は、単独の格子ではなく、複数の格子であり、例えばサンプリングしたピクセルの距離に対応する格子群である。

【0060】

上記のステップＳ２０５〜Ｓ２１３の処理は、観測データ空間全体（あるいは、観測データ空間のうちの関心のある空間全体）内のすべてのピクセルに関して実行される（ステップＳ２１５）。そして、すべてのピクセルのサンプリング及び計算処理が完了したとき、各格子ｉに占有確率の値が関連付けられた占有格子地図が作成される。

【0061】

以下、ステップＳ２０７における対数尤度比の計算（強度値尤度計算部１０３における計算処理）について説明する。センサー５０における受信信号は、物体から反射してきた信号に加えて、何らかのノイズの影響による信号も存在する。受信信号の強度値から、物体の反射によるものか、ただのノイズなのかを判別するため、本発明の第１の実施の形態では強度値（ｄＢ）を確率変数とみなす。すなわち、観測データ空間内のあるピクセルの強度値（ｄＢ）を考えた場合、物体であった場合（占有状態ｍである）の強度（ｄＢ）の確率分布と、ノイズであった場合（占有状態ｍではない）の強度（ｄＢ）の確率分布を考え、これらの確率分布に基づいて、物体又はノイズのどちらから得られた強度値（ｄＢ）であるかを示す確率変数を計算する。それを図７（ａ）及び（ｂ）に示す。

【0062】

図７（ａ）には、本発明の第１の実施の形態において用いられる、物体であった場合（占有状態ｍである）の強度（ｄＢ）の確率分布のグラフの一例と、ノイズであった場合（占有状態ｍではない）の強度（ｄＢ）の確率分布のグラフの一例とが示されている。これらの分布は、本発明の第１の実施の形態における処理の前に得られるものである。例えば、物体の確率（以下、ターゲットの確率と記載することもある）ついては、既知の環境下で、物体から反射してきた場合の強度分布を収集することによって、物体から反射してきた場合の強度の確率分布（ターゲットの分布）を事前に求めることができる。なお、物体から反射してきた場合の強度は物体の材質や形状などに依存するが、ここでは、ノイズの強度未満の値をとることはない、ノイズの平均強度以上の強度をとる場合は物体の強度は一定の確率になる、として、図７（ａ）に示すターゲットの分布を求めている。

【0063】

また、ノイズの確率については、周辺に物体が何も存在しない場所でのレーダーからの距離ｒにおける強度分布を収集することによって、ノイズの確率分布（ノイズの分布）を事前に求めることができる。ここでは、ガウス分布であると仮定している。

【0064】

なお、上記のターゲットの分布及びノイズの分布は、環境や使うセンサー５０などによって異なり、設計者が自由に設定可能である。また、本発明ではその分布の形状については限定しない。また、ターゲットの分布及びノイズの分布は表形式（数値の対応が事前に設定されたデータベースなど）で保持されてもよく、あるいは、関数（入力された数値に基づいて計算を行い、その計算結果を出力する）として保持されていてもよい。

【0065】

図７（ａ）に示すグラフの３０〜４０（ｄＢ）付近を拡大した図７（ｂ）を参照しながら、より詳細に説明する。例えば、ある観測データのピクセルｊから強度ｚ_ｊ（ｄＢ）が得られた場合、この強度ｚ_ｊが図７（ａ）及び図７（ｂ）の横軸の値となり、図７（ｂ）に示すように、この強度ｚ_ｊに対応するターゲットの確率及びノイズの確率が得られる。そして、この強度ｚ_ｊに対応するターゲットの確率及びノイズの確率から、物理的空間の格子内に物体が存在するかどうか（すなわち、受信信号が被検出物から反射してきたものであるか、あるいはノイズに由来するものであるか）を表す占有確率ｐ（ｍ｜ｚ_ｊ）を求めることができる。さらに、非占有確率１−ｐ（ｍ｜ｚ_ｊ）についても求められ、下記の式（４）で表されるような占有確率の対数オッズを対数尤度比として定義する。

【0066】

【数4】

【0067】

なお、式（４）に含まれる状態ｍにはインデックスｉがなく、ある単独のピクセルｊに対して、どの占有格子ｉが更新されるべきか特定されていない。上記の理由から、式（４）で表される対数尤度比を、格子不定の対数尤度比と呼ぶことにする。

【0068】

次に、ステップＳ２１１における重みの計算（分布関数計算部１０５における計算処理）について説明する。上述のように、本発明の第１の実施の形態では、観測データ空間内のある１つのピクセルに含まれる観測データを、複数の格子にまたがって分布させて計算することに特徴がある。すなわち、式（４）で定義された格子不定の対数尤度比は、複数の格子の更新計算に用いられる点に特徴がある。例えば、ある単独のピクセルｊの強度から式（４）を用いて得られた格子不定の対数尤度比は、そのピクセルｊに対応する物理的空間の格子と等しい距離ｒを持つが異なる角度φをとる複数の格子群（例えば、図４（ｂ）において、車両前方のセンサー５０からの距離が等しい複数の格子）の更新計算に用いられる。

【0069】

前述のように、観測データ空間内のあるピクセルjは、物理的空間内においてどの方位角φに対応しているかを、受信位相アレイを構成する複数の受信アンテナ１２同士の位相差によって推定できる。しかしながら、ピクセルjからサンプリングした強度が弱い場合、一般に、推定された方位角φの分散が広がることが知られている。それはつまり、受信アンテナ１２同士の位相差から推定した方位角φ_{ｃｅｎｔｅｒ}を中心として、その左右方向（円周方向）の複数の格子群にまたがって不確定となることを意味する。

【0070】

言い換えれば、従来の技術では、強度に依存して推定位置が不確定となることから、所定の閾値よりも強度の強いピクセルのみをサンプリングし、方位角φを１つの格子に限定できるように計算していた。しかしながら、本発明では、強度の強弱にかかわらず、すべてのピクセルをサンプリングし、広がった分布をもって物理的空間内の複数の格子を更新していく。

【0071】

ここで、非特許文献３（ページ９９４、Figure 8.19を参照）に記載されている、アレイセンサーにおいて観測された信号のＳＮＲ（信号対雑音比：signal-to-noise ratio）と物体の推定角度の分散の変化を表す図を図８に示す。これは、一般にＣＲＢ（クラメール・ラオ限界：Cramer Rao Bound）と呼ばれる。

【0072】

図８には、縦軸を角度方向とし横軸をＳＮＲ値とした座標系が示されている。角度方向を示す縦軸は、中心をゼロ角度であり、上下プラスマイナスに広がりを持つ。中心であるゼロ角度は正しい方位角であり、ゼロ角度から上下方向に離れれば離れるほど、方位角に大きな誤差があることを表している。また、ＳＮＲ値を示す横軸は、左から右にＳＮＲ値が大きい状態を表している。すなわち、横軸方向の左側ほどＳＮＲが弱く、観測された信号の強度が弱いことを表し、右側ほどＳＮＲが強く、観測された信号の強度が強いことを表している。

【0073】

さらに、図８には、上記の座標系に観測データをプロットした状態が示されている。観測データは、どの程度の強度を持ち（ＳＮＲ値）、かつ、どの角度方向を表しているか（方位角）に応じて、図８に示す座標系に『＊』でプロットされている。すなわち、図８中に存在する『＊』は、１つ１つが個々の観測結果を表しており、具体的には、ＳＮＲ値（ｄＢ）が−１２、−９、−６、−３、０、３、６の観測データがそれぞれ１００個ずつプロットされている。

【0074】

図８を参照すると理解できるように、ＳＮＲが強いものほど、角度方向の分布幅が狭く、推定角度φが正しい角度（ゼロ角度）に分散なく分布する状態となる。また、ＳＮＲが弱いものほど、角度方向の分布幅が広く、推定角度φが正しい角度（ゼロ角度）から離れた分散された状態となる。すなわち、これは、ＳＮＲが弱いほど、推定角度φからずれた場所に物体が存在する可能性が大きいことを示している。なお、図８には、ＳＮＲ値（ｄＢ）が−６の場合（分布幅が広い状態）、ＳＮＲ値（ｄＢ）が０の場合（分布幅が狭い状態）のそれぞれに関して、分散の様子をより明確に表す分布曲線が模式的に図示されている。

【0075】

このように、強度が弱い観測ほど、確率的に本来の方位角（正しい方位角）からプラス又はマイナスにずれた角度となる場合が多いことが分かる。このことは、観測データ空間内の１つのピクセルｊの観測データに対し、そのピクセルｊの観測データの影響を受ける格子ｉは１つに限定されるものではなく、強度が弱いほど、強度が左右方向（方位角φの周方向）に広がって、分布の広がりが大きくなることを示している。その様子を図９（ａ）に示す。

【0076】

図９（ａ）は、図８を用いて説明した分散の様子に関して、角度方向を示す横軸を、角度方向に配置された格子の配列と重ね合わせて図示したものである。図９（ａ）では、受信アンテナ１２同士の位相差から推定した方位角φ_{ｃｅｎｔｅｒ}（図８のゼロ角度に相当）を中心として、物体の存在する確率が広がりを持つことが示されている。また、図９（ａ）では、物体の存在する確率が複数の格子ｉにわたって広がっていることが、グラデーションによって示されている。なお、グラデーションが濃いほど、物体の存在する確率が大きいことを表している。このように、確率の広がりにより、方位角φ_{ｃｅｎｔｅｒ}から離れた位置に存在する格子ｉにおいても物体の存在する確率はゼロではない。

【0077】

また、図９（ｂ）及び図９（ｃ）には、物体の存在する確率が広がりを持つことを考慮して、観測データ空間内のピクセルｊが物理的空間内の複数の格子ｉにマッピングされる様子が模式的に示されている。図９（ｃ）においても、受信アンテナ１２同士の位相差から推定した方位角φ_{ｃｅｎｔｅｒ}を中心として、物体の存在する確率が広がりを持つことが示されている。また、物体の存在する確率が複数の格子ｉにわたって広がっていることが、グラデーションによって示されている。

【0078】

以上の知見に基づき、ステップＳ２１１における重みの計算では、サンプリングしたピクセルｊに対し、受信アンテナ１２同士の位相差から方位角φ_{ｃｅｎｔｅｒ}を計算し、その方位角φ_{ｃｅｎｔｅｒ}を中心に左右に離れるにつれて存在確率がどのように分散するかを示す分布関数を、ピクセルｊの強度からＣＲＢを使って計算する。さらに、このようにして得られた分布関数を規格化した関数Ｗ（Φ_ｉ，ｊ）を計算する。なお、後述するように、この関数Ｗ（Φ_ｉ，ｊ）は、対数尤度比に対する重みとして用いられることから、重みＷ（Φ_ｉ，ｊ）と呼ぶ。

【0079】

次に、ステップＳ２１３における計算（占有格子地図計算部１０７における計算処理）について説明する。上述したように、重みＷ（Φ_ｉ，ｊ）は、あるピクセルｊにおける強度が弱いほど、左右方向（方位角φの周方向）に広がって分布の広がりが大きくなることを表す規格化された分布関数である。この重みＷ（Φ_ｉ，ｊ）によって表現された分散の考えを、式（４）で表される格子不定の対数尤度比に適用し、さらに、観測による更新を加算する従来の更新式（式（１）〜式（３）の式を参照）の概念を拡張して、下記の式（５）で表される更新式を定義する。

【0080】

【数5】

【0081】

なお、式（５）は、サンプリングしたピクセルｊに関して、そのピクセルｊに対応する物理的空間の格子と等しい距離ｒを持つすべての格子群に対して適用される。これは、重みＷ（Φ_ｉ，ｊ）が、格子ｉの左右方向（周方向）の分散であることを考慮している。

【0082】

式（５）の左辺は、格子ｉにおいて、ピクセルｊまでに得られた占有確率ｐ（ｍ_ｉ｜ｚ_１，．．．，ｚ_ｊ）の対数オッズである。一方、式（５）の右辺第３項は、式（５）の左辺は、格子ｉにおいて、ピクセルｊ−１までに得られた占有確率ｐ（ｍ_ｉ｜ｚ_１，．．．，ｚ_ｊ−１）の対数オッズであり、直前の計算によって得られた計算結果（直前の計算で得られた式（５）の左辺）の値が、この右辺第３項に代入される。

【0083】

すなわち、例えば図６に示す処理では、観測データ空間内からピクセルを順次サンプリングして計算処理を行う際、１つ前にサンプリングしたピクセル（例えば、ピクセルｊ−１）における式（５）の計算結果を記憶しておき、次にサンプリングするピクセル（例えば、ピクセルｊ−１）における式（５）の左辺第３項に代入する処理が行われる。

【0084】

なお、式（５）の右辺第２項は任意に設定される初期値であり、例えばｐ（ｍ_ｉ）＝０．５とした場合には、右辺第２項はゼロとなる。したがって、最も単純化して右辺第２項をゼロに設定すると、式（５）は、式（５）の右辺第１項の総和として表される。

【0085】

式（５）の右辺第１項は、本発明の特徴を表す項であり、格子不定の対数尤度比に対して、格子ｉごとに異なる重みＷ（Φ_ｉ，ｊ）を掛けた項である。この重みＷ（Φ_ｉ，ｊ）は、図９（ａ）又は図９（ｃ）に模式的に図示されている分散を示した分布のある角度φに対する密度とも言える。

【0086】

観測データから求められた対数尤度比に対して、この分布の密度を重みとした積を考えればよいことは、従来の技術における式（３）との対比から、以下のように導き出せる。すなわち、図８での分布が、個別の観測（図８中の『＊』）が複数回行われた頻度数によって表現されていることに着目する。ここで、観測データ空間内のあるピクセルｊからの影響が、式（３）で表現されているように、多数のＫ回の観測『＊』として、Ｗ（Φ_ｉ，ｊ）という分布をなして、複数の格子にまたがっていると考える。ある格子ｉに対して、あるピクセルｊのうちの、重みＷ（Φ_ｉ，ｊ）に比例した観測『＊』の数の分だけ足されたものが更新分の値となる。最終的に、Ｋ個を無限数と考えて全体を規格化することは、重みＷ（Φ_ｉ，ｊ）を格子不定の対数尤度比（式（４））に掛けることに等しい。

【0087】

以上の手法により、観測データ空間内の各ピクセルの値をサンプリングして、物理的空間内の複数の格子における占有確率を更新していくことで、センサー５０によってセンシングした観測データを占有格子に反映し、占有格子地図を作成することが可能となる。

【0088】

また、図１０（ａ）及び（ｂ）並びに図１１（ａ）及び（ｂ）に、上述した手法を用いて実際に占有格子地図を作成した例を示す。図１０（ａ）及び図１１（ａ）には、センサー５０が搭載された車両から前方を見た風景が示されており、図１０（ｂ）及び図１１（ｂ）には、このときセンサー５０によってセンシングした観測データから作成された占有格子地図が示されている。

【0089】

図１０（ｂ）は、静止しているとみなすことができる反射点を対象として計算を行った結果であり、図１０（ａ）に写っている建物や停止している車両が精度良く検知されている。また、図１１（ｂ）は、移動しているとみなすことができる反射点を対象として計算を行った結果であり、図１１（ａ）に写っている先行車両が精度良く検知されている。

【0090】

以上説明したように、本発明の第１の実施の形態によれば、確率変数に基づく新たなセンサーモデルが提案され、観測データから得られる強度の強弱に応じて占有確率（格子不定の対数尤度比）を計算するとともに、強度の強弱に応じて異なる分散の幅（分布関数）で複数の格子にわたって影響を与えることを考慮して、占有格子地図を作成することが可能となる。

【0091】

＜第２の実施の形態＞
次に、本発明の第２の実施の形態について説明する。上述した本発明の第１の実施の形態では、観測データ空間内のあるピクセルｊに関して、そのピクセルｊからサンプリングした強度に基づいて格子不定の対数尤度比を計算するとともに、そのピクセルｊの距離に対応する複数の格子（観測点からの距離が等しい格子）に対して、対数尤度比に重みを掛けて加算する計算を行っている。この計算は、格子ｉに関する式（５）によって表されており、この計算はすべてのピクセルｊに関して行われる。すなわち、１つのピクセルｊに含まれる値が、複数の格子ｉの占有確率に影響を与えている。

【0092】

これを、逆に格子側から見れば、ある格子ｉは、複数のピクセルｊに含まれる値の加算を受けることになると言える。しかしながら、実際には複数のピクセルｊの観測データから、重複した同一の格子ｉに物体が重なって存在することはないため、ある格子ｉは、複数のピクセルｊに由来する情報のうち、一番確率の高いものだけを加算させるという考え方もできる。この考え方を反映したのが、下記の式（６）である。

【0093】

【数6】

【0094】

上記式（６）の右辺第１項に記載されているｍａｘは、配列の数値（ｋ＝１，・・・,ｊ）のうち最大のものだけをとるｍａｘ演算子である。右辺第１項は、ｍａｘ演算子は、格子不定の対数尤度比と重みＷ（Φ_ｉ，ｊ）との積に関して最大のものをとることを表しており、最終的に、格子不定の対数尤度比と重みＷ（Φ_ｉ，ｊ）との積が最大の値となるもののみが、格子ｉの占有確率を表す値として採用されることになる。

【0095】

また、図１２（ａ）及び（ｂ）に、上述した手法を用いて実際に占有格子地図を作成した例を示す。なお、図１２（ａ）には、本発明の第１の実施の形態における計算によって作成された占有格子地図が比較例として示されており、図１２（ｂ）には、本発明の第１の実施の形態における計算によって作成された占有格子地図が示されている。図１２（ａ）では、実際には物体が存在しない位置において、占有状態ｍであるとの結果が得られてしまっているが、一方、図１２（ｂ）では、占有状態ｍではないとの結果が得られている。

【0096】

なお、本発明の第１の実施の形態では占有格子地図計算部１０７における計算に式（５）を用いる一方、本発明の第２の実施の形態では式（６）を用いる点を除いて、本発明の第１及び第２の実施の形態は基本的に同一である。すなわち、本発明の第２の実施の形態における周辺環境推定装置１０の構成は、図２に示す構成と基本的に同一であり、占有格子地図計算部１０７で用いられる数式及び計算処理のみが異なっている。また、本発明の第２の実施の形態における周辺環境推定装置１０で行われる処理は、図６に示す処理と基本的に同一であり、ステップＳ２１３で用いられる数式及び計算処理のみが異なっている。

【0097】

以上説明したように、本発明の第２の実施の形態によれば、本発明の第１の実施の形態と同様に、確率変数に基づく新たなセンサーモデルが提案され、観測データから得られる強度の強弱に応じて占有確率（格子不定の対数尤度比）を計算するとともに、強度の強弱に応じて異なる分散の幅（分布関数）で複数の格子にわたって影響を与えることを考慮して、占有格子地図を作成することが可能となる。

【0098】

＜第３の実施の形態＞
次に、本発明の第３の実施の形態について説明する。上述した本発明の第１又は第２の実施の形態における手法では、静止している物体及び移動している物体に関して、その存在確率を計算し、図１０（ｂ）、図１１（ｂ）、図１２（ａ）及び（ｂ）のように占有格子地図を作成している。このとき、物体が占有格子地図上の占有確率は、式（５）や式（６）のように、加算することで更新できる。

【0099】

移動している物体と移動していない物体を選別して、例えば図１０（ｂ）に示すように、移動していない物体のみを選別して占有格子地図を作成することも可能であり、例えば図１１（ｂ）に示すように、移動している物体のみを選別して占有格子地図を作成することも可能である。

【0100】

ここで、車両に搭載されたセンサー５０を用いてセンシングした観測データから、道路のガードレールや壁を占有状態ｍである位置として検知し、車両が走行できる道路部分を占有状態ｍではない位置として検知した占有格子地図を作成したい場合を考える。この場合、例えば図１０（ａ）に示すように、移動していない物体のみを選別することで、移動している物体を排除した道路部分を作成することが可能である。

【0101】

しかしながら、移動している物体を排除するのではなく、例えば先行車などの移動している物体が検知された位置も、走行可能な領域（道路部分）であることを示す有力な証拠であると考えることもできる。つまり、ある観測の瞬間では、移動する先行車がその位置を占有している位置であっても、移動している先行車は走行しており、その位置に存在していることから、その位置は車両が走行できる道路部分であるという証拠と考えることができる。

【0102】

本発明の第３の実施の形態では、上述のように、移動している物体が存在している位置は、道路部分であるという意味で『非占有』の位置である証拠として考え、上記の式（５）式（６）を拡張して、移動している物体の位置については占有確率を下げるように計算する下記の式（７）及び式（８）を定義する。

【0103】

【数7】

【数8】

【0104】

式（７）は、上述した式（５）の式を拡張したものであり、式（８）は、上述した式（６）を拡張したものである。式（７）及び式（８）は両方共、それぞれ式（５）及び式（６）の観測項（左辺第１項）に対してｓｇｎ（ｊ）がさらに積算されている。ｓｇｎ（ｊ）は、ピクセルｊが、静止物の時は『１』、移動体の場合は『−１』となる符号関数である。なお、ピクセルｊから得られる相対速度（図４（ａ）の観測データ空間の横軸の値）の値を参照することで、ピクセルｊに含まれる観測データが静止物に対応しているか移動体に対応しているかを判別することができる。また、ここでは、式（７）及び式（８）の観測項（左辺第１項）において、単純に『１』又は『−１』の値をとるｓｇｎ（ｊ）を掛けているが、ｓｇｎ（ｊ）を物体の速度に依存する関数（例えば、物体の速度が小さいほど『１』に近づき、物体の速度が大きいほど『−１』に近づく関数）としてもよい。

【0105】

式（７）及び式（８）では、移動体の場合には、符号関数ｓｇｎ（ｊ）＝『−１』の積によって観測項（左辺第１項）が負の値となる。つまり、移動体の場合には、観測項（左辺第１項）が加算されるのではなく除算されるようにすることで、移動体が存在する場合は占有状態である確率を下げる（すなわち、占有状態ではない確率を上げる）ように更新計算を行う。

【0106】

以下、図１３を参照しながら、本発明の第３の実施の形態における周辺環境推定装置の構成について説明する。図１３には、本発明の第３の実施の形態における周辺環境推定装置の構成の一例が図示されている。図１３に示す周辺環境推定装置３０は、図２に示す構成に加えて、さらに静止物／移動体判別部３０１を有する。なお、観測データ空間処理部１０１、強度値尤度計算部１０３、分布関数計算部１０５は、図２と同一であり、ここでは説明を省略する。

【0107】

静止物／移動体判別部３０１は、観測データ空間に存在する特定のピクセルについて、ピクセルに関連する相対速度（ドップラー）を取得し、その相対速度から、特定のピクセルに存在する物体が静止物か移動体かを判別する機能を有する。静止物／移動体判別部３０１における判別結果は、占有格子地図計算部３０７に供給される。

【0108】

また、占有格子地図計算部３０７は、図２に示す占有格子地図計算部１０７と比較すると、静止物／移動体判別部３０１から供給される判別結果を受け取る機能、及び、この情報を用いて、式（７）及び式（８）で表されている計算を行う機能を更に有している。占有格子地図計算部３０７は、静止物／移動体判別部３０１から供給された判別結果に基づいて、例えば式（７）及び式（８）の観測項（左辺第１項）に含まれるｓｇｎ（ｊ）の値を決定して計算を行う。

【0109】

また、本発明の第３の実施の形態における周辺環境推定装置３０で行われる処理は、図６に示す処理と基本的に同一であるが、静止物／移動体判別部３０１において物体が静止物か移動体かを判別する処理が更に行われる点と、ステップＳ２１３で用いられる数式及び計算処理が異なっている。

【0110】

以下、図１４を参照しながら、本発明の第３の実施の形態に特徴的な処理について説明する。図１４は、本発明の第３の実施の形態に特徴的な処理の一例を示すフローチャートである。なお、図１４は、図６のステップＳ２１３内で行われる処理（占有格子地図計算部３０７による処理）に対応している。

【0111】

図１４において、占有格子地図計算部３０７は、静止物／移動体判別部３０１から判別結果を取得し（ステップＳ２１３１）、その判別結果から、観測された物体（処理対象のピクセルｊに含まれる物体）が静止物か移動体かを判断する（ステップＳ２１３３）。
観測された物体が静止物である場合には、式（７）又は式（８）の観測項（左辺第１項）に含まれるｓｇｎ（ｊ）の値を『１』に設定する（ステップＳ２１３５）。一方、観測された物体が移動体である場合には、式（７）又は式（８）の観測項（左辺第１項）に含まれるｓｇｎ（ｊ）の値を『−１』に設定する（ステップＳ２１３７）。そして、占有格子地図計算部３０７は、強度値尤度計算部１０３による計算結果（対数尤度比）、分布関数計算部１０５による計算結果（重み）を利用して、複数の格子に対して式（７）又は式（８）を用いて更新計算を行う（ステップＳ２１３９）。

【0112】

以上説明したように、本発明の第３の実施の形態によれば、本発明の第１及び第２の実施の形態と同様に、確率変数に基づく新たなセンサーモデルが提案され、観測データから得られる強度の強弱に応じて占有確率（格子不定の対数尤度比）を計算するとともに、強度の強弱に応じて異なる分散の幅（分布関数）で複数の格子にわたって影響を与えることを考慮して、占有格子地図を作成することが可能となる。さらに、本発明の第３の実施の形態によれば、移動体が存在している位置については非占有であると判断されるよう計算を行うことによって、より精度の高い占有格子地図を作成することが可能となる。

【0113】

＜第４の実施の形態＞
以下、本発明の第４の実施の形態について説明する。上述の本発明の第１〜第３の実施の形態では、ある瞬間の観測データのスナップショットから占有格子地図を作成する場合について説明したが、さらに、これらの占有格子地図を組み合わせた占有格子地図を作成することが可能である。なお、以下では、ある瞬間の観測データのスナップショットから作成された占有格子地図をローカル占有格子地図と呼び、複数のローカル占有格子地図を重ね合わせることを目的として座標系が変換された占有格子地図をグローバル占有格子地図と呼ぶ。

【0114】

グローバル占有格子地図は、円座標系ではなく直交座標系の格子で表される。各格子には、占有確率に関する値と、その格子の観測回数（何回その格子が観測されたかを示す情報）が格納される。

【0115】

以下、図１５（ａ）〜（ｃ）を参照しながら、ローカル占有格子地図に含まれている情報を、グローバル占有格子地図に反映させる手法について説明する。

【0116】

ローカル占有格子地図の各格子（円座標系）には、センサー５０から得られた観測データに基づく占有確率の値が格納されている。この情報をグローバル占有格子地図（直交座標系）に反映させる場合、ローカル及びグローバル占有格子地図は異なる格子座標系であることから、格子間の補間が必要となる。具体的には、グローバル占有格子地図の各格子の中心座標の位置が、円座標系のローカル占有格子地図のどの位置に該当するかを計算し、その該当位置が、周辺の円座標の各格子までｒ軸、Φ軸に沿ってどれだけずれているのかを求め、そのずれに応じた重みを掛けた周辺の各格子の値による補間の値を、グローバル占有格子地図の格子の値とする。

【0117】

以下、複数のローカル占有格子地図を組み合わせて、グローバル占有格子地図を作成する手法について説明する。例として、グローバル占有格子地図上において、ある瞬間のローカル占有格子地図（図１５（ｂ））と、観測者が移動した後の次の瞬間のローカル占有格子地図（図１５（ａ））を示す。また、図１５（ｃ）には、図１５（ａ）及び図１５（ｂ）のそれぞれのローカル占有格子地図を重ね合わせた状態を示す。

【0118】

この２つの時間フレームを組み合わせ、複数のローカル占有格子地図のそれぞれにおける占有確率に関する値（上述のように、直交座標系に合わせるように補間された値）を足すことで、グローバル占有格子地図の各格子の値を計算する。また、ローカル占有格子地図の値が足し合わせた回数、つまり、観測された回数は各格子によって異なるため、各格子で観測された回数を記録しておく。そうすることで、最終的に、グローバル占有格子地図から確率地図を表示させる際、各格子の占有確率に関する値の和を観測された回数で割ることで、占有確率に関する値の「平均値」を計算することができる。この平均値を最終的な占有確率の値として表示させることができる。

【0119】

上述の手法により作成したグローバル占有格子地図を図１６に示す。図１６に示すグローバル占有格子地図は、各格子には０から１の値をとる占有確率が含まれており、占有確率が占有状態を示す『１』に近いほど白く、非占有状態を示す『０』に近いほど黒く表示されている。

【0120】

図１６に示すグローバル占有格子地図は、上述した本発明の第３の実施の形態における手法に基づいて式（８）を用いて作成されたものであり、複数のローカル占有格子地図を組み合わせてグローバル占有格子地図が作成できることを示すとともに、本発明の第３の実施の形態における手法から有意な結果が得られることを示している。

【0121】

本発明の第３の実施の形態では、道路の壁などの静止物を占有領域とし、道路部分の走行できる領域を非占有領域と考えている。式（８）を用いて、道路の両側の防音壁などの静止物が検知した位置はｓｇｎ（ｊ）が正の符号『１』で加算された結果、図１６では占有領域として白く表示されている。一方、先行車が検知された位置はｓｇｎ（ｊ）が負の符号『−１』で加算された結果（すなわち、除算の項となる）、図１６では非占有領域として黒く表示されている。このように、本発明の第３の実施の形態における手法では、道路部分の走行できる領域が非占有領域として明確に判断され、図１６に示すように黒く表示されるので、占有領域及び非占有領域のコントラストがより明確になり、車両が走行できる領域、フリースペースがより推測しやすくなる。

【0122】

以上説明したように、本発明の第４の実施の形態によれば、複数のローカル占有格子地図を組み合わせて、グローバル占有格子地図を作成することが可能となる。さらに、例えば、移動体が存在している位置については非占有であるとした計算結果（本発明の第３の実施の形態における計算結果）から得られるローカル占有格子地図を重ね合わせて、グローバル占有格子地図を作成した場合に、占有領域及び非占有領域のコントラストがより明確に表示されるようになる。

【0123】

なお、本発明の各実施の形態において、装置構成の説明に用いられているブロック図は、本発明に関連した機能を表しているにすぎず、実際の実装では、ハードウェア、ソフトウェア、ファームウェア、又はそれらの任意の組み合わせによって実現されてもよい。ソフトウェアで実装される機能は、１つ又は複数の命令若しくはコードとして任意のコンピュータ可読媒体に記憶され、これらの命令又はコードは、ＣＰＵ（Central Processing Unit：中央処理ユニット）又はＥＣＵ（Engine Control Unit：エンジンコントロールユニット）などのハードウェアベースの処理ユニットによって実行可能である。また、本発明に関連した機能は、ＩＣ（Integrated Circuit：集積回路）やＩＣチップセットなどを含む様々なデバイスによって実現されてもよい。

【0124】

また、本発明の各実施の形態において、処理の説明に用いられているフローチャートは、ブロック図は、本発明を実現するための各処理の一例を表しているにすぎず、本発明の基本的な概念を実現するための処理を限定するものではない。また、フローチャート内の各処理の流れは、各処理が実行される順序を限定するものではなく、適宜順序の入れ替えや並行処理が行われてもよく、また、処理の省略や更なる処理の追加が行われてもよい。

【0125】

さらに、本明細書に開示されている内容は、たとえ明示されていない組み合わせであっても、任意の内容を適宜組み合わせることが可能である。さらに、本明細書に開示されている内容と既に知られている従来技術との組み合わせについても同様に、適宜組み合わせることが可能である。こうした組み合わせについても、本発明の開示範囲の一部に含まれる。

【産業上の利用可能性】

【0126】

本発明は、センシング機器（センサー）を用いて周辺をセンシングすることで得られた観測データから、周辺の環境を推定する技術に適用可能である。また本発明は、特に、車両に搭載されたセンシング機器を用いて得られた観測データから、車両周辺の環境を推定する技術に適用可能である。

【符号の説明】

【0127】

１０、３０ 周辺環境推定装置
５０ センサー
１０１ 観測データ空間処理部
１０３ 強度値尤度計算部
１０５ 分布関数計算部
１０７、３０７ 占有格子地図計算部
３０１ 静止物／移動体判別部

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】

【図9】

【図10】

【図11】

【図12】

【図13】

【図14】

【図15】

【図16】