(57)【要約】【課題】  ビデオカメラを使って手動操作により撮影した画像から外界の物体の空間情報を復元するときにそれぞれ焦点距離が変化するため、従来の因子分解法を利用することができず、ピント調整のために変化する焦点距離を同時に復元する必要があるという問題を解決する。【解決手段】  因子分解処理部４３で各画像の特徴点の時間的変化量を示す画像座標値から行列データを作成し、この行列データを特異値分解し雑音除去を行って運動情報を表す行列データと３次元情報を表す行列データとを得て、運動情報の成分において運動を規定するために設定した条件を満足する変換行列を求め、運動情報となる行列データに変換行列を作用させてカメラ視点運動を復元し、３次元情報を表す行列データに変換行列の逆行列を作用させて物体形状を構成する３次元情報を復元し、光軸並進運動復元・焦点距離推定部４４でカメラ視点の光軸方向の並進運動の再現、焦点距離の推定を行う。【選択図】  図４
2006106930.tif 000002
【特許請求の範囲】
【請求項１】
カメラ手段によって取得した対象物の時系列画像に設定した画像座標系において、各画像の特徴点の時間的変化量を示す画像座標値を観測する特徴点観測手段と、
前記画像座標値から行列データを作成し、当該行列データから前記カメラ手段のカメラ視点運動と前記対象物の外観形状を構成する３次元形状とを復元する因子分解処理手段と、
前記復元したカメラ視点運動と３次元形状とから各特徴点に関する画像座標系でのカメラ視点運動の再投影誤差である再投影座標値を求め、当該再投影座標値と前記観測した画像座標値とから、前記カメラ視点の光軸方向の並進運動を再現し、かつ、前記カメラ手段の焦点距離を推定する光軸並進運動復元・焦点距離推定手段と、
前記再投影誤差と、前記観測した画像座標値に前記カメラ視点運動および前記３次元形状から構成される係数並びにカメラ内部パラメータから構成される係数を掛けた座標値との間の誤差が規定値以下に収束しているかを判断し、収束していない場合には、前記カメラ視点運動および前記３次元形状から構成される係数並びにカメラ内部パラメータから構成される係数を掛けた画像座標値を新たな行列要素として行列データを作成して前記因子分解処理手段からの処理を反復させる反復処理手段と、を備えることを特徴とするパラメータ推定装置。
【請求項２】
前記光軸並進運動復元・焦点距離推定手段では、前記復元したカメラ視点運動と3次元形状とから各特徴点に関する画像座標系でのカメラ視点運動の再投影誤差である再投影座標値を求め、当該再投影座標値と観測した画像座標値とから、前記カメラ視点の光軸方向の並進運動を再現し、かつ、スケールファクタ比と焦点距離とを推定し、
前記反復処理手段では、前記再投影誤差と前記観測した画像座標値に前記カメラ視点運動および前記3次元形状から構成される係数並びにカメラ内部パラメータから構成される係数を掛けた座標値との間の誤差が規定値以下に収束しているかを判断し、収束していない場合には、前記カメラ視点運動および前記3次元形状から構成される係数並びにカメラ内部パラメータから構成される係数を掛けた画像座標値を新たな行列要素として行列データを作成して前記因子分解処理手段からの処理を反復させ、誤差が前記規定値以下に収束するまで当該反復処理を行い、
前記反復処理が終了した時に、前記推定した焦点距離およびスケールファクタ比を用いてスケールファクタ比を算出するスケールファクタ比算出手段を備える、ことを特徴とする請求項１に記載のパラメータ推定装置。
【請求項３】
前記因子分解処理手段では、
前記作成した行列データを特異値分解し雑音除去を行って運動情報を表す行列データと３次元情報を表す行列データとを得て、
前記運動情報の成分において、運動を規定するために設定した条件を満足する変換行列を求め、運動情報となる行列データにこの変換行列を作用させてカメラ視点運動を復元し、
前記３次元情報を表す行列データにこの変換行列の逆行列を作用させて物体形状を構成する３次元情報を復元することを特徴とする請求項１または２に記載のパラメータ推定装置。
【請求項４】
特徴点観測手段がカメラ手段によって取得した対象物の時系列画像に設定した画像座標系において、各画像の特徴点の時間的変化量を示す画像座標値を観測する特徴点観測ステップと、
因子分解処理手段が前記画像座標値から行列データを作成し、当該行列データから前記カメラ手段のカメラ視点運動と前記対象物の外観形状を構成する3次元形状とを復元する因子分解処理ステップと、
光軸並進運動復元・焦点距離推定手段が前記復元したカメラ視点運動と3次元形状とから各特徴点に関する画像座標系でのカメラ視点運動の再投影誤差である再投影座標値を求め、当該再投影座標値と前記観測した画像座標値とから、前記カメラ視点の光軸方向の並進運動を再現し、かつ、前記カメラ手段の焦点距離を推定する光軸並進運動復元・焦点距離推定ステップと、
反復処理手段が前記再投影誤差と、前記観測した画像座標値に前記カメラ視点運動および前記3次元形状から構成される係数並びにカメラ内部パラメータから構成される係数を掛けた座標値との間の誤差が規定値以下に収束しているかを判断し、収束していない場合には、前記カメラ視点運動および前記3次元形状から構成される係数並びにカメラ内部パラメータから構成される係数を掛けた画像座標値を新たな行列要素として行列データを作成して前記因子分解処理手段からの処理を反復させる反復処理ステップと、を備えることを特徴とするパラメータ推定方法。
【請求項５】
前記光軸並進運動復元・焦点距離推定ステップでは、前記復元したカメラ視点運動と3次元形状とから各特徴点に関する画像座標系でのカメラ視点運動の再投影誤差である再投影座標値を求め、当該再投影座標値と観測した画像座標値とから、前記カメラ視点の光軸方向の並進運動を再現し、かつ、スケールファクタ比と焦点距離とを推定し、
前記反復処理ステップでは、前記再投影誤差と前記観測した画像座標値に前記カメラ視点運動および前記3次元形状から構成される係数並びにカメラ内部パラメータから構成される係数を掛けた座標値との間の誤差が規定値以下に収束しているかを判断し、収束していない場合には、前記カメラ視点運動および前記3次元形状から構成される係数並びにカメラ内部パラメータから構成される係数を掛けた画像座標値を新たな行列要素として行列データを作成して前記因子分解処理手段からの処理を反復させ、誤差が前記規定値以下に収束するまで当該反復処理を行い、
前記反復処理が終了した時に、スケールファクタ比算出手段が前記推定した焦点距離およびスケールファクタ比を用いてスケールファクタ比を算出するスケールファクタ比算出ステップを備える、ことを特徴とする請求項４に記載のパラメータ推定方法。
【請求項６】
前記因子分解処理ステップでは、
前記作成した行列データを特異値分解し雑音除去を行って運動情報を表す行列データと３次元情報を表す行列データとを得て、
前記運動情報の成分において、運動を規定するために設定した条件を満足する変換行列を求め、運動情報となる行列データにこの変換行列を作用させてカメラ視点運動を復元し、
前記３次元情報を表す行列データにこの変換行列の逆行列を作用させて物体形状を構成する３次元情報を復元することを特徴とする請求項４または５に記載のパラメータ推定装置。
【請求項７】
上記の請求項１〜６のいずれか１項に記載のパラメータ推定装置またはパラメータ推定方法を、コンピュータプログラムで記載してそれを実行可能にしたことを特徴とするパラメータ推定プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、画像入力装置（例えば、カメラ等）を使って取得した車載画像、室内画像、歩行撮影画像、船上からの海上画像、空撮画像などの時系列画像全般に利用可能で、カメラで取得した時系列画像から、カメラ視点の並進運動と光軸周りの回転（ヨー回転）、及び、時系列画像に映っている外界の３次元形状、すなわち、被写体（物体）の外観形状を構成する３次元情報を復元し、同時に、カメラ内部パラメータを推定することに関する。
【背景技術】
【０００２】
コンピュータビジョン分野では、時系列画像データから、対象物の形状を計測または獲得する手法には、ステレオ計測やエピポーラ面解析を用いた３次元解析手法がある。この手法によれば、物体が撮影されている複数の時系列画像から、空間形状または空間構造に関する３次元位置情報、及び、カメラ視点に関する運動を復元することができる。しかし、カメラを動かしながら撮影した時系列画像にはランダム性の雑音が混入するため、カメラ運動と物体形状を同時に、かつ、高精度に復元することが困難である。
【０００３】
これに対して、従来の因子分解法はランダム性の雑音が混入するようなビデオ画像などの時系列画像において、カメラ運動と物体形状をロバスト、かつ、高精度に復元することができる(例えば、特許文献１参照。)。
【０００４】
しかし、この手法では、カメラの内部パラメータ（焦点距離など）を既知としており、事前にカメラキャリブレーションによりカメラ内部パラメータを校正しておく必要がある。
【０００５】
一方、ディジタル・スチルカメラやハンディカメラなどにはオートフォーカス機能が搭載されており、被写体との距離に合わせてピント調節を行い、ピントの合った画像の取得が可能である。これにより、オートフォーカスで撮影された画像は各画像において焦点距離が微妙に調整される。従来法の因子分解法ではカメラの内部パラメータが固定（画像において焦点距離が変わらず一定）であることを前提としてカメラ運動と３次元形状を復元するため、オートフォーカスで撮影された画像には応用することができなかった。
【０００６】
これに対して、簡易な撮影ターゲット（紙に印刷した白黒の市松模様パターン）をハンドフリーで撮影した画像から、カメラの内部パラメータを推定する技術も開示されている(例えば、非特許文献１参照。)。しかし、時間的に変化する内部パラメータに対しては、この技術は適用できなかった。
【特許文献１】特開２００３−２７１９２５号公報
【非特許文献１】Ｚ．Ｚｈａｎｇ，「Ａ  Ｆｌｅｘｉｂｌｅ  Ｎｅｗ  Ｔｅｃｈｎｉｑｕｅ  ｆｏｒ  Ｃａｍｅｒａ  Ｃａｌｉｂｒａｔｉｏｎ」，ＩＥＥＥ  Ｔｒａｎｓ．Ｐａｔｔｅｒｎ  Ａｎａｌ．＆  Ｍａｃｈ．Ｉｎｔｅｌｌ，Ｖｏｌ．２２，Ｎｏ．１１，ｐｐ．１３３０−１３３４，２０００．
【発明の開示】
【発明が解決しようとする課題】
【０００７】
Ｓｈａｐｅ  Ｆｒｏｍ  Ｍｏｔｉｏｎの原理（カメラが運動することで３次元形状を復元する考え方）を応用した計測方法により、ディジタル・スチルカメラやハンディタイプのビデオカメラを使って撮影した時系列画像から外界の物体の空間情報を獲得、復元するとき、各画像でのピント合わせは必要であり、市販のカメラには被写体との距離に応じてフォーカスを調整するというオートフォーカス機能が装備されており、この機能により画像ごとにピントが調整された画像が得られる。しかし、オートフォーカス機能により各画像を取得した際の焦点距離が変化するため、従来の因子分解法を利用することができず、ピント調整のために変化する焦点距離を同時に復元しなければならないという問題がある。
【０００８】
本発明は、かかる事情に鑑みてなされたもので、その目的は上記課題を解決したカメラパラメータ推定技術を提供することにある。
【課題を解決するための手段】
【０００９】
本発明では、焦点距離を推定するための投影モデルを
【００１０】
【数１】
2006106930.tif 000003
【００１１】
と仮定し、この投影モデルに基づいてカメラ運動、３次元形状、焦点距離を復元する。式１では、第ｉフレームの画像において観測される第ｊ番目の特徴点の画像座標値を（ｘij，ｙij）、そのときの焦点距離をｆi、ＸＹ並進運動を（Ｔｘi，Ｔｙi，Ｔｚi）、光軸周りの回転をθiとする。また、画像座標値（ｘij，ｙij）は主点（ｐｒｉｎｃｉｐｌｅ  ｐｏｉｎｔ）を原点としており、主点は既知とし、画像の横縦のスケールファクタ比は１：１とする。
【００１２】
上記を全フレーム（ｉ＝１，２，・・・，Ｆ）、全特徴点（ｊ＝１，２，・・・，Ｐ）に対して行列表記すると、
【００１３】
【数２】
2006106930.tif 000004
【００１４】
となる。ここで、ｆi＝１（ｉ＝１，２，・・・，Ｆ）とした場合は、（式２）右辺はカメラ運動を平面運動と仮定し、因子分解法を利用してカメラ運動と３次元形状を復元することができる。ところが、カメラ運動により光軸（Ｚ軸）方向の並進運動が変化し、これに伴い、オートフォーカス機能により焦点距離ｆiが変化してフレームごとに異なるため、式２左辺のＴｚi（光軸方向の並進運動）と焦点距離ｆiを同時に推定しなければならない。そこで、式２左辺に示す係数（１−Ｔｚi／Ｚｊ）／ｆiを掛けた画像座標値（（１−Ｔｚi／Ｚｊ）ｘij／ｆi，（１−Ｔｚi／Ｚｊ）ｙji／ｆi）から平面運動と３次元形状を復元する。
【００１５】
因子分解法により復元した平面運動と３次元形状から再投影座標値（Δｕij，Δｖij）は、
【００１６】
【数３】
2006106930.tif 000005
【００１７】
と表される。式３全特徴点、全フレームに展開すると、
【００１８】
【数４】
2006106930.tif 000006
【００１９】
となり、光軸並進運動Ｔｚiと焦点距離ｆiは、
【００２０】
【数５】
2006106930.tif 000007
【００２１】
【数６】
2006106930.tif 000008
【００２２】
で得られる。
【００２３】
さらに、画素の横縦のスケールファクタ比が１：βの場合は、
【００２４】
【数７】
2006106930.tif 000009
【００２５】
が投影モデルであり、各フレームでの焦点距離ｆiとして、αi＝ｆi，βi＝βｆiとする。ここで、カメラ運動と３次元形状を復元されると、再投影座標値（Δｕij，Δｖij）を算出し、式７から、
【００２６】
【数８】
2006106930.tif 000010
【００２７】
【数９】
2006106930.tif 000011
【００２８】
【数１０】
2006106930.tif 000012
【００２９】
となる。上式を全特徴点について連立して、
【００３０】
【数１１】
2006106930.tif 000013
【００３１】
となるため、
【００３２】
【数１２】
2006106930.tif 000014
【００３３】
の計算により、第ｉフレームでの光軸並進運動Ｔｚi，αi，βiを同時に求めることができる。
【００３４】
本発明では、平面運動に光軸並進を加えた運動で撮像した時系列画像からカメラ運動と３次元形状を復元し、カメラ内部パラメータを推定するため、カメラを図１、図２に示すように平面板上に支柱で固定されたカメラを半径Ｒの円周上で移動させることにより、平面運動に光軸並進運動を加えた運動での時系列画像を取得でき、その時系列画像中の特徴点の画像座標値（ｘij，ｙij）を観測することで、式２、式７の投影モデルに基づき、オートフォーカスによる焦点距離の変化を考慮した式５、式６、式１２の計算によりカメラ内部パラメータを推定することができる。
【発明の効果】
【００３５】
本発明により、カメラ等を使って取得した時系列画像全般（移動手段を利用して撮影した車載画像、海上画像、空撮画像、屋内画像など）から、オートフォーカス機能が作用した状況でも、カメラ運動と３次元形状を同時に復元し、オートフォーカス機能により変化したカメラ内部パラメータを推定することができる。
【００３６】
本発明ではカメラが図１、図２の撮像手段により時系列画像を取得することで、特徴点観測手段での特徴点の画像座標値を高精度に観測することに効果がある。これにより、カメラ運動と３次元形状、カメラ内部パラメータをより高精度に復元、推定することができる。なお、本発明で使用する計算は、大半が線形演算で構成されるため、コンピュータ言語での実装が容易である。
【発明を実施するための最良の形態】
【００３７】
以下、図面を用いて本発明の実施形態を説明する。
【００３８】
(実施形態１)
図１では水平面板１１に回転中心から半径Ｒの距離を隔て、垂直に立てた支柱１２にカメラ１３を固定し、カメラ１３を半径Ｒの円周上を回転させながら、被写体を撮像する。
【００３９】
また、図２では垂直に立てた回転平面板２１に回転中心から半径Ｒの距離を隔て、水平に立てた支柱２２にカメラ２３を固定し、カメラ２３を半径Ｒの円周上を回転させながら、被写体を撮像する。
【００４０】
これらの撮像手段では、カメラを全周回転させる必要はなく、被写体がフレームアウトしない程度で撮影を行い、時系列画像を得る。これにより、平面運動に光軸並進運動のある自由度４のカメラ運動で撮像することができる。一方、カメラを手動で動かしながら被写体を撮影した時系列画像にも本発明を適用することができるが、図１、図２の撮影手段で得ることで、より高精度にカメラ内部パラメータを推定することができる。
【００４１】
図３は図１、図２の平面運動に光軸方向の並進運動があるカメラ運動での撮影手段により、時系列画像を取得する状況のカメラと対象物の位置関係を表している。運動の中心は視点としており、視点を原点としたカメラ座標系ＸＹＺ，原点Ｏとした世界座標系ＸｗＹｗＺｗを設定する。この座標系において、カメラ運動として、Ｔi（Ｔｘi，Ｔｙi，Ｔｚi）の並進運動と光軸（Ｚ軸）周りのヨー回転θiとする。対象物の点Ｐj（Ｘj，Ｙj，Ｚj）はカメラにより画像面において投影中心（主点）を原点とした画像座標値（ｘij，ｙij）へ投影されるとする。なお、説明の都合上、初期フレームでの視点とＯは一致しているとし、光軸はＺｗ軸と平行関係にあり、θiはＸとＸｗ軸との成す角とするが、一般性を損うものではない。
【００４２】
図４は本発明のパラメータ推定装置の基本構成図であり、本実施形態を図４にそって説明する。本発明は、時系列画像を蓄積した時系列画像データベース４１、各画像における特徴点の座標値を観測する特徴点観測部４２、画像座標値からなる行列データからカメラ運動と３次元形状を復元する因子分解処理部４３、光軸方向の並進運動を復元し、かつ、カメラ内部パラメータの焦点距離を推定する光軸並進運動復元・焦点距離推定部４４、カメラ運動と３次元形状、焦点距離が得られるまで処理を繰り返すための反復処理部４５から構成される。
【００４３】
この構成において、時系列画像データベース４１には、ハードディスク、ＲＡＩＤ装置、ＣＤ−ＲＯＭなどの記録媒体を利用する形態、または、ネットワークを介したリモートなデータ資源を利用する形態でもどちらでも構わない。
【００４４】
さらに、図５はリアルタイムで処理する場合の処理構成図であり、本発明は必ずしも時系列画像データベース部４１などの記憶装置を必要としない。
【００４５】
まず、図４の特徴点観測部４２において、対象物を撮影した時系列画像として時系列画像データベース４１からフレームＦ分の画像系列を取り出す。この取り出した画像系列において特徴点追跡を行う。特徴点は従来から用いられているような以下の手順により抽出する。
【００４６】
図６に示すように画像１の領域１において、
（１）各画素に対する２×２のヘッセ行列を求める。次に、
（２）各点の３×３近傍領域において極大点かどうか判定し、極大点以外の点を削除する（ｎｏｎ−ｍａｘｉｍａ  ｓｕｐｐｒｅｓｓｉｏｎ）。さらに、
（３）得られた各点のヘッセ行列の固定値σl，σs（σs≦σl）を求め、σsがσp以上となる点を抽出する。最後に、
（４）抽出した点のσsの大きさの順にソートし、上位の点から順番にその点（ｐl）より上位の点（ｐh）がσd画素以内の距離に存在するかどうかを判定し、もし、存在する場合は下位の点ｐlを削除する。さらに、抽出した特徴点（ｊ＝１，２，・・・，Ｐ）をＫＬＴ法（Ｋａｎａｄｅ−Ｌｕｃａｓ−Ｔｏｍａｓｉ）により画像ｉ（ｉ＝２，・・・，Ｆ）にわたって追跡し、画像座標値（ｘij，Ｙij）を観測する。このようにして得られた特徴の画像座標値を
【００４７】
【数１３】
2006106930.tif 000015
【００４８】
に示す配列に並べた２Ｆ×Ｐの行列データ（行列データ［Ａ］）を用意する。
【００４９】
次に、因子分解処理部４３、光軸並進運動復元・焦点距離推定部４４、及び、反復処理部４５の処理内容について、図７の処理フロー図にそって説明する。
【００５０】
まず、因子分解処理部４３が行列データの入力を行う(Ｓ１)。
【００５１】
次に、因子分解処理部４３が行列分解データの生成を行う(Ｓ２)。すなわち、係数δijを
【００５２】
【数１４】
2006106930.tif 000016
【００５３】
のように初期化し、保持しておいた式１３の２Ｆ×Ｐの行列データ［Ａ］から
【００５４】
【数１５】
2006106930.tif 000017
【００５５】
に示すデータを求めておき、これを行列分解データ［Ｂ］として保持しておく。
【００５６】
次に、因子分解処理部４３が、特異値分解を行う（Ｓ３）。すなわち、行列分解データ［Ｂ］を入力データとして特異値分解により
【００５７】
【数１６】
2006106930.tif 000018
【００５８】
に示す３つの行列［Ｕ］，［Ｗ］，［Ｖ］に行列分解する。ここで、［Ｕ］は２Ｆ×Ｐサイズの行列、［Ｗ］はＰ×Ｐサイズの対角行列、［Ｖ］はＰ×Ｐサイズの行列である。
【００５９】
次に、因子分解処理部４３が雑音除去を行う（Ｓ４）。すなわち、
【００６０】
【数１７】
2006106930.tif 000019
【００６１】
の第二項に示すように、ランク４以降の各行列の成分を削除する。この削除のときは、行列［Ｕ］を取り出し、この行列の要素において第４から第Ｐ列目までを削除し、残りの成分からなる行列を保持する。次に、行列［Ｗ］を取り出し、この行列の要素において第４から第Ｐ行目並びに第４から第Ｐ列目までを削除し、残りの成分からなる行列を保持する。さらに、行列［Ｖ］を取り出し、この行列の要素において第４から第Ｐ行目までを削除し、残りの成分からなる行列を保持する。
【００６２】
【数１８】
2006106930.tif 000020
【００６３】
雑音除去に続き、第４から第Ｐ行目並びに第４から第Ｐ列目までを削除した行列［Ｗ］の対角要素の平方をとった行列から、
【００６４】
【数１９】
2006106930.tif 000021
【００６５】
【数２０】
2006106930.tif 000022
【００６６】
に示す行列［Ｕ’］と行列［Ｖ’］を得る。
【００６７】
次に、因子分解処理部４３が変換行列算出を行う(Ｓ５)。すなわち、保持してある行列［Ｕ’］を取り出し、
【００６８】
【数２１】
2006106930.tif 000023
【００６９】
【数２２】
2006106930.tif 000024
【００７０】
【数２３】
2006106930.tif 000025
【００７１】
に示す連立する条件式における対称行列［Ｃ］の各要素にかかる係数を計算する。なお、対称行列［Ｃ］は
【００７２】
【数２４】
2006106930.tif 000026
【００７３】
で表されている。
【００７４】
これらの係数計算は、行列演算により容易に得られるものであり、この条件式を全フレームに対して計算する。次に、全フレームの式２１〜２３に示す連立する条件式に対して、最小二乗法などの数値計算を利用して３×３サイズの行列［Ｃ］の各要素を決定する。求めた行列［Ｃ］を
【００７５】
【数２５】
2006106930.tif 000027
【００７６】
に示すように固有値分解する。ここで、固有値行列の平方と固有値行列から、
【００７７】
【数２６】
2006106930.tif 000028
【００７８】
の行列［Ｃ’］を生成し、この行列要素を成分にもつ行列［Ｑ］を
【００７９】
【数２７】
2006106930.tif 000029
【００８０】
に従って算出する。
【００８１】
次に、因子分解処理部４３が平面運動復元を行う(Ｓ６)。すなわち、求めた行列［Ｑ］と、保持しておいた行列［Ｕ’］から、
【００８２】
【数２８】
2006106930.tif 000030
【００８３】
の行列演算により行列［Ｍ’］を計算する。行列［Ｍ’］から各フレーム（第ｉフレーム）の行列要素（ｍix，ｎix）または（ｍiy，ｎiy）を取り出し、
【００８４】
【数２９】
2006106930.tif 000031
【００８５】
を使って、ヨー回転θiを復元する。次に、行列［Ｍ’］から各フレーム（第ｉフレーム）の行列要素（Ｔiu，Ｔiv）を取り出す。この（Ｔiu，Ｔiv）から、
【００８６】
【数３０】
2006106930.tif 000032
【００８７】
を使って第ｉフレームにおけるユークリッド空間でのＸＹ並進運動（Ｔｘi，Ｔｙi）を計算する。
【００８８】
次に、因子分解処理部４３が３次元情報復元を行う(Ｓ７)。すなわち、先に保持しておいた行列［Ｖ’］と、変換行列算出で得られた行列［Ｑ］から、
【００８９】
【数３１】
2006106930.tif 000033
【００９０】
に示す行列演算を行い、行列［Ｓ’］を求める。次に、行列［Ｓ’］の要素に対して、
【００９１】
【数３２】
2006106930.tif 000034
【００９２】
に示す変換を行い行列［Ｐ］を得る。行列［Ｐ］の列ベクトルは、それぞれ第ｊ番目の特徴点のユークリッド空間での３次元座標値（Ｘj，Ｙj，Ｚj）になっている。
【００９３】
次に、光軸並進運動復元・焦点距離推定部４４が光軸並進運動復元・焦点距離推定を行う(Ｓ８)。すなわち、平面運動復元と３次元情報復元で復元したカメラ運動と３次元形状から
【００９４】
【数３３】
2006106930.tif 000035
【００９５】
に示す（Δｕij，Δｖij）を得て、各フレームでの画像座標値（ｘij，ｙij）から、式３に示すΔｗijを全フレーム、全特徴点について求め、
【００９６】
【数３４】
2006106930.tif 000036
【００９７】
のＦ×Ｐの行列［Ｄ1］を準備する。次に、３次元情報復元で得たＺj（ｊ＝１，２，・・・，Ｐ）からなる
【００９８】
【数３５】
2006106930.tif 000037
【００９９】
の行列［Ｄ2］を準備する。
【０１００】
【数３６】
2006106930.tif 000038
【０１０１】
の行列演算を行い、
【０１０２】
【数３７】
2006106930.tif 000039
【０１０３】
のＦ×２の行列［Ｄ］を得る。さらに、行列［Ｄ］の第ｉ行第１列のｄi1，第ｉ行第２列のｄi2と取り出し、
【０１０４】
【数３８】
2006106930.tif 000040
【０１０５】
の計算を全フレームにわたり計算して、焦点距離ｆi、光軸並進運動Ｔｚiを得る。この計算は、式５、式６に基づいて行列計算を行っている。
【０１０６】
次に、反復処理部４５が、
【０１０７】
【数３９】
2006106930.tif 000041
【０１０８】
【数４０】
2006106930.tif 000042
【０１０９】
の計算で得られる誤差△Ｅを算出する(Ｓ９)。この誤差が、ある一定値以下に収束しているかを判断し(Ｓ１０)、ある一定値以下に収束している場合は処理を終了する。もし、そうでない場合は、
【０１１０】
【数４１】
2006106930.tif 000043
【０１１１】
に従って係数δijを更新する(Ｓ１１)。さらに、この係数更新に伴い式１５により行列分解データ［Ｂ］を更新する。これを誤差△Ｅが一定値以下に収束するまで反復を繰り返す。
【０１１２】
以上、本発明の実施例により、時系列画像中の特徴点の画像座標値から、カメラ視点の運動：光軸周りの回転とＸＹＺ軸方向の並進運動、及び、物体形状を構成する３次元情報を復元し、同時に、カメラ内部パラメータの焦点距離を推定することができる。
【０１１３】
(実施形態２)
図８は本発明のパラメータ推定装置の基本構成図であり、本実施形態を図８にそって説明する。本実施形態は、時系列画像を蓄積した時系列画像データベース８１、各画像における特徴点の座標値を観測する特徴点観測部８２、画像座標値からなる行列データからカメラ運動と３次元形状を復元する因子分解処理部８３、光軸方向の並進運動を復元し、カメラ内部パラメータの焦点距離を推定する光軸並進運動復元・焦点距離推定部８４、カメラ運動と３次元形状、焦点距離が得られるまで処理を繰り返すための反復処理部８５、カメラ内部パラメータの画素の横縦のスケールファクタ比を求めるスケールファクタ比算出部８６から構成される。
【０１１４】
この構成において、時系列画像データベース８１には、ハードディスク、ＲＡＩＤ装置、ＣＤ−ＲＯＭなどの記録媒体を利用する形態、または、ネットワークを介したリモートなデータ資源を利用する形態でもどちらでも構わない。
【０１１５】
さらに、図９はリアルタイムで処理する場合の処理構成図であり、本実施形態は必ずしも時系列画像データベース部８１などの記憶装置を必要としない。
【０１１６】
本実施形態は、実施形態１に対して光軸並進運動復元・焦点距離推定部８４、スケールファクタ比算出部８６の処理が異なっているため、以下の説明では、この部分についてだけ記載する。
【０１１７】
図１０は図８に示すパラメータ推定装置における因子分解処理部８３、光軸並進運動復元・焦点距離推定部８４、スケールファクタ比算出部８６の処理フローである。
【０１１８】
最初の反復では、
【０１１９】
【数４２】
2006106930.tif 000044
【０１２０】
のように係数δij，γijを初期化し、
【０１２１】
【数４３】
2006106930.tif 000045
【０１２２】
の行列［Ｂ］を得る。光軸並進運動復元・焦点距離推定部８４での処理は、式３３に示す（Δｕij，Δｖij）を得て、各フレームでの画像座標値（ｘij，ｙij）から、全特徴点について式８、式９、式１０のＡj，Ｂjを得て、
【０１２３】
【数４４】
2006106930.tif 000046
【０１２４】
の２×Ｐの行列［Ｄ1］を準備する。式３６の行列演算で
【０１２５】
【数４５】
2006106930.tif 000047
【０１２６】
の各フレームでの２×２の行列［Ｄ］を得る。行列［Ｄ］から第１行第１列のｄi1、第２行第１列のｄ’i1を取り出し
【０１２７】
【数４６】
2006106930.tif 000048
【０１２８】
によりαi，βiを求める。また、光軸並進運動Ｔｚiは行列［Ｄ］の各要素から
【０１２９】
【数４７】
2006106930.tif 000049
【０１３０】
により計算する。これらの計算は式１２に基づいて計算を行っている。
【０１３１】
反復処理部８５では、式３９、式４０の計算で得られる△Ｅを算出する。この誤差が、ある一定値以下に収束している場合は処理を終了する。もし、そうでない場合は、係数δij，γijを
【０１３２】
【数４８】
2006106930.tif 000050
【０１３３】
に従って更新し、この更新に伴い式４３により行列分解データ［Ｂ］を更新する。これを誤差△Ｅが一定値以下に収束するまで反復を繰り返す。
【０１３４】
反復処理が終了したとき、スケールファクタ比算出部８６ではβi／αiの平均値を求めて画素の横縦のスケールファクタ比１：βを得る。
【０１３５】
以上、実施形態２により、時系列画像中の特徴点の画像座標値から、カメラ視点の運動：光軸周りの回転とＸＹＺ軸方向の並進運動、並びに、物体形状を構成する３次元情報を復元し、同時に、カメラ内部パラメータの焦点距離、画素の横縦の量子化比（スケールファクタ比）を推定することができる。
【０１３６】
なお、実施形態１及び２に記載のパラメータ推定装置は、例えば、パラメータ推定装置を構成するコンピュータ装置が有するＣＰＵによって実現され、必要とする特徴点観測、因子分解処理、光軸並進運動復元・焦点距離推定、反復処理、スケールファクタ比算出などをアプリケーションプログラムとして搭載することができる。
【０１３７】
また、コンピュータ装置には上記説明した特徴点観測、因子分解処理、光軸並進運動復元・焦点距離推定、反復処理、スケールファクタ比算出など行った処理結果や計算結果等のデータを内部メモリや外部記憶装置等に書き込み・読み出しができるようにしてもよい。
【０１３８】
また、本発明は、前述した実施形態の機能を実現するソフトウェアのプログラムコードを記録した記録媒体を、システムあるいは装置に供給し、そのシステムあるいは装置のＣＰＵ(ＭＰＵ)が記録媒体に格納されたプログラムコードを読み出し実行することによっても実現できる。その場合、記憶媒体から読み出されたプログラムコード自体が上述した実施形態の機能を実現することになり、そのプログラムコードを記憶した記憶媒体、例えば、ＣＤ−ＲＯＭ、ＤＶＤ−ＲＯＭ、ＣＤ−Ｒ、ＣＤ−ＲＷ、ＭＯ、ＨＤＤ、等は本発明を構成する。
【図面の簡単な説明】
【０１３９】
【図１】カメラが水平面板上で半径Ｒの回転運動で撮像する場合の図。
【図２】カメラが垂直面板上で半径Ｒの回転運動で撮像する場合の図。
【図３】カメラ視点と３次元位置、並びに、本発明の説明で使用する座標系を示す図。
【図４】パラメータ推定装置の基本構成図。
【図５】パラメータ推定装置の基本構成図。
【図６】時系列画像での特徴点の観測を示す図。
【図７】因子分解処理部、光軸並進運動復元・焦点距離推定部、及び、反復処理部に関する処理フロー図。
【図８】パラメータ推定装置の基本構成図。
【図９】パラメータ推定装置の基本構成図。
【図１０】因子分解処理部、光軸並進運動復元・焦点距離推定部、反復処理部、及び、スケールファクタ比算出部に関する処理フロー図。
【符号の説明】
【０１４０】
１１…平面版
１２…支柱
１３…カメラ
２１…回転平面版
２２…支柱
２３…カメラ
４１…時系列画像データベース
４２…特徴点観測部
４３…因子分解処理部
４４…光軸並進運動復元・焦点距離推定部
４５…反復処理部
５１…画像入力部
５２…特徴点観測部
５３…因子分解処理部
５４…光軸並進運動復元・焦点距離推定部
５５…反復処理部
８１…時系列画像データベース
８２…特徴点観測部
８３…因子分解処理部
８４…光軸並進運動復元・焦点距離推定部
８５…反復処理部
８６…スケールファクタ比算出部
９１…画像入力部
９２…特徴点観測部
９３…因子分解処理部
９４…光軸並進運動復元・焦点距離推定部
９５…反復処理部
９６…スケールファクタ比算出部
【図１】
2006106930.tif 000051
【図２】
2006106930.tif 000052
【図３】
2006106930.tif 000053
【図４】
2006106930.tif 000054
【図５】
2006106930.tif 000055
【図６】
2006106930.tif 000056
【図７】
2006106930.tif 000057
【図８】
2006106930.tif 000058
【図９】
2006106930.tif 000059
【図１０】
2006106930.tif 000060