(57)【要約】      （修正有）【課題】移動体の周囲の物体と前記移動体とが接触するまでの時間を算出する装置を提供する。【解決手段】移動体に設けられ、移動体の周囲の時間変化する画像を連続的に取得する画像センサと、連続する前記画像間での画素の輝度変化量を算出する輝度変化量算出手段と、前記輝度変化量算出手段によって算出された輝度変化量に基づき最適化輝度勾配を算出し、前記最適化輝度勾配に基づき前記物体と前記移動体とが接触するまでの時間を算出する接触時間算出手段と、を備える。【選択図】図１
2008059569.tif 000002
【特許請求の範囲】
【請求項１】
移動体の周囲の物体と前記移動体とが接触するまでの時間を算出する装置であって、
前記移動体に設けられ、前記移動体の周囲の時間変化する画像を連続的に取得する画像センサと、
連続する前記画像間での画素の輝度変化量を算出する輝度変化量算出手段と、
前記輝度変化量算出手段によって算出された輝度変化量に基づき最適化輝度勾配を算出し、前記最適化輝度勾配に基づき前記物体と前記移動体とが接触するまでの時間を算出する接触時間算出手段と、を備えた接触時間算出装置。
【請求項２】
請求項１に記載の接触時間算出装置であって、
前記接触時間算出手段は、前記時間変化する画像の少なくとも一部の画素における前記輝度変化量の総和を最小化する解を求めることにより前記最適化輝度勾配を算出する。
【請求項３】
請求項１または２に記載の接触時間算出装置であって、
前記接触時間算出手段は、前記物体の自由度及び前記移動体の自由度のうち少なくとも一方が拘束される場合といずれも拘束されない場合とに分けて、前記物体と前記移動体とが接触するまでの時間を算出する。
【請求項４】
請求項１から３のいずれか一項に記載の接触時間算出装置であって、
前記接触時間算出手段は、前記物体の自由度及び前記移動体の自由度のいずれもが拘束されない場合に、各計算においては前記物体の姿勢または前記移動体の並進経路のうちいずれか一方が固定されていると仮定して繰り返し計算を行うことにより、前記物体と前記移動体とが接触するまでの時間を算出する。
【請求項５】
請求項３または４に記載の接触時間算出装置であって、
前記接触時間算出手段は、前記物体が固定物体か移動物体かに基づいて、前記物体の自由度及び前記移動体の自由度のうち少なくとも一方が拘束される場合であるか、いずれも拘束されない場合であるかを決定する。
【請求項６】
請求項１から５のいずれか一項に記載の接触時間算出装置であって、
前記接触時間算出手段は、前記物体の平面の向き及び前記移動体が前記物体に向って進む並進経路に応じ、
前記平面が前記画像センサの光軸に対し垂直であり前記移動体が前記光軸に沿って前記物体に近づく第１ケースと、
前記平面が前記画像センサの光軸に対し垂直であり前記移動体が任意の方向から前記物体に近づく第２ケースと、
前記平面が任意の方向を向き前記移動体が前記光軸に沿って前記物体に近づく第３ケースと、
前記平面が任意の方向を向き前記移動体が任意の方向から前記物体に近づく第４ケースとに分けて、前記物体と前記移動体とが接触するまでの時間を算出する。
【請求項７】
移動体の周囲の物体と前記移動体とが接触するまでの時間を算出する方法であって、
前記移動体の周囲の時間変化する画像を前記移動体から連続的に取得する画像取得工程と、
連続する前記画像間での画素の輝度変化量を算出する輝度変化量算出工程と、
前記輝度変化量算出工程で算出された輝度変化量に基づき最適化輝度勾配を算出し、前記最適化輝度勾配に基づき前記物体と前記移動体とが接触するまでの時間を算出する接触時間算出工程と、を有する接触時間算出方法。
【請求項８】
請求項７に記載の接触時間算出方法であって、
前記接触時間算出工程は、前記時間変化する画像の少なくとも一部の画素における前記輝度変化量の総和を最小化する解を求めることにより前記最適化輝度勾配を算出する工程を含む。
【請求項９】
請求項７または８に記載の接触時間算出方法であって、
前記接触時間算出工程では、前記物体の自由度及び前記移動体の自由度のうち少なくとも一方が拘束される場合といずれも拘束されない場合とに分けて、前記物体と前記移動体とが接触するまでの時間を算出する。
【請求項１０】
請求項７から９のいずれか一項に記載の接触時間算出方法であって、
前記接触時間算出工程では、前記物体の自由度及び前記移動体の自由度のいずれもが拘束されない場合に、各計算においては前記物体の姿勢または前記移動体の並進経路のうちいずれか一方が固定されていると仮定して繰り返し計算を行うことにより、前記物体と前記移動体とが接触するまでの時間を算出する。
【請求項１１】
請求項９または１０に記載の接触時間算出方法であって、
前記接触時間算出工程では、前記物体が固定物体か移動物体かに基づいて、前記物体の自由度及び前記移動体の自由度のうち少なくとも一方が拘束される場合であるか、いずれも拘束されない場合であるかを決定する。
【請求項１２】
請求項７から１１のいずれか一項に記載の接触時間算出方法であって、
前記接触時間算出工程では、前記物体の平面の向き及び前記移動体が前記物体に向って進む並進経路に応じ、
前記平面が前記画像センサの光軸に対し垂直であり前記移動体が前記光軸に沿って前記物体に近づく第１ケースと、
前記平面が前記画像センサの光軸に対し垂直であり前記移動体が任意の方向から前記物体に近づく第２ケースと、
前記平面が任意の方向を向き前記移動体が前記光軸に沿って前記物体に近づく第３ケースと、
前記平面が任意の方向を向き前記移動体が任意の方向から前記物体に近づく第４ケースとに分けて、前記物体と前記移動体とが接触するまでの時間を算出する。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、２つの物体または面が接触するまでに残された時間を求めるための画像システムの使用に関するものである。画像センサは、手動操縦、遠隔操縦、あるいは自動操縦される陸上車両、水中車両、または航空機に、あるいは、自動組立装置または自動化工場における運動チェーンまたはロボットアームに装着することができる。
【背景技術】
【０００２】
物体までの、または、物体上の一点までの距離や間隔を推定するために、ソナー、レーザー、レーダを用いた能動的なシステムが知られている。このようなシステムは、物体の幅を微分したり、あるいはドップラー効果を利用することにより、速度を算出することができる。次いで、距離と速度の割合を利用して接触までの時間を推定することができる。このようなシステムは自ら電磁波を発信し、物体に反射した電磁波を検出する。発信波と受信波との間の時間遅延、発信波と受信波との間の位相ずれから距離が推定される。電磁波の受信センサが発信源と並列に設置されている場合には、三角測量法を用いて距離が推定される。 
【０００３】
このような「能動的な」システムは、複雑かつ高価であり、干渉の影響を受け、空間における分解能が限られ、また、自ら電磁波を発信し物体に反射させて受信することが必要であるため用途が限られる。これに対し、いわゆる「受動的な」システムは、電磁波を発信せず、物体からの自然発光、自然反射光、輻射電磁波により機能する。
【０００４】
受動型光学システムの一つに、２つの画像センサを用いて二眼ステレオ視の原理によって距離を推定するものがある。接触までの時間は、速度がわかっていれば計算することができる。二眼ステレオシステムは高価であり、慎重なキャリブレーションを必要とし、ステレオ照合処理に膨大な計算が要求される（非特許文献１）。
【０００５】
図１６は、自車両１０が、自車両１０の走行する道路に交差する道路を走行中の他車両２０に接近する状況を図示している。図１６に概略的に示すように、接触までの時間（ＴＴＣ）は、グレイ値コーナー部のような、ある「特徴点」の移動量から推定することができる。このような特徴点は、いわゆる「インテレスト・オペレータ（interest operators）」によって特定されるエッジ、点、または微小領域であってもよい。しかしながら、「特徴点追跡」法は、例えば、特徴点同士を誤って対応付けてしまった場合、あるいは、特徴点を誤って抽出したり抽出に失敗したりして特徴点同士の対応付けに失敗した場合のように、特徴または特徴点の抽出結果にエラーが含まれている場合、及び特徴点を用いて行われる計算にエラーが含まれている場合、ＴＴＣの推定における信頼性は低い（図１７に概略を示す）。また、特徴点の抽出には膨大な計算量が要求されるため、ＴＴＣの高速推定を実現することは困難である。その結果、「特徴点追跡」法は、ＴＴＣの推定において成功裏に用いられたことはない。  
【０００６】
接触までの時間は、「オプティカルフロー」に基づいた二段階で構成される手法で推定することができる。「オプティカルフロー」は、画像内の各点における画像の明らかな動きを特定するベクトル場である（非特許文献２）。物体に対する画像システムの動きは、オプティカルフローから推定することができる（非特許文献３）。次いで、衝突までの時間は、物体との距離が既知であれば、相対的な動きから推定することができる。
【０００７】
オプティカルフローの推定手法は、１９８０年代から知られている（非特許文献２）。画像輝度の、空間および時間による微分は、オプティカルフローの速度成分に制約を持たせる。空間および時間による微分に基づく手法は、「勾配法」としても知られている。しかしながら、オプティカルフローは、通常、二つの速度成分の二次偏微分方程式を数値解析的に解いて求められるため、オプティカルフロー手法は計算負荷が大きいという問題点がある。
【０００８】
粗いグリッドまたは散在点に基づく結果のみが要求されるのであれば、オプティカルフローは、種々の方法で求めることができる。そのうちの一つは、全体画像を分割して得られる部分画像ブロック内においてオプティカルフロー速度は一定または固定値であるという仮定に基づく最小自乗法である（非特許文献３のセクション4.3）。このような「固定フロー」手法は、上記フルオプティカルフロー手法に比べて計算量は少ないが、各ブロック内においてフローが一定であると仮定しているので精度が低いという問題があり、また、解像度と精度とを適度に両立させるための画像ブロックサイズを決定する手段が要求される。
【０００９】
ある種の光学式マウスは、画像輝度の空間的時間的微分を用いて「固定フロー」を推定することにより、面上におけるマウスの左右方向の動きを推定している。例えば、特許文献１、特許文献２、特許文献３を参照。
【００１０】
オプティカルフローは、画像フレームにおける各ブロックの、前画像フレームから最良の整合性を与える移動量を求めることによっても推定することができる。画像ブロック同士がどの程度整合するかを測る手段としては、例えば、相関、正規化相関、差分絶対値総和、あるいは、差分自乗総和がある（例えば、非特許文献５参照）。これらの手法は、前述の手法と同様に精度が低いという欠点を有しており、加えて、最良の整合性を与える移動量を探し出す必要があるため計算負荷が大きい。
【００１１】
オプティカルフローに基づく全ての手法は、オプティカルフローの推定精度を制約する、いわゆる「アパーチャ効果」の影響を受ける。「アパーチャ効果」は、局所画像の、オプティカルフローの等照線方向（すなわち、局所輝度勾配に直交する方向）の成分を拘束する能力の欠如によって生じる。
【００１２】
「衝突までの時間」を推定する他の手法として、生物の視覚の中立構造を模擬したものがある（非特許文献６参照）。これらの手法では、計算量が膨大であり、精度が低く、著しい遅延の後にようやく結果が得られる。回避行動をとるには遅過ぎるタイミングで警告を発するので、当然ながら、衝突警告システムとしては現実性がない。
【００１３】
剛体の動きを推定するのに、オプティカルフローによる推定を迂回することにより二段階手法を回避し、時間的かつ空間的な輝度の微分を直接利用する、いわゆる「直接法」と呼ばれる手法が提案されている。固定環境下におけるカメラの回転運動と共に、固定環境下におけるカメラの並進運動を推定する手法は既に記述されている。本明細書は、非特許文献７及び非特許文献４の内容を参照し、ここに援用する。しかしながら、これらの手法は、接触までの時間あるいは衝突までの時間を求めるまでには至っていない。
【００１４】
接触までの時間を求める手法を現実的なものとするためには、下級の汎用コンピュータあるいは廉価な特殊目的回路（例えば、ＤＰＳ，ＦＰＧＡ，ＰＬＤ，ＡＳＩＣ）によって必要な計算を実行する必要がある。さらに、計算結果は直ちに、例えば、解析される連続画像の１フレーム時間内に、利用可能となる必要がある。このことは、潜在期間を増長するような計算の継続の克服が必要であることを意味する。
【特許文献１】米国特許第５，７９３，３５７号明細書
【特許文献２】米国特許第６，０８４，５７４号明細書
【特許文献３】米国特許第６，１２４，５８７号明細書
【非特許文献１】ビー・ディー・ルーカス（B.D. Lucas）、ティー・カネイド（T. Kanade）共著、「イタラティブ・イメージ・レジストレーション・テクニーク・ウィズ・アプリケーション・トゥー・ステレオ・ビジョン（An Iterative Image Registration Technique with an Application to Stereo Vision）」、プロシーディング・オブ・イメージング・アンダースタンディング・ワークショップ（Proceedings of Imaging Understanding Workshop）、１９８１年、p.121-130
【非特許文献２】ビー・ケー・ピー・ホーン（B.K.P. Horn）、ビー・ジー・シャンク（B.G. Schunck）共著、「デターミニング・オプティカル・フロー（Determining Optical Flow）」、アーティフィシャル・インテリジェンス（Artificial Intelligence）、１９８１年８月、Vol. 16、No. 1-3、p.185-203
【非特許文献３】エー・アール・ブラス（A.R. Bruss）、ビー・ケー・ピー・ホーン（B.K.P. Horn）共著、「パッシブ・ナビゲーション」（Passive Navigation）、コンピュータ・ビジョン・グラフィクス・アンド・イメージ・プロセッシング（Computer Vision, Graphics, and Image Processing）、1983年1月、Vol. 21、No. 1、p. 3-20
【非特許文献４】ビー・ケー・ピー・ホーン（B.K.P. Horn）、イー・ジェー・ウェルダン・ジュニア（E.J. Weldon, Jr.）共著、「ディレクト・メソッド・フォー・リカバリング・モーション（Direct Methods for Recovering Motion）、インターナショナル・ジャーナル・オブ・コンピュータ・ビジョン（International Journal of Computer Vision）１９８８年６月、Vol. 2、No. 1、p. 51-76
【非特許文献５】ジェー・イー・タナー（J.E. Tanner）、シー・エー・ミード（C.A. Mead）共著、「コリレーティング・オプティカル・モーション・ディテクタ（A Correlating Optical Motion Detector）、プロシーディングス・コンフェレンス・オン・アドバンスト・リサーチ・イン・ヴィ・エル・エス・アイ（Proceedings, Conference On Advanced Research In VLSI）、１９８４年１月２３−２５日、p. 57-64
【非特許文献６】ジェー・エム・ガルブレイス（J.M. Galbraith）、ジー・ティー・ケニオン（G.T. Kenyon）、アール・ダブリュー・ジオルコスキー（R.W. Ziolkowski）共著、「タイム・トゥー・コリジョン・エスティメーション・フロム・モーション・ベイスト・オン・プライメート・ビジュアル・プロセッシング（Time-to-Collision Estimation from Motion Based on Primate Visual Processing）、アイ・イー・イー・イー・トランザクション・オン・パターン・アナリシス・アンド・マシン・インテリジェンス（IEEE Transactions on Pattern Analysis and Machine Intelligence）、２００５年８月、Vol. 27、No. 8、p. 1279-1291
【非特許文献７】ビー・ケー・ピー・ホーン（B.K.P. Horn）、 エス・ネガダリポア（S. Negahdaripour）共著、「ディレクト・パッシブ・ナビゲーション（Direct Passive Navigation）、アイ・イー・イー・イー・トランザクション・オン・パターン・アナリシス・アンド・マシン・インテリジェンス（IEEE Transactions on Pattern Analysis and Machine Intelligence、１９８７年１月、Vol. PAMI-9、No. 1、p. 168-176
【非特許文献８】アイ・エス・マククワーク（I.S. McQuirk）、ビー・ケー・ピー・ホーン（B.K.P. Horn）、 エイチ・エス・リー（H.-S. Lee）、ジェー・エル・ワイアット（J.L. Wyatt）共著、「エスティメーティング・ザ・フォーカス・オブ・エクスパンション・イン・アナログ・ヴィ・エル・エス・アイ（Estimating the Focus of Expansion in Analog VLSI）、インターナショナル・ジャーナル・オブ・コンピュータ・ビジョン（International Journal of Computer Vision）、１９９８年、Vol. 28、No. 3、p. 261-277
【発明の開示】
【発明が解決しようとする課題】
【００１５】
本発明の課題は、時間変化する画像あるいは画像システムからの連続画像に基づき、効率的な計算及び短い潜在期間を実現できる接触時間推定方法を提供すると共に、接触時間推定装置を提供することである。
【課題を解決するための手段】
【００１６】
上記課題を解決するために本発明は、移動体の周囲の物体と前記移動体とが接触するまでの時間を算出する装置であって、前記移動体に設けられ、前記移動体の周囲の時間変化する画像を連続的に取得する画像センサと、連続する前記画像間での画素の輝度変化量を算出する輝度変化量算出手段と、前記輝度変化量算出手段によって算出された輝度変化量に基づき最適化輝度勾配を算出し、前記最適化輝度勾配に基づき前記物体と前記移動体とが接触するまでの時間を算出する接触時間算出手段と、を備えた接触時間算出装置を提供する。
【００１７】
上記接触時間算出装置によれば、移動体が物体に接触するまでの時間を、１つの画像センサ（例えば、カメラ）を用いるだけで推定することができる。
【００１８】
さらに、従来の２レンズ式ステレオ視装置あるいはレーザー装置とは異なり、本発明の接触時間算出装置は、使用に際して正確な位置調整を必要としない。従って、本発明の接触時間算出装置を車両の衝突警告システムに適用した場合には、２レンズ式ステレオ視装置あるいはレーザー装置のように他のタイプのセンサを用いる場合と比較して、廉価なシステムを構成することができる。
【００１９】
さらに、本発明の接触時間算出装置では、ハードウェアとして必要になるのは１つの画像センサと計算機のみであるため、装置を小型にすることができる。従って、センサの設置に小さい空間しか利用できない小型車、あるいは、さらに小さい自動車にも、接触時間算出装置を容易に搭載することができる。
【００２０】
上記接触時間算出装置において、前記接触時間算出手段は、前記時間変化する画像の少なくとも一部の画素における前記輝度変化量の総和を最小化する解を求めることにより前記最適化輝度勾配を算出してもよい。
【００２１】
上記接触時間算出装置において、前記接触時間算出手段は、前記物体の自由度及び前記移動体の自由度のうち少なくとも一方が拘束される場合といずれも拘束されない場合とに分けて、前記物体と前記移動体とが接触するまでの時間を算出してもよい。
【００２２】
上記接触時間算出装置において、前記接触時間算出手段は、前記物体の自由度及び前記移動体の自由度のいずれもが拘束されない場合に、各計算においては前記物体の姿勢または前記移動体の並進経路のうちいずれか一方が固定されていると仮定して繰り返し計算を行うことにより、前記物体と前記移動体とが接触するまでの時間を算出してもよい。
【００２３】
上記接触時間算出装置において、前記接触時間算出手段は、前記物体が固定物体か移動物体かに基づいて、前記物体の自由度及び前記移動体の自由度のうち少なくとも一方が拘束される場合であるか、いずれも拘束されない場合であるかを決定してもよい。
【００２４】
上記接触時間算出装置において、前記接触時間算出手段は、前記物体の平面の向き及び前記移動体が前記物体に向って進む並進経路に応じ、前記平面が前記画像センサの光軸に対し垂直であり前記移動体が前記光軸に沿って前記物体に近づく第１ケースと、前記平面が前記画像センサの光軸に対し垂直であり前記移動体が任意の方向から前記物体に近づく第２ケースと、前記平面が任意の方向を向き前記移動体が前記光軸に沿って前記物体に近づく第３ケースと、前記平面が任意の方向を向き前記移動体が任意の方向から前記物体に近づく第４ケースとに分けて、前記物体と前記移動体とが接触するまでの時間を算出してもよい。
【００２５】
本発明はさらに、移動体の周囲の物体と前記移動体とが接触するまでの時間を算出する方法であって、前記移動体の周囲の時間変化する画像を前記移動体から連続的に取得する画像取得工程と、連続する前記画像間での画素の輝度変化量を算出する輝度変化量算出工程と、前記輝度変化量算出工程で算出された輝度変化量に基づき最適化輝度勾配を算出し、前記最適化輝度勾配に基づき前記物体と前記移動体とが接触するまでの時間を算出する接触時間算出工程と、を有する接触時間算出方法を提供する。
【００２６】
上記接触時間算出方法において、前記接触時間算出工程は、前記時間変化する画像の少なくとも一部の画素における前記輝度変化量の総和を最小化する解を求めることにより前記最適化輝度勾配を算出する工程を含んでいてもよい。
【００２７】
上記接触時間算出方法において、前記接触時間算出工程では、前記物体の自由度及び前記移動体の自由度のうち少なくとも一方が拘束される場合といずれも拘束されない場合とに分けて、前記物体と前記移動体とが接触するまでの時間を算出してもよい。
【００２８】
上記接触時間算出方法において、前記接触時間算出工程では、前記物体の自由度及び前記移動体の自由度のいずれもが拘束されない場合に、各計算においては前記物体の姿勢または前記移動体の並進経路のうちいずれか一方が固定されていると仮定して繰り返し計算を行うことにより、前記物体と前記移動体とが接触するまでの時間を算出してもよい。
【００２９】
上記接触時間算出方法において、前記接触時間算出工程では、前記物体が固定物体か移動物体かに基づいて、前記物体の自由度及び前記移動体の自由度のうち少なくとも一方が拘束される場合であるか、いずれも拘束されない場合であるかを決定してもよい。
【００３０】
上記接触時間算出方法において、前記接触時間算出工程では、前記物体の平面の向き及び前記移動体が前記物体に向って進む並進経路に応じ、前記平面が前記画像センサの光軸に対し垂直であり前記移動体が前記光軸に沿って前記物体に近づく第１ケースと、前記平面が前記画像センサの光軸に対し垂直であり前記移動体が任意の方向から前記物体に近づく第２ケースと、前記平面が任意の方向を向き前記移動体が前記光軸に沿って前記物体に近づく第３ケースと、前記平面が任意の方向を向き前記移動体が任意の方向から前記物体に近づく第４ケースとに分けて、前記物体と前記移動体とが接触するまでの時間を算出してもよい。
【００３１】
本発明の手法は、画像システムからの時間変化する画像を用い、光学システムの慎重なキャリブレーションを要求せず、計算効率が高い。接触時間を求める新規な手法は、中間結果としてのオプティカルフローの計算に依存しないという意味で「直接法」である。新規な手法は、輝度の空間的かつ時間的微分を直接利用する。その意味で、本手法は、「勾配に基づく」手法である。本手法は、誤差項の最小自乗法に基づくものであり、簡単なケースでは閉形式で実現される。計算には、注意領域における、あるいは画像全体における輝度の空間的時間的微分値と画像座標との積の特定和の算出が含まれる。
【００３２】
概して、本手法は、新しい各画像に対して以下の実工程を含む。
（１）累積器一式をゼロにリセットする。
（２）連続画像あるいは時間変化画像における画像を必要に応じて処理する。処理は、ローパスフィルタリング、ブロック平均化、及び／または、サブサンプリングを含む。
（３）各画像セル（あるいは、サブサンプリングされた画像の各セル）における画像輝度の空間的時間的微分値を求める。
（４）各画像セル（あるいは、サブサンプリングされた画像の各セル）における空間的時間的微分値と画像座標との積の特定和を算出する。算出されるべき特定の項は、相対運動及び表面形状に関する特定の仮定に依存する。
（５）累積器に保持されている総和に、これらの項を加算する。
（６）全ての画像セル（あるいはサブサンプリングされた画像のセル）の寄与度を累積総和に加算した後、相対運動のパラメータ及び表面形状のパラメータを算出する。相対運動のパラメータは、３個までの並進運動成分及び３個までの回転運動成分を含み得る。表面形状のパラメータは、表面距離、表面傾斜、あるいは、より一般的には、表面形状の級数展開係数、例えば、累乗級数の係数、フーリエ級数の係数を含む。
（７）運動及び形状のパラメータから接触までの時間を算出する。最も簡単なケースにおける接触までの時間は、表面までの距離と、表面に接近しつつある速度との比である。
【００３３】
上記ステップは、必ずしも順序どおり実行されなくてもよい。例えば、ステップ２，３，４，５は、全画像セルに対してステップ２を実行してから全画像セルに対してステップ３を実行し、さらにその次のステップに移るのではなく、１つの画像セルに対して各ステップを実行し、次の画像セルに移るといった交錯順序でもよい。
【発明を実施するための最良の形態】
【００３４】
以下、添付図面を参照して本発明の好ましい実施形態を説明する。
【００３５】
図１は、輝度勾配を用いて接触までの時間を推定する手法を示す概略図である。図１は、自車１０（移動体）が、自車１０の走行する道路に交差する道路を走行中の他車２０（物体）に接近するケースを示している。
【００３６】
図２は、座標系の定義を示す概略図であり、画像座標系における物体の位置ｐは座標（ｘ,ｙ）で、世界座標系における物体の位置Ｐは座標（Ｘ,Ｙ,Ｚ）で表現される。
【００３７】
図３は、オプティカルフローの概念を示す概略図である。図３に示すように、オプティカルフローは、物体の連続する画像間での移動ベクトルであり、物体の動きの測定または検出に用いられる。実際には、オプティカルフローは、第１画像における第１の点が、後続の第２画像においては前記第１の点と同じ輝度を有しかつ前記第１の点から僅かに移動した第２の点として現れるという仮定のもとに、２点間で引かれる移動ベクトルとして算出される。図４に示すように、画像内の点Ｐ（ｘ,ｙ）の時刻ｔにおける輝度値をＥ（ｘ,ｙ,ｔ）と表現する。時刻ｔから時刻ｔ＋Δｔの間に点Ｐが、輝度値を維持しながら（Δｘ,Δｙ）だけ移動したと仮定すると、下記の式（１）が成り立つ。
【００３８】
【数１】
2008059569.tif 000003
【００３９】
式（１）の右辺をテーラー展開し、２次以降の項を無視すると下記の式（２）を得る。
【００４０】
【数２】
2008059569.tif 000004
【００４１】
両辺をΔｔで割ると下記の式（３）を得る。
【００４２】
【数３】
2008059569.tif 000005
【００４３】
画像上の移動ベクトルは下記の式（４）で定義される。
【００４４】
【数４】
2008059569.tif 000006
【００４５】
輝度値の空間的勾配は下記の式（５）で定義される。
【００４６】
【数５】
2008059569.tif 000007
【００４７】
輝度値の時間的勾配は下記の式（６）で定義される。
【００４８】
【数６】
2008059569.tif 000008
【００４９】
式（４），（６），（６）を用いて式（３）を書き直し、下記の式（７）を得る。
【００５０】
【数７】
2008059569.tif 000009
【００５１】
式（７）は輝度勾配の拘束式である。式（７）を満たすパラメータｕ，ｖがオプティカルフロー、すなわち求めるべき移動ベクトルである。
【００５２】
実際のＥｘ,Ｅｙ,Ｅｔは、図５に示すように、時刻ｔおよび時刻ｔ＋Δｔの画像について、隣接する画素の輝度値を用いて以下のように算出する。ここで、ｉ，ｊは画像内の任意の座標である。
【００５３】
【数８】
2008059569.tif 000010
【００５４】
Ｅｘ,Ｅｙ,Ｅｔは、図６に示すように、キューブピクセルという概念でも表すことができる（非特許文献２）。
【００５５】
本発明では、接触までの時間を精度良く算出するために、物体（図中における平面１）と画像センサ（図中におけるカメラ２）の相対位置関係を、図７から１０に示す４ケースに分類する。
【００５６】
以下の説明において、「カメラ２の自由度」とは、カメラ２がその光軸２ａに沿って移動するか否かによって定義される。すなわち、カメラ２が光軸２ａに沿って移動する場合、カメラ２の自由度は拘束されているとみなされ、カメラ２の移動方向が必ずしも光軸２ａに一致しない場合、カメラ２の自由度は拘束されていないとみなされる。さらに、「平面１の自由度」とは、平面１が光軸２ａに対して垂直であるか否かによって定義される。すなわち、平面１が光軸２ａに対して垂直である場合、平面１の自由度は拘束されているとみなされ、平面１が光軸２ａに対して必ずしも垂直でない場合、平面１の自由度は拘束されていないとみなされる。
【００５７】
（ケース１）
図７に示すように、ケース１では、カメラ２の光軸２ａに平行な並進経路３に沿い、光軸２ａに垂直な平面１（すなわち、平面１の法線１ａが光軸２ａ及び並進経路３に対して平行）に向かうカメラ２の並進運動を考察する。この場合、平面１及びカメラ２の自由度が拘束され、すなわち、カメラ２は回転運動を伴わずに、光軸２ａに沿って直進するのみであるため、輝度勾配の拘束式は最も簡単な方程式で表現できる。さらに、対象とする物体を、カメラ２の光軸２ａに対して垂直な平面であると仮定すると、画像上の移動ベクトル（ｕ,ｖ）（ｕ＝ｄｘ／ｄｔ，ｖ＝ｄｙ／ｄｔ）は、単純に以下のように表現できる。
【００５８】
【数９】
2008059569.tif 000011
【００５９】
そして、接触までの時間ＴＴＣは、次のように表現される。
【００６０】
【数１０】
2008059569.tif 000012
【００６１】
ここで、Ｗは図２に示すＺ軸に沿ったカメラ２の速度、Ｚは平面１までの距離である。
【００６２】
式（７）より、以下の輝度勾配の拘束式が成り立つ。
【００６３】
【数１１】
2008059569.tif 000013
【００６４】
ここで、Ｃ（＝Ｗ／Ｚ）はＴＴＣの逆数である。
【００６５】
さらに、下記の各式Ｃ１−（４）が成り立つ。ここで、ｒは原点からの径方向距離である。
【００６６】
【数１２】
2008059569.tif 000014
【００６７】
未知のパラメータＣは、式Ｃ１−（３）に基づき、最小自乗法を用いることにより求めることができる。式Ｃ１−（３）の両辺を自乗して和を求めることにより、下記の総和項Ｃ１−（５）が求められる。
【００６８】
【数１３】
2008059569.tif 000015
【００６９】
総和項Ｃ１−（５）は、注意領域（ＲＯＩ）の画素、または全画像の画素から算出する。Ｃの最適解を得るために、総和項Ｃ１−（５）をＣで微分し、結果をゼロとする。
【００７０】
【数１４】
2008059569.tif 000016
【００７１】
式Ｃ１−（６）は、式Ｃ１−（７）に変換できる。
【００７２】
【数１５】
2008059569.tif 000017
【００７３】
ＴＴＣの逆数である未知のＣは、下記の式Ｃ１−（８）から求めることができる。
【００７４】
【数１６】
2008059569.tif 000018
【００７５】
図１１を参照し、ケース１におけるＴＴＣ算出手順を以下に示す。
（１）連続する２枚（または、３枚以上）の画像を用い、一画像及び前画像における画素の輝度と隣接する画素の輝度とから、注目領域の各画素における微分値Ｅｘ,Ｅｙ,Ｅｔを算出する（ステップＣ１−Ｓ１）。
（２）各画素において， Ｇ＝ｒＥｒ＝ｘＥｘ＋ｙＥｙを算出し、その結果から、ＧＥｔ及びＧ２を算出し（ステップＣ１−Ｓ２）、 これら２項をそれぞれ累積加算する。
（３）式Ｃ１−（８）を用いて未知のＣを推定する（ステップＣ１−Ｓ３）。
（４）Ｃの逆数としてＴＴＣを算出する（ステップＣ１−Ｓ４）。
【００７６】
空間的時間的微分値計算には、２×２×２キューブピクセル（画素）輝度値が優先的に用いられる（非特許文献２）。しかしながら、中心差分法など、連続する３以上の画像フレームを用いる手法（この手法は潜在時間が長くなる可能性がある）を含む他の微分値推定手法を用いることもできる。
【００７７】
ステップ（１）で算出される微分値は直ちにステップ（２）で用いられ、ステップ（２）で算出される項は直後のステップにおいて累積器に加算可能であり保存の必要がないため、中間結果のために大きな保存領域を確保する必要はない。
【００７８】
累積総和に対する寄与度は重み付けすることができる。これによって、例えば、より信頼のおける情報を、最終結果に、より大きく反映させることができる。これが行われる場合、両累積器の寄与度に対して同じ重み付け係数が付与されるべきである。例えば、画像上において顕著な特徴を有さない領域からの情報は、累積総和に対して有効に寄与しない。また、鏡面のような光輝部を有する領域は、鏡面反射に起因するエラーを低減するために圧縮してもよい。
【００７９】
注意領域の画素数に比例する算術ステップを有し、調査や繰り返し計算を必要とせず、複雑な方程式を解く必要もない計算手法が極めて直接的であることは明らかである。
【００８０】
さらに、エッジあるいはグレイ値コーナー部のような画像の「特徴部分」を検出する必要もない。このことにより、本手法は、小さいマイクロプロセッサ（アナログ・デバイス・ブラック・フィン・マイクロプロセッサなど）に、あるいは、ＦＰＧＡ、ＰＬＤ、ＡＳＩＣといった特殊目的回路にさえも搭載することができる。実際、計算を実行する回路を画像センサと一体化することもできる。同じ考察が、以下に述べる他のケースにも当てはまる。
【００８１】
（ケース２）
図８に示すように、ケース２では、任意に決定される並進経路３に沿い、光軸２ａに垂直な平面１（すなわち、平面１の法線１ａが光軸２ａに対して平行）に向かうカメラ２の並進運動を考察する。この場合、平面１の自由度が拘束され、カメラ２の自由度は拘束されていない。すなわち、カメラ２は、必ずしも光軸２ａに沿って移動するのではなく、また、必ずしも平面１の法線１ａに平行な直線に沿って移動するのではない。図２に示すＵ，Ｖ，Ｗは、それぞれ、世界座標系における速度のＸ，Ｙ，Ｚ成分を表す。これらの項に関して、下記の式Ｃ２−（１）（すなわち、輝度勾配の拘束式）が成り立つ。ここで、ｆは、カメラ２の焦点距離である。  
【００８２】
【数１７】
2008059569.tif 000019
【００８３】
上記考察におけるケース１は、単に、ケース２においてＵ＝Ｖ＝０とした特別な場合であることが明らかである。ここで、下記輝度勾配の拘束式が成り立ち、Ａ＝ｆ・Ｕ／Ｚ，Ｂ＝ｆ・Ｖ／Ｚ，Ｃ＝Ｗ／Ｚ（ケース１同様、ＣはＴＴＣの逆数）である。
【００８４】
【数１８】
2008059569.tif 000020
【００８５】
また、ケース１同様、下記の各式Ｃ２−（３）が成り立つ。
【００８６】
【数１９】
2008059569.tif 000021
【００８７】
未知のパラメータＡ，Ｂ，Ｃは、最小自乗法を用いることにより求めることができる。ここで、総和項は下記のとおりである。
【００８８】
【数２０】
2008059569.tif 000022
【００８９】
ケース１同様、総和項Ｃ２−（４）は、注意領域（ＲＯＩ）の画素、または全画像の画素から算出する。Ａ，Ｂ，Ｃの最適解を得るために、総和項Ｃ２−（４）をＡ，Ｂ，Ｃそれぞれについて微分し、３つの結果をゼロとする。
【００９０】
【数２１】
2008059569.tif 000023
【００９１】
３個の未知パラメータＡ，Ｂ，Ｃを含むこれら３つの線形方程式は、ガウスの消去法など公知の方法を用いて容易に解くことができる。
【００９２】
図１２を参照し、ケース２におけるＴＴＣ算出手順を以下に示す。
（１）連続する２枚（または、３枚以上）の画像を用い、一画像及び前画像における画素の輝度と隣接する画素の輝度とから、注目領域の各画素における微分値Ｅｘ,Ｅｙ,Ｅｔを算出する（ステップＣ２−Ｓ１）。
（２）各画素において， Ｇ＝ｒＥｒ＝ｘＥｘ＋ｙＥｙを算出する（ステップＣ２−Ｓ２）。
（３）式Ｃ２−（５）、Ｃ２−（６）、Ｃ２−（７）に基づき、未知のパラメータＡ，Ｂ，Ｃを推定する（ステップＣ２−Ｓ３）。
（４）Ｃの逆数としてＴＴＣを算出する（ステップＣ２−Ｓ４）。
【００９３】
係数の３×３対称行列において識別される６個の量、及び右辺のベクトルにおける３個の量は全て、輝度微分値Ｅｘ,Ｅｙ,Ｅｔ及び画像座標ｘ,ｙの種々の組合せ積和である。
【００９４】
接触までの時間は、Ｃ（＝Ｗ／Ｚ）の逆数である。接触時間を算出するために焦点距離ｆ（基準距離）を知る必要はない。前述したように、上記考察におけるケース１は、ケース２においてＡ＝Ｂ＝０とした特別な場合である。この特別な場合では、３つの方程式のうち第３の方程式のみを考えればよく、２つの項は不要となる。
【００９５】
接触時間に加えて、もし必要であれば、Ｕ／Ｗ＝（Ａ／Ｃ）／ｆ、Ｖ／Ｗ＝（Ｂ／Ｃ）／ｆであることに着目して、並進運動の方向を算出することができる。しかしながら、並進運動の大きさは、スケールファクタの多義性により、この手法では求めることができない。
【００９６】
２個の替わりに９個（６＋３）の累積器が必要であること、各画像セルにつき２項ではなく９項の計算が必要であることを除けば、基本的には上記ケース１と同様に８ステップの手順を用いて結果を算出することができる。
【００９７】
もし必要であれば、上記係数の３×３対称行列の固有値及び固有ベクトルを求めることにより、解法の安定性を評価することができる。特に小さい固有値が無いのであれば、結果は良好である。これに対し、少なくとも１つの固有値が小さい場合には、精度が低いと見込まれる。
【００９８】
（ケース３）
図９に示すように、ケース３では、カメラ２の光軸２ａに平行な並進経路３に沿い、任意方向を向く平面１（すなわち、平面１の法線１ａが必ずしも光軸２ａに対して平行ではない）に向かうカメラ２の並進運動を考察する。この場合、カメラ２の自由度が拘束され、平面１の自由度は拘束されていない。  
【００９９】
画像座標系で平面１のＸ，Ｙ方向における傾斜をｐ、ｑで定義した場合、平面１は下記の式Ｃ３−（１）で表現される。
【０１００】
【数２２】
2008059569.tif 000024
【０１０１】
式Ｃ３−（１）にＸ＝Ｚ・ｘ／ｆ、Ｙ＝Ｚ・ｙ／ｆを代入して下記の式Ｃ３−（２）を得る。
【０１０２】
【数２３】
2008059569.tif 000025
【０１０３】
輝度間の関係を用いると、これら移動ベクトルの拘束式は下記の式Ｃ３−（３）で表される。
【０１０４】
【数２４】
2008059569.tif 000026
【０１０５】
式Ｃ３−（３）を変換して下記の式Ｃ３−（４）を得る。ここで、Ｃ＝Ｗ／Ｚ０、Ｐ＝（ｐ／ｆ）・Ｃ、Ｑ＝（ｑ／ｆ）・Ｃ  である。
【０１０６】
【数２５】
2008059569.tif 000027
【０１０７】
上記考察におけるケース１は、単に、ケース３においてＰ＝Ｑ＝０とした特別な場合であることが明らかである。ケース１，２同様、未知のパラメータＰ，Ｑ，Ｃは、最小自乗法を用いることにより求めることができる。ここで、総和項は下記のとおりである。
【０１０８】
【数２６】
2008059569.tif 000028
【０１０９】
ケース１，２同様、総和項Ｃ３−（５）は、注意領域（ＲＯＩ）の画素、または全画像の画素から算出する。Ｐ，Ｑ，Ｃの最適解を得るために、総和項Ｃ３−（５）をＰ，Ｑ，Ｃそれぞれについて微分し、３つの結果をゼロとする。
【０１１０】
【数２７】
2008059569.tif 000029
【０１１１】
３個の未知パラメータＰ，Ｑ，Ｃを含むこれら３つの線形方程式は、ガウスの消去法など公知の方法を用いて容易に解くことができる。
【０１１２】
図１３を参照し、ケース３におけるＴＴＣ算出手順を以下に示す。
（１）連続する２枚（または、３枚以上）の画像を用い、一画像及び前画像における画素の輝度と隣接する画素の輝度とから、注目領域の各画素における微分値Ｅｘ,Ｅｙ,Ｅｔを算出する（ステップＣ３−Ｓ１）。
（２）各画素において， Ｇ＝ｒＥｒ＝ｘＥｘ＋ｙＥｙを算出する（ステップＣ３−Ｓ２）。
（３）式Ｃ３−（６）、Ｃ３−（７）、Ｃ３−（８）に基づき、未知のパラメータＰ，Ｑ，Ｃを推定する（ステップＣ３−Ｓ３）。
（４）Ｃの逆数としてＴＴＣを算出する（ステップＣ３−Ｓ４）。
【０１１３】
係数の３×３対称行列における９個の量、及び右辺のベクトルにおける３個の量は全て、輝度微分値Ｅｘ,Ｅｙ,Ｅｔ及び画像座標ｘ,ｙの積和である。
【０１１４】
接触までの時間は、Ｃ（＝Ｗ／Ｚ０）の逆数である。接触時間を算出するために焦点距離ｆ（基準距離）を知る必要はない。前述したように、上記考察におけるケース１は、本ケースにおいてＰ＝Ｑ＝０とした特別な場合である。この特別な場合では、３つの方程式のうち第３の方程式のみを考えればよく、２つの項は不要となる。
【０１１５】
もし必要であれば、ｐ＝ｆ（Ｐ／Ｃ）、ｑ＝ｆ（Ｑ／Ｃ）であることに着目して、平面の方向を算出することができる。
【０１１６】
２個の替わりに９個（６＋３）の累積器が必要であること、各画像セルにつき２項ではなく９項の計算が必要であることを除けば、基本的には上記ケース１と同様に８ステップの手順を用いて結果を算出することができる。
【０１１７】
もし必要であれば、上記係数の３×３対称行列の固有値及び固有ベクトルを求めることにより、解法の安定性を評価することができる。特に小さい固有値が無いのであれば、結果は良好である。これに対し、少なくとも１つの固有値が小さい場合には、精度が低いと見込まれる。
【０１１８】
基本的な解析手法は、詳細に説明した上記３ケースから明らかである。以下、より一般的なケースについてさらに説明する。
【０１１９】
（ケース４）
図１０に示すように、ケース４では、任意の並進経路３に沿い、任意方向を向く平面１（すなわち、平面１の法線１ａが必ずしも光軸２ａに対して平行ではない）に向かうカメラ２の並進運動を考察する。この場合、カメラ２及び平面１いずれの自由度も拘束されていない。
【０１２０】
ケース２における移動に関する方程式と、ケース３における平面に関する方程式とを組合せることにより、下記の式Ｃ４−（１）を得る。
【０１２１】
【数２８】
2008059569.tif 000030
【０１２２】
ここで、Ａ’＝ｆ・（Ｕ／Ｗ）、Ｂ’＝ｆ・（Ｖ／Ｗ）、Ｐ＝（ｐ／ｆ）・（Ｗ／Ｚ０）、Ｑ＝（ｑ／ｆ）・（Ｗ／Ｚ０）である。
【０１２３】
前述同様、Ｃ（＝Ｗ／Ｚ０）は、接触時間の逆数である。これは、下記のように書き換えることができる。
【０１２４】
【数２９】
2008059569.tif 000031
【０１２５】
ここで、Ａ＝ｆ・（Ｕ／Ｚ０）、Ｂ＝ｆ・（Ｖ／Ｚ０）、Ｐ’＝ｐ／ｆ、Ｑ’＝ｑ／ｆである。
【０１２６】
ケース１は、ケース４においてＡ＝Ｂ＝Ｕ＝Ｖ＝０とした単なる特別なケースであり、ケース２は、この一般ケースにおいてＡ＝Ｂ＝０とした特別なケースであり、ケース３は、この一般ケースにおいてＵ＝Ｖ＝０とした特別なケースである。
【０１２７】
５個の未知パラメータＡ’，Ｂ’，Ｃ，Ｐ，Ｑを求めるために最小自乗法を適用する。ここで、総和項は下記のとおりである。
【０１２８】
【数３０】
2008059569.tif 000032
【０１２９】
同様に、５個の未知パラメータＡ，Ｂ，Ｃ，Ｐ’，Ｑ’を求めるために最小自乗法を適用する。ここで、総和項は下記のとおりである。
【０１３０】
【数３１】
2008059569.tif 000033
【０１３１】
総和項Ｃ４−（３），Ｃ４−（４）は、注意領域（ＲＯＩ）の画素、または全画像の画素から算出する。Ａ，Ｂ，Ｐ，Ｑと、Ａ’，Ｂ’，Ｐ’，Ｑ’との間には次の関係がある。
Ａ＝Ａ’Ｃ，Ｂ＝Ｂ’Ｃ，Ｐ＝Ｐ’Ｃ，Ｑ＝Ｑ’Ｃ
【０１３２】
５個の未知パラメータの最適解を求めるために、５個のパラメータで２つの総和項のいずれかを式を微分し，その結果をゼロとする。これにより、５個の未知パラメータに対して５個の方程式を得ることができる。上記ケース同様、これら方程式における全ての係数は、輝度微分値及び画像座標の積和である。しかしながら、上記ケース１〜３とは異なり、これらの方程式は非線形なので、容易には解くことができない。以下、２つの解法について説明する。
【０１３３】
（解法１）
解法１では、接触時間（ＴＴＣ）を規定する未知のパラメータを求めるために繰返し演算処理を適用する。初めに、パラメータＰ’，Ｑ’を定数として扱うことで、総和項Ｃ４−（４）から導出される方程式は、残りの未知パラメータＡ，Ｂ，Ｃに関する線形方程式となり、Ａ，Ｂ，Ｃの最適解を求めるために、線形解法を適用することができる。次に、パラメータＡ’，Ｂ’を定数として扱うことで、総和項Ｃ４−（３）から導出される方程式は、残りの未知パラメータＰ，Ｑ，Ｃに関する線形方程式となり、Ｐ，Ｑ，Ｃの最適解を求めるために、線形解法を適用することができる。
【０１３４】
具体的な手順を説明する。初めに、パラメータＰ’，Ｑ’を定数として扱う。Ａ，Ｂ，Ｃの最適解を得るために、総和項Ｃ４−（４）をＡ，Ｂ，Ｃそれぞれについて微分し、３つの結果をゼロとする。
【０１３５】
【数３２】
2008059569.tif 000034
【０１３６】
ここで、Ｆ＝１−（ｐ／ｆ）ｘ−（ｑ／ｆ）ｙ  である。
【０１３７】
次に、パラメータＡ’，Ｂ’を定数として扱う。Ｐ，Ｑ，Ｃの最適解を得るために、総和項Ｃ４−（３）をＰ，Ｑ，Ｃそれぞれについて微分し、３つの結果をゼロとする。
【０１３８】
【数３３】
2008059569.tif 000035
【０１３９】
ここで、Ｄ＝Ａ’Ｅｘ＋Ｂ’Ｅｙ−ｒＥｒ  である。
【０１４０】
このようにして初期値が与えられると、Ｐ，Ｑを固定してＡ，Ｂ，Ｃについて解き、次いで、新たなＡ，Ｂの推定値を用い、Ａ，Ｂを固定してＰ，Ｑ，Ｃについて解くことができる。対をなすこれらのステップを数回繰り返すことで、厳密解に十分近い近似解を得ることができる。
【０１４１】
図１４を参照し、解法１におけるＴＴＣ算出手順を以下に示す。
（１）連続する２枚（または、３枚以上）の画像を用い、一画像及び前画像における画素の輝度と隣接する画素の輝度とから、注目領域の各画素における微分値（輝度勾配）Ｅｘ,Ｅｙ,Ｅｔを算出する（ステップＣ４−１−Ｓ１）。
（２）各画素において， Ｇ＝ｒＥｒ＝ｘＥｘ＋ｙＥｙを算出する（ステップＣ４−１−Ｓ２）。
（３）初期値を算出する（ステップＣ４−１−Ｓ３）。
（４）未知のパラメータＰ，Ｑ，Ｃを推定する（ステップＣ４−１−Ｓ４）。
（５）解が収束したか判定する（ステップＣ４−１−Ｓ５）。解が十分に収束した場合は処理を終了し、そうでない場合は、以下のステップを実行する。
（６）ｐ，ｑを算出する（ステップＣ４−１−Ｓ６）。
（７）Ｆ（＝１−（ｐ／ｆ）ｘ−（ｑ／ｆ）ｙ）を算出する（ステップＣ４−１−Ｓ７）。
（８）未知のパラメータＡ，Ｂ，Ｃを推定し（ステップＣ４−１−Ｓ８）、ステップ（３）に戻る。
【０１４２】
（解法２）
解法２を用いることによりＴＴＣの近似解を求めることができる。
【０１４３】
初めに、式Ｃ４−（１），Ｃ４−（２）から下記の項Ｃ４−（１１），Ｃ４−（１２）を導出する。
【０１４４】
【数３４】
2008059569.tif 000036
【０１４５】
添字Ｔは、行ベクトルを列ベクトルに（また、その逆に）入れ替える転置を意味する。これらの表記は、ベクトルと、３×３行列と、他のベクトルとの積を用いて以下のように書き換えることができる。
【０１４６】
【数３５】
2008059569.tif 000037
【０１４７】
ここで、Ｍは、ベクトル（Ａ  Ｂ  Ｃ）Ｔとベクトル（−Ｐ’  −Ｑ’  １）Ｔとのダイアディック積であり、ベクトル（Ａ’  Ｂ’  １）Ｔとベクトル（−Ｐ  −Ｑ  Ｃ）Ｔとのダイアディック積と等価である。すなわち、下記式が成り立つ。
【０１４８】
【数３６】
2008059569.tif 000038
【０１４９】
３×３行列Ｍの９個の未知パラメータが瞬間において独立の変数として扱われるならば、これら９個の未知パラメータを求めるために線形最小自乗法を用いることができる。上記項は（ｅＴｍ＋Ｅｔ）の形態で表記することができ、最小化されるべき総和は下記のようになる。
【０１５０】
【数３７】
2008059569.tif 000039
【０１５１】
ここで、ｍは、行列Ｍを「フラット化（flattening）」して得られるベクトル、すなわち、行列の９個の係数を並べたもの、ｅはベクトル（−Ｅｘ −Ｅｙ（ｒＥｒ））Ｔと、ベクトル（ｘ ｙ １）Ｔとのダイアディック積をフラット化したものである。  
【０１５２】
【数３８】
2008059569.tif 000040
【０１５３】
式Ｃ４−（１６），Ｃ４−（１７），Ｃ４−（１８）から下記の式Ｃ４−（１９），Ｃ４−（２０），Ｃ４−（２１）を得る。
【０１５４】
【数３９】
2008059569.tif 000041
【０１５５】
項Ｃ４−（１５）をベクトルｍについて微分し結果をゼロとして、下記の式Ｃ４−（２２），Ｃ４−（２３）を得る。
【０１５６】
【数４０】
2008059569.tif 000042
【０１５７】
これは、行列における９個の係数に関する９つの線形方程式であり、ガウスの消去法など公知の方法を用いて解くことができる。９×９対称行列及び方程式の右辺のベクトルの係数は、全て、輝度微分値及び画像座標の積和である。これらの項は、５４個（４５＋９）の累積器が存在し、各画像セルにつき５４項の計算が累積に必要であることを除けば、上記ケース同様、累積することができる。
【０１５８】
最後に、５個の未知パラメータＡ，Ｂ，Ｃ，Ｐ’，Ｑ’、あるいは等価なＡ’，Ｂ’，Ｃ，Ｐ，Ｑは、行列Ｍを特異値分解（singular value decomposition）することで求めることができる。
【０１５９】
別の方法として、この特殊なケースでは、行列における各列が（Ａ，Ｂ，Ｃ）に比例し、行列における各行が（−Ｐ，−Ｑ，Ｃ）に比例することに着目して未知パラメータを求めることもできる。実際、４個の未知パラメータＡ，Ｂ，Ｐ，Ｑは、フラット化されたベクトルｍを再編成することにより、行列Ｍの最終行及び最終列から読み取ることができる。
【０１６０】
しかしながら、この計算の結果として得られる３×３行列は、３個の対角項に加えられた任意の係数λを含む可能性がある。なぜなら、
−ｘＥｘ−ｙＥｙ＋（ｒＥｒ）＝０ 
であるからである。正しい対角項は、特異値分解を用いて求めることができる。あるいは、特にこのケースでは、より簡単に、非対角項の関係から求めることができる。
【０１６１】
【数４１】
2008059569.tif 000043
【０１６２】
式Ｃ４−（２６）から、ＴＴＣの逆数であるＣが求められる。接触までの時間を求めることが目的であるならば、実際には、行列における下段右の係数（式Ｃ４−（２６））のみが必要である。
【０１６３】
解法２は、上記解法１によって解いた最小自乗問題を全く同様に解くわけではないが、便宜的な近似形態の結果を与える。測定ノイズ（誤差）が無ければ結果は同じになる。画像輝度の測定にノイズがある場合、すなわち、Ｅｘ，Ｅｙ，Ｅｔにノイズがある場合、解法２による結果は真値の近似解となる。これは、行列Ｍの９個の要素が独立ではなく、５個の未知パラメータに依存するからである。この近似解法は、直接用いても場合によっては十分であり、あるいは、最小自乗問題を繰り返し演算を用いて解く際（すなわち解法１における）の非線形数値解法の初期値推定には十分である。
【０１６４】
前述同様、接触までの時間はパラメータＣの逆数である。この計算は、基準距離ｆが未知であっても実行することができる。
【０１６５】
平面の向きが
ｐ＝ｆ（Ｐ／Ｃ），ｑ＝ｆ（Ｑ／Ｃ）
で規定されるように、必要であれば、並進運動の方向を下記式で求めることができる。
【０１６６】
【数４２】
2008059569.tif 000044
【０１６７】
図１５を参照し、解法２におけるＴＴＣ算出手順を以下に示す。
（１）連続する２枚（または、３枚以上）の画像を用い、一画像及び前画像における画素の輝度と隣接する画素の輝度とから、注目領域の各画素における微分値（輝度勾配）Ｅｘ,Ｅｙ,Ｅｔを算出する（ステップＣ４−２−Ｓ１）。
（２）ベクトルｅを定義する（ステップＣ４−２−Ｓ２）。
（３）最小自乗法を適用する（ステップＣ４−２−Ｓ３）。
（４）行列Ｍの対角要素を算出する（ステップＣ４−２−Ｓ４）。
（５）ベクトルｍの各要素を算出する（ステップＣ４−２−Ｓ５）。
（６）ＴＴＣを推定する（ステップＣ４−２−Ｓ６）。
【０１６８】
以上説明した実施形態は例示に過ぎず、輝度微分値を用いて接触までの時間を求める手法の詳細については種々の代替手法が存在することは当業者には明らかである。さらに、並進運動と回転運動との他の組合せ、及び、他の表面形状についても、本手法を用いて取り扱うことが可能である。
【０１６９】
もし複数の移動体が存在するならば、個々の移動体の画像に対応する注意領域毎に画像を分割することができる。
【０１７０】
接触までの時間は、センサの突出部の中心が、視認されている物体の表面と交差するまでの時間と考えられることが多い。形状既知の突出物体にセンサが装着されているのであれば、その形状が、視認されている物体と最初に交差するまでの時間を求めることができる。
【０１７１】
接触までの時間を算出する際には、通常、運動パラメータは変化しないものと仮定する。センサまたはセンサが装着された物体の運動力学が既知であり、その運動に対する外乱が既知である場合には、運動力学モデルを構築し、例えば、カルマンフィルタ（Kalman filter）を用いて更新を行ってもよい。
【０１７２】
画像センサは、可視光用のものである必要はなく、例えば、赤外線センサ、紫外線センサ、あるいは、その他の電磁波用センサであってもよい。
【０１７３】
画像システムは、識別可能なスナップショットあるいは「フレーム」を撮像するものでなくてもよい。例えば、光ダイオードのように、連続的に動作する検出器であってもよい。時間微分値は、このような連続信号から、公知の電気回路を用いて求めることができる。
【０１７４】
画像輝度自体は直接的に計算に用いられるわけではないので、画像システムは、必ずしも実際の画像を出力する必要はなく、接触時間の算出に必要な画像輝度微分値のみを生成してもよい。
【０１７５】
画像センサは、必ずしも矩形グリッド形状の画像セルを有している必要はなく、例えば、六角形グリッド形状、三角形グリッド形状、あるいは、空間的に連続形状であってもよい。このような形状は、輻射の絶対量を測定するのではなく、輻射の変化にのみ反応する画像センサの場合には重要である可能性がある。
【０１７６】
輝度微分値及び微分値と画像座標との積の算出は、必ずしもデジタル方式で行う必要はなく、アナログ方式でも、ハイブリッド方式であってもよい。例えば、乗算は、４象限アナログ乗算器を用いて、あるいは、デジタル−アナログ（Ｄ／Ａ）変換乗算器を用いて実行可能である（非特許文献８）。
【０１７７】
微分値と座標の積和の加算は、例えば、抵抗を流れる電流を加算したり、キャパシタに蓄積される電荷を加算したりすることにより、アナログ方式でも実行することができる（非特許文献８）。演算にアナログ方式、ハイブリッド方式を採用した場合、構成要素の数、または製造コストの点で有利である。
【図面の簡単な説明】
【０１７８】
【図１】輝度勾配を用いて接触までの時間を推定する手法を示す概略図である。
【図２】座標系の定義を示す概略図である。
【図３】オプティカルフローの概念を示す概略図である。
【図４】輝度取得の概念を示す概略図である。
【図５】輝度微分値の算出方法を示す概略図である。
【図６】キューブピクセル（画素）による輝度の微分表現を示す概略図である。
【図７】平面及びカメラいずれの自由度も拘束されたケース１における平面とカメラとの関係を示す概略図である。
【図８】平面の自由度が拘束されたケース２における平面とカメラとの関係を示す概略図である。
【図９】カメラの自由度が拘束されたケース３における平面とカメラとの関係を示す概略図である。
【図１０】平面及びカメラいずれの自由度も拘束されていないケース４における平面とカメラとの関係を示す概略図である。
【図１１】ケース１における接触時間推定フローを示す概略図である。
【図１２】ケース２における接触時間推定フローを示す概略図である。
【図１３】ケース３における接触時間推定フローを示す概略図である。
【図１４】ケース４における解法１の接触時間推定フローを示す概略図である。
【図１５】ケース４における解法２の接触時間推定フローを示す概略図である。
【図１６】特徴点追跡法による接触時間推定を示す概略図である。
【図１７】特徴点追跡法における問題点を示す概略図である。
【符号の説明】
【０１７９】
１…平面、１ａ…法線、２…カメラ、２ａ…光軸、３…並進経路、１０…自車（移動体）、他車（物体）
【図１】
2008059569.tif 000045
【図２】
2008059569.tif 000046
【図３】
2008059569.tif 000047
【図４】
2008059569.tif 000048
【図５】
2008059569.tif 000049
【図６】
2008059569.tif 000050
【図７】
2008059569.tif 000051
【図８】
2008059569.tif 000052
【図９】
2008059569.tif 000053
【図１０】
2008059569.tif 000054
【図１１】
2008059569.tif 000055
【図１２】
2008059569.tif 000056
【図１３】
2008059569.tif 000057
【図１４】
2008059569.tif 000058
【図１５】
2008059569.tif 000059
【図１６】
2008059569.tif 000060
【図１７】
2008059569.tif 000061
【外国語明細書】
2008059569.tif 000062
2008059569.tif 000063
2008059569.tif 000064
2008059569.tif 000065
2008059569.tif 000066
2008059569.tif 000067
2008059569.tif 000068
2008059569.tif 000069
2008059569.tif 000070
2008059569.tif 000071
2008059569.tif 000072
2008059569.tif 000073
2008059569.tif 000074
2008059569.tif 000075
2008059569.tif 000076
2008059569.tif 000077
2008059569.tif 000078
2008059569.tif 000079
2008059569.tif 000080
2008059569.tif 000081
2008059569.tif 000082
2008059569.tif 000083
2008059569.tif 000084
2008059569.tif 000085
2008059569.tif 000086
2008059569.tif 000087
2008059569.tif 000088
2008059569.tif 000089
2008059569.tif 000090
2008059569.tif 000091
2008059569.tif 000092
2008059569.tif 000093
2008059569.tif 000094
2008059569.tif 000095
2008059569.tif 000096
2008059569.tif 000097
2008059569.tif 000098
2008059569.tif 000099
2008059569.tif 000100
2008059569.tif 000101
2008059569.tif 000102
2008059569.tif 000103