JP7052788B2 - カメラパラメータ推定装置、カメラパラメータ推定方法、及びプログラム - Google Patents

カメラパラメータ推定装置、カメラパラメータ推定方法、及びプログラム Download PDF

Info

Publication number
JP7052788B2
JP7052788B2 JP2019505764A JP2019505764A JP7052788B2 JP 7052788 B2 JP7052788 B2 JP 7052788B2 JP 2019505764 A JP2019505764 A JP 2019505764A JP 2019505764 A JP2019505764 A JP 2019505764A JP 7052788 B2 JP7052788 B2 JP 7052788B2
Authority
JP
Japan
Prior art keywords
camera
optimum solution
parameters
parameter
calculated
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2019505764A
Other languages
English (en)
Other versions
JPWO2018168255A1 (ja
Inventor
学 中野
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
NEC Corp
Original Assignee
NEC Corp
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by NEC Corp filed Critical NEC Corp
Publication of JPWO2018168255A1 publication Critical patent/JPWO2018168255A1/ja
Application granted granted Critical
Publication of JP7052788B2 publication Critical patent/JP7052788B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/80—Analysis of captured images to determine intrinsic or extrinsic camera parameters, i.e. camera calibration
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
    • G01B11/00—Measuring arrangements characterised by the use of optical techniques
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/70—Determining position or orientation of objects or cameras
    • H—ELECTRICITY
    • H04—ELECTRIC COMMUNICATION TECHNIQUE
    • H04N—PICTORIAL COMMUNICATION, e.g. TELEVISION
    • H04N17/00—Diagnosis, testing or measuring for television systems or their details
    • H04N17/002—Diagnosis, testing or measuring for television systems or their details for television cameras
    • H—ELECTRICITY
    • H04—ELECTRIC COMMUNICATION TECHNIQUE
    • H04N—PICTORIAL COMMUNICATION, e.g. TELEVISION
    • H04N23/00—Cameras or camera modules comprising electronic image sensors; Control thereof
    • H04N23/60—Control of cameras or camera modules
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10004—Still image; Photographic image
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10016—Video; Image sequence
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30244—Camera pose

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Theoretical Computer Science (AREA)
  • Multimedia (AREA)
  • Signal Processing (AREA)
  • Health & Medical Sciences (AREA)
  • Biomedical Technology (AREA)
  • General Health & Medical Sciences (AREA)
  • Studio Devices (AREA)

Description

本発明は、カメラの外部パラメータと内部パラメータとを推定する、カメラパラメータ推定装置及びカメラパラメータ推定方法に関し、更に、これらを実現するためのプログラムに関する。
既知の３次元点を観測した画像からカメラパラメータを推定する技術は、ロボットの自己位置推定、および画像合成において、重要な要素技術である（例えば、非特許文献１～５参照）。また、カメラパラメータには、外部パラメータと、内部パラメータとの２種類が存在する。外部パラメータは、カメラの３次元空間における位置と姿勢とを表すパラメータであり、内部パラメータは、レンズの焦点距離、光学中心、アスペクト比、せん断係数、レンズ歪み係数などのパラメータである。以下では、既知の３次元点と、この３次元点を画像上で観測することによって得られた２次元点とのペアを、単に「対応点」と表記する。
非特許文献１と非特許文献２とは、共に、表面模様の３次元点が既知の平面、または立体のパターンを複数の角度から撮影して得られた画像を用いて、外部パラメータと内部パラメータとを計算する方法を開示している。これらの非特許文献に開示された方法では、まずレンズ歪み係数がすべてゼロに仮定され、姿勢に関する拘束条件を無視した上で内部パラメータと外部パラメータとが計算される。そして、計算された内部パラメータおよび外部パラメータを初期値として、ニュートン法などの非線形最適化が行われる。
また、非特許文献３は、非特許文献１または非特許文献２に開示された方法を利用して、予め内部パラメータを計算しておき、３組以上の対応点から外部パラメータのみを計算する方法を開示している。具体的には、非特許文献３に開示された方法は、３次元点と２次元点との射影関係に基づく誤差関数を多項式で表現し、誤差関数の停留点である勾配のゼロ点をすべて計算することで、局所解を回避して大域的に最適な唯一の解を推定する。
更に、非特許文献４は、内部パラメータのうち焦点距離のみを未知として、４組以上の対応点を用いて、非特許文献３に開示された方法と同様にして、多項式表現による大域的最適解を直接計算する方法を開示している。
また、非特許文献５および非特許文献６は、それぞれ、４組または５組の対応点を入力として、内部パラメータのうち焦点距離とレンズ歪み係数を未知として、これらの内部パラメータと外部パラメータとを共に算出する方法を開示している。
更に、非特許文献７は、内部パラメータのうち焦点距離とレンズ歪み係数とを未知として、５組以上の対応点を用いて内部パラメータと外部パラメータとを算出する方法を開示している。具体的には、誤差関数をカメラパラメータの依存関係に基づいて複数の部分問題に分割し、それぞれの部分問題において、非特許文献３に開示された方法と同様にして、多項式表現による大域的最適解を計算する。
Z. Zhang,"Flexible camera calibration by viewing a plane from unknown orientations,"ICCV 1999. R. Tsai, "Aversatile camera calibration technique for high-accuracy 3D machine visionmetrology using off-the-shelf TV cameras and lenses,"Robotics and Automation, IEEEJournal of 3.4 (1987): 323-344. G. Nakano, "Globally Optimal DLS Method for PnP Problem with Cayleyparameterization," BMVC 2015. Y. Zheng, etal., "A General and Simple Method for Camera Pose andFocal Length Determination," CVPR 2014. M. Bujnak, etal., "New efficient solution to the absolute poseproblem for camera with unknown focal length and radial distortion," ACCV 2010. Z. Kukelova, etal., "Real-time solution to the absolute pose problemwith unknown radial distortion and focal length," ICCV2013. G. Nakano, "A Versatile Approach for Solving PnP, PnPf, and PnPfr Problems," ECCV 2016.
しかしながら、上述の非特許文献１～７に開示された方法には、以下のような問題がある。まず、非特許文献１および２に開示された方法では、安定的な計算には数十枚の画像が必要であり、また、非線形最適化の際に望ましくない局所解に収束することがあり、内部パラメータおよび外部パラメータの推定精度が低いという問題がある。
また、非特許文献３に開示された方法では、固定焦点距離のカメラが想定されており、ズームレンズを搭載したカメラが用いられた場合に、内部パラメータおよび外部パラメータを推定できないという問題がある。一方、非特許文献４に開示された方法では、ズームレンズを搭載したカメラのズームを変化させた際の焦点距離は推定できるものの、ズームの変化と同時に生じるレンズ歪み係数の変化に対応できない。よって、非特許文献４に開示された方法にも、内部パラメータおよび外部パラメータの推定精度が低いという問題がある。
更に、非特許文献５及び非特許文献６に開示された方法では、一般に複数の解が得られ、一意に解を選択できず、加えて、４組以上または５組以上の対応点を入力として最小二乗的な意味で最適な解を推定することもできない。従って、非特許文献５または非特許文献６に開示された方法にも、内部パラメータおよび外部パラメータの推定精度が低いという問題がある。
加えて、非特許文献７に開示された方法では、５組以上の対応点を入力とするため最小二乗的な意味で最適な解を算出できる。しかし、非特許文献７に開示された方法では、多項式表現に基づく数値計算方法がノイズに敏感であるため、本来実数である望ましい解が複素数になってしまい、誤った大域的最適解が選択されてしまう問題がある。
このように、上述した非特許文献１～６に開示された方法には、入力である複数の対応点から、レンズ歪み係数を未知として含む内部パラメータと外部パラメータとを高精度に推定することは難しいという問題がある。
本発明の目的の一例は、上記問題を解消し、レンズ歪み係数を含む内部パラメータと外部パラメータとを、数値的に安定に、且つ、高精度に推定し得る、カメラパラメータ推定装置、カメラパラメータ推定方法、およびプログラムを提供することにある。
上記目的を達成するため、本発明の一側面におけるカメラパラメータ推定装置は、カメラによる対象物の撮影時におけるカメラパラメータを推定するための装置であって、
前記対象物に関する６組以上の３次元点と、前記３次元点それぞれに対応する画像上の２次元点との変換を表す誤差関数を、前記カメラパラメータの依存関係と前記３次元点の分布の平面性とに基づいて、複数の部分問題に分割し、前記複数の部分問題それぞれにおいて最適解候補を算出する、パラメータ算出部と、
前記最適解候補を初期値として、前記誤差関数によって得られる誤差を最小とする最適解を求め、求めた前記最適解を最適なカメラパラメータとして出力する、パラメータ最適化部と、
を備えている、ことを特徴とする。
また、上記目的を達成するため、本発明の一側面におけるカメラパラメータ推定方法は、カメラによる対象物の撮影時におけるカメラパラメータを推定するための方法であって、
（ａ）前記対象物に関する６組以上の３次元点と、前記３次元点それぞれに対応する画像上の２次元点との変換を表す誤差関数を、前記カメラパラメータの依存関係と前記３次元点の分布の平面性とに基づいて、複数の部分問題に分割し、前記複数の部分問題それぞれにおいて最適解候補を算出する、ステップと、
（ｂ）前記最適解候補を初期値として、前記誤差関数によって得られる誤差を最小とする最適解を求め、求めた前記最適解を最適なカメラパラメータとして出力する、ステップと、
を有する、ことを特徴とする。
更に、上記目的を達成するため、本発明の一側面におけるプログラムは、コンピュータによって、カメラによる対象物の撮影時におけるカメラパラメータを推定するためのプログラムであって、
前記コンピュータに、
（ａ）前記対象物に関する６組以上の３次元点と、前記３次元点それぞれに対応する画像上の２次元点との変換を表す誤差関数を、前記カメラパラメータの依存関係と前記３次元点の分布の平面性とに基づいて、複数の部分問題に分割し、前記複数の部分問題それぞれにおいて最適解候補を算出する、ステップと、
（ｂ）前記最適解候補を初期値として、前記誤差関数によって得られる誤差を最小とする最適解を求め、求めた前記最適解を最適なカメラパラメータとして出力する、ステップと、
を実行させる、プログラム。
以上のように、本発明によれば、レンズ歪み係数を含む内部パラメータと外部パラメータとを、数値的に安定に、且つ、高精度に推定することができる。
図１は、本発明の実施の形態におけるカメラパラメータ推定装置の概略構成を示すブロック図である。 図２は、本発明の実施の形態におけるカメラパラメータ推定装置の具体的構成を示すブロック図である。 図３は、本発明の実施の形態におけるカメラパラメータ推定装置の動作を示すフロー図である。 図４は、本発明の実施の形態におけるカメラパラメータ推定装置を実現するコンピュータの一例を示すブロック図である。
（実施の形態１）
以下、本発明の実施の形態における、カメラパラメータ推定装置、カメラパラメータ推定方法、及びプログラムについて、図１～図４を参照しながら説明する。
［装置構成］
最初に、本発明の実施の形態におけるカメラパラメータ推定装置の概略構成及び機能について図１を用いて説明する。図１は、本発明の実施の形態におけるカメラパラメータ推定装置の概略構成を示すブロック図である。
図１に示す本実施の形態におけるカメラパラメータ推定装置１０は、カメラによる対象物の撮影時におけるカメラパラメータを推定する装置である。図１に示すように、カメラパラメータ推定装置１０は、パラメータ算出部２０と、パラメータ最適化部３０とを備えている。
パラメータ算出部２０は、まず、対象物に関する６組以上の３次元点と、３次元点それぞれに対応する画像上の２次元点との変換を表す誤差関数を、カメラパラメータの依存関係と３次元点の分布の平面性に基づいて、複数の部分問題に分割する。次いで、パラメータ算出部２０は、複数の部分問題それぞれを線形最小二乗法により解き、最適解候補を算出する。
パラメータ最適化部３０は、まず、最適解候補を初期値として誤差関数によって得られる誤差を最小とする最適解を求める。次いで、パラメータ最適化部３０は、求めた最適解を最適なカメラパラメータとして出力する。
このように、本実施の形態では、誤差関数は複数の部分問題に分割され、部分問題毎に線形最小二乗法に基づき最適パラメータ候補が計算され、そして、計算された最適パラメータ候補を用いて最適解が計算される。よって、本実施の形態によれば、カメラパラメータ、即ち、レンズ歪み係数を含む内部パラメータと外部パラメータとが、数値的に安定に、且つ、高精度に推定される。
続いて、図２を用いて、本発明の実施の形態におけるカメラパラメータ推定装置の構成及び機能について、更に具体的に説明する。図２は、本発明の実施の形態におけるカメラパラメータ推定装置の具体的構成を示すブロック図である。
図２に示すように、本実施の形態では、カメラパラメータ推定装置１０には、撮影の対象物に関する６組以上の３次元点と、３次元点それぞれに対応する画像上の２次元点とが、上述した対応点として入力される。
本実施の形態において、カメラパラメータには、カメラの外部に関する外部パラメータと、カメラの内部に関する内部パラメータとが含まれる。具体的には、外部パラメータは、従来と同様に、カメラの３次元区間における位置及び姿勢のいずれか一方または両方を含む。内部パラメータは、カメラにおける、レンズの焦点距離、光学中心、アスペクト比、せん断係数、及びレンズ歪み係数のうち、少なくとも１つを含む。なお、本実施の形態では、カメラパラメータには、外部パラメータのみが含まれていても良い。
本実施の形態では、カメラパラメータ推定装置１０は、６組以上の対応点に基づいて、カメラの外部パラメータ、更にはカメラの外部パラメータに加えて内部パラメータも推定する情報処理装置である。また、本実施の形態では、図２に示すように、パラメータ算出部２０には、分割処理部２１と、部分パラメータ算出部２２と、残存パラメータ算出部２３とが備えられている。
分割処理部２１は、対応点に含まれる３次元点と２次元点との所定の変換を表す誤差関数を、カメラパラメータの依存関係と３次元点の分布の平面性とに基づいて、複数の部分問題に分割する。
ここで「カメラパラメータ依存関係に基づく分割」とは、誤差関数を構成する複数の項のうち、例えば外部パラメータのみが変数として現れる項、または、内部パラメータのみが変数として現れる項のように、未知のパラメータ同士の影響を受けない、又は影響が少ない項を選出することを表す。
また、「３次元点の分布の平面性に基づく分割」とは、入力に含まれる３次元点が、１つの平面上に分布しているか、２つ以上の平面、つまりは非平面に分布しているかに応じて、部分問題の定式化を切り替えることを表す。
部分パラメータ算出部２２は、対応点を入力として、複数の部分問題のうち、外部パラメータの一部のみに依存する１つの部分問題について、線形最小二乗法に基づき、外部パラメータの一部の最適解候補を算出する。
残存パラメータ算出部２３は、部分パラメータ算出部２２で算出された外部パラメータの一部の最適解候補と、対応点とを入力として、未だ算出されていないパラメータの部分問題について線形最小二乗法に基づき解を算出する。
パラメータ最適化部３０は、本実施の形態では、部分パラメータ算出部２２と残存パラメータ算出部２３とで算出された最適解候補と対応点とを入力として、これらに基づいて、分割される前の元の誤差関数を最適化することによって、誤差を最小とする最適解を求める。そして、パラメータ最適化部３０は、求めた最適解を最適なカメラパラメータとして出力する。
［装置動作］
次に、本発明の実施の形態におけるカメラパラメータ推定装置１０の動作について図３を用いて説明する。図３は、本発明の実施の形態におけるカメラパラメータ推定装置の動作を示すフロー図である。以下の説明においては、適宜図１を参酌する。また、本実施の形態では、カメラパラメータ推定装置１０を動作させることによって、カメラパラメータ推定方法が実施される。よって、本実施の形態におけるカメラパラメータ推定方法の説明は、以下のカメラパラメータ推定装置１０の動作説明に代える。
図３に示すように、まず、カメラパラメータ推定装置１０に６組以上の対応点が入力されると、分割処理部２１は、対応点に含まれる３次元点と２次元点との所定の変換を表す誤差関数を、カメラパラメータの依存関係と３次元点の分布の平面性とに基づいて、複数の部分問題に分割する（ステップＳ１１）。
本実施の形態において、入力される対応点は、推定すべきカメラパラメータにおける、いわゆる退化配置（Degenerate ConfigurationまたはCriticalConfigurationとも呼ばれる）ではないとする。なぜならば、そのような対応点に対しては理論的にカメラパラメータを推定できないためである。また、座標値に関して、３次元点も２次元点も共にすべて異なるとする。これは、例えば、１つの３次元点が複数の異なる２次元点に対応することは現実的に起こりえないためである。カメラパラメータ推定装置１０において、これらの不適切な対応点は以降の処理から除外される。また、本実施の形態は、カメラパラメータ推定装置１０の代わりに、対応点を外部からカメラパラメータ推定装置に入力する別の装置（図示せず）が、対応点のエラー判定等を行う態様であっても良い。
次に、部分パラメータ算出部２２は、ステップＳ１１の分割によって得られた部分問題の中で外部パラメータの一部のみに依存する１つの部分問題について、線形最小二乗法に基づき最適解候補を算出する（ステップＳ１２）。
ここで、最適解候補とは、最小化したい誤差関数の最小の停留点に近い値、つまり、誤差関数を最小化するための初期として適切なことを意味する。また、停留点とは、最適化問題における目的関数（より正確にはラグランジェ関数であるが、以下単に「目的関数」と表記する）の勾配ベクトルがゼロとなる未知数の点のことを意味する。すなわち、停留点は、勾配ベクトルを連立方程式とみなすと、その連立方程式の解のことである。ステップＳ１２においては、目的関数は外部パラメータの一部のみに依存する１つの部分問題を指し、停留点は当該部分問題の勾配ベクトルをゼロとする外部パラメータを指す。また、最適解候補は、当該停留点ではないが、それに近しい値を指す。
次いで、残存パラメータ算出部２３は、対応点と、ステップＳ１２で算出された外部パラメータの一部の最適解候補とを入力として、残りのカメラパラメータを未知数とするその他の部分問題、即ち、最適解候補が未だ算出されていない部分問題における、すべての最適解候補を算出する（ステップＳ１３）。
ここで、残りのカメラパラメータとは、推定したい全カメラパラメータから、ステップＳ１２で最適解候補が算出された部分問題に依存する、外部パラメータの一部を除いた、カメラパラメータのことである。元の分割される前の誤差関数においては、すべてのカメラパラメータが未知数であるが、ステップＳ１３における部分問題においては、ステップＳ１２で算出された最適解候補（カメラパラメータ）が既知の定数として扱われる。すなわち、ステップＳ１３では、ステップＳ１２で算出されなかった残りのカメラパラメータが、ステップＳ１２では対象とならなかった部分問題の停留点として算出される。
また、ステップＳ１２で出力され最適解候補が複数の場合、残存パラメータ算出部２３は、複数の停留点それぞれ毎に、その他の部分問題の停留点の算出を行う。すなわち、例えば、ステップＳ１２でＫ個の最適解候補が算出された場合は、ステップＳ１３は、その他の部分問題をＫ回解くことになる。
最後に、パラメータ最適化部３０は、ステップＳ１２とステップＳ１３とで出力されたすべての最適解候補と対応点とを入力として、元の誤差関数の最適化を行い、誤差が最小となる解（最適解）を求め、求めた解をカメラパラメータとして出力する（ステップＳ１４）。
ステップＳ１３では、パラメータ最適化部３０は、更に、入力または出力されたカメラパラメータが物理的または理論的に意味のあるかどうかの判定を行うこともできる。物理的な意味とは、例えば、すべての３次元点はカメラの前面に存在しなくてはならないことをいう。従って、パラメータ最適化部３０は、一部またはすべての３次元点がカメラの裏側に存在するような不適切なカメラパラメータを除外したり、カメラパラメータの正負符号を逆にして3次元点がカメラの前面に存在するように補正したりしてもよい。また、一部のカメラパラメータが取りうる値の範囲が予めわかっている場合、パラメータ最適化部３０は、それを閾値として不適切な入出力を除外してもよい。例えば、ズームレンズの仕様がわかっていれば、推定される焦点距離の上限と下限とを事前に決めることができるので、パラメータ最適化部３０は、この上限及び下限から外れる入出力を除外する。
理論的な意味とは、例えば、計算結果が誤差関数の２次最適性を満たさないことをいう。パラメータ最適化部３０は、計算されたカメラパラメータが誤差関数の２次最適性を満たさない場合は、このカメラパラメータを除外する。
このように不適切なカメラパラメータが除外された場合、パラメータ最適化部３０は、物理的または理論的に意味のある値の範囲内での、誤差最小となる解を出力する。
もし上記の条件を満たす解が複数存在し、それらの誤差が同一とみなせる場合には、パラメータ最適化部３０は、すべての解を出力してもよい。この場合、大域的最適解が複数存在するとみなせるため、出力されるカメラパラメータも複数である。
［具体例］
ここで、本実施の形態の具体例について以下に説明する。なお、以下の具体例では、外部パラメータは、それぞれ３自由度の位置と姿勢とする。また、ズームレンズにおいてズームを変更した際に影響を受けるのは焦点距離とレンズ歪み係数であるため、内部パラメータは１自由度の焦点距離と３自由度のレンズ歪み係数とする。
レンズ歪みの種類は樽型歪（radialdistortion）とし、計算モデルはAndrew Fitzgibbonが提案したdivision modelとする。それ以外の内部パラメータは、光学中心については画像の中心（つまり画像の大きさの半分）と仮定し、アスペクト比については１と仮定し、せん断係数については０と仮定する。これらの仮定は、上述した非特許文献４、５、６、７でも採用されている。昨今のデジタルカメラは、高精度に製造されているため、上記の仮定は妥当な仮定である。
但し、どのカメラパラメータを未知または既知とするかはこの限りではなく、利用シーンに応じて当業者が自由に変更することができる。また、レンズ歪みの種類および計算モデルについても同様である。
また、以下の説明では、上添字の「Ｔ」は行列およびベクトルの転置を表し、「０」はゼロ行列およびゼロベクトルを表し、「Ｉ」は単位行列を表し、「|| ||」はベクトルのＬ２ノルムを表す。外部パラメータにおいては、カメラ姿勢を３×３の回転行列Ｒで表し、カメラ位置を３×１のベクトルｔで表す。内部パラメータにおいては、焦点距離を「ｆ」で表し、焦点距離の逆数１／ｆを「ｗ」で表記し、レンズ歪み係数を３×１のベクトル「ｋ」で表す。また、ベクトルｋを下記の数１で表す。
Figure 0007052788000001
対応点は、入力点数をｎ組（ｎ≧５）とし、i番目の３次元座標を「Ｘｉ」、division modelにより斉次化した２次元座標を「ｍｉ」と表記する。２次元座標ｍｉは下記の数２で表される。
Figure 0007052788000002
また、誤差関数の定義を行う。誤差関数は、Ｘｉとのｍｉとの変換関係を表しており、以下の数３で定義される。
Figure 0007052788000003
上記数３において、「～」は定数倍に不定なことを表している。「Ｋ」は、［１，１，ｗ］を対角成分とする３×３の行列である。
また、上記数３より、ｎ組の対応点に対する誤差関数は以下の数４で定義される。
Figure 0007052788000004
ここで、［ｍｉ］ｘはベクトルの外積を示す３×３の歪対称行列であり、具体的な成分は下記の数５で表される。
Figure 0007052788000005
上記数４は、代数学的誤差と呼ばれ、上記数３の両辺のベクトル方向が互いに近しくなることを要求している。
次に、部分問題の定義を行う。［ｍｉ］ｘの行ベクトルを下記の数６～数８に示す通りとすると上記数４は以下の数９と等価である。なお、制約条件の記述は省略している。
Figure 0007052788000006
Figure 0007052788000007
Figure 0007052788000008
Figure 0007052788000009
ここで、上記数９の各項を独立の最適化問題とみなすと、以下の数１０～数１２（制約条件の記述は省略）で表される３つの部分問題が得られる。
Figure 0007052788000010
Figure 0007052788000011
Figure 0007052788000012
次に、各部分問題の具体的な解法について説明する。まず、上記数１２に着目すると、下記の数１３が成立するため、上記数１２を下記の数１４で表すことができる。
Figure 0007052788000013
Figure 0007052788000014
上記数１４において、ｒｊ ＴはＲの第j行、ｔｊはｔの第j成分である。また、上記数１４は、内部パラメータを含まず、ｔ３を除く外部パラメータのみで表現されている。また、Ｒに関する制約条件を考慮すると、｛ｒ１，ｒ２，ｔ１，ｔ２｝の自由度は５であるため、対応点数nが６組以上であれば、上記数１２を解くことができる。また、上記数１２は上記数１４を用いて、下記の数１５で表すことができる。
Figure 0007052788000015
ここで、Ａ，Ｂ，Ｍは、ｕｉ、ｖｉ、Ｘｉからなるｎ×６、ｎ×２、６×６の行列である。これらについての詳しい導出方法は、上記非特許文献７に開示されているため、ここでは記述を省略する。
上記数１５を線形最小二乗法により解く場合、すべての３次元点が１つの平面上に分布しているか否か、および点の数で解法が異なる。なぜならば、入力にノイズを含む一般的な場合、行列Ｍの階数は、３次元点が１つの平面上に分布しており、かつ、６点以上の場合は４となり，２つ以上の平面に分布しており、かつ、６点の場合は５、７点以上の場合６となるからである。例えば、ｚ＝０の平面上に分布する場合を考えると、上記数１４より階数が４であることが証明できる。任意の平面は、剛体変換によりｚ＝０の平面へ変換できるため、すべての平面についてＭの階数は４である。
すべての３次元点が１つの平面上に分布しているか否かは、例えば、下記の数１６を用いて判定できる。
Figure 0007052788000016
ここで、Ｙはすべての３次元点の重心であり、ｄｅｔ（）は行列式を計算する関数である。すべての３次元点が1つの平面上に分布している場合、３点以上の３次元点は互いに従属であるため、上記数１６は、理論的に０となる。実際には、計算精度及び３次元点には、微小のノイズが含まれるため、上記数１６は完全な０にはならない。そのため、事前に決定したしきい値以下であれば平面と判定する。
次に、上記数１５を線形最小二乗法により解く具体的な方法について説明する。まず、３次元点が２つ以上の平面に分布している場合、つまり上記数１６で非平面と判定された場合について説明する。このとき、変数ベクトル［ｒ１ Ｔ，ｒ２ Ｔ］について、回転行列の拘束条件を無視して、０ベクトルでないという非自明な拘束条件を導入すると、その自由度は５である。
また、対応点の数が７点以上の場合、行列Ｍの階数は６であるから、上記数１５を最小化する変数ベクトル［ｒ１ Ｔ，ｒ２ Ｔ］は、一般的な線形最小二乗法として定式化でき、Ｍの最小固有値に対応する固有ベクトルとして算出できる。また、非平面かつ対応点の数が６点の場合、行列Ｍの階数は５であり、上記数１５を最小化する変数ベクトル［ｒ１ Ｔ，ｒ２ Ｔ］は、行列Ｍのゼロ空間を表す２つの直交基底ベクトルの線形和で表される。例えば、２つの直交基底ベクトルをｖ１、ｖ２、係数をaとすると、［ｒ１ Ｔ，ｒ２ Ｔ］＝ａ＊ｖ１＋ｖ２である。ｒ１とｒ２のノルムが等しく直交するという制約条件を用いればaを算出でき、変数ベクトル［ｒ１ Ｔ，ｒ２ Ｔ］を決定できる。そして、ｒ３は、ｒ１とｒ２との外積として、ｒ１とｒ２との推定後に算出すればよい。こうして得られた回転行列は、上記数４で表された拘束条件を一般に満たさない。そこで、特異値分解を利用して拘束条件を満たす回転行列へ補正する手法が広く知られている。
次に、３次元点が一つの平面上に分布している場合の上記数１５を、線形最小二乗法によって解く具体的な解法について説明する。ここでは簡単のため、すべての３次元点はｚ＝０の平面上に分布しているものとする。任意の平面は剛体変換によりｚ＝０の平面へ変換できるため、上記数１６で行列式が０と判定された場合に、ｚ＝０の平面へと変換させてもよい。上述したように平面の場合、行列Ｍの階数は４であり、さらに、ｚ＝０の平面である場合、Ｍの第３列、第６列、第３行、第６行の成分は、すべて０となる。行列Ｍの非ゼロ成分を抽出した４×４の行列をＮ、ｒ１とｒ２の第１成分及び第２成分を、それぞれｐ１、ｐ２とすると、上記数１５の第１項は下記数１７で表される。また、下記の数１７に含まれる、変数ベクトルを下記の数１８に示すように定義する。
Figure 0007052788000017
Figure 0007052788000018
上記数１７において、変数ベクトルＲは４次元であり、Ｍの階数は４であるから、上記数１７を最小化する変数ベクトルＲは、３次元点が非平面上に分布する場合と同様に、線形最小二乗法により算出することができる。推定されていないｒ１とｒ２それぞれの第３成分を算出するには、回転行列の直交性を利用すればよい。例えば、以下のようにしてもよい。回転行列の第i行j列成分をrijと書くと、第３行３列成分は、一般にｒ３３＝ｒ１１ｒ２２－ｒ１２ｒ２１と計算できる。ここで、残る未知数ｒ１３、ｒ２３、ｒ３１、ｒ３２は、下記数１９で表される直交性の拘束条件を満たす。
Figure 0007052788000019
上記数１９の第１式と第２式との第１行３列、第２行３列の成分は、ｒ１３、ｒ２３、ｒ３１、ｒ３２についての線形１次方程式であるから、それを解けばスケール不定性を持ったｒ１３、ｒ２３、ｒ３１、ｒ３２が算出でき、第１式と第２式との対角成分を用いて、全rijについてのスケールを線形最小二乗法により算出でき、スケール不定性を除いた回転行列が推定できる。
例えば、下記数２０で表される連立方程式を解いてもよい。
Figure 0007052788000020
上記数２０は、ｒ１３またはｒ２３について、４次式となるので、最大４個の解が得られる。ｒ１とｒ２のすべての成分を求めたあと、ノルムが１になるように正規化し、ｒ１とｒ２との外積としてｒ3を計算すれば、回転行列の成分をすべて算出することができる。上記数２０は、行方向の直交性を元にしているが、列方向の直交性を用いても同様に計算できる。
そこで、３次元点が非平面に分布する場合と同様に、特異値分解を利用して拘束条件を満たした回転行列へ補正することができる。
以上より、３次元点が、平面及び非平面のいずれかに分布する場合において、上記数１５のｒ１，ｒ２を算出する方法を説明した。上記数１５のｔ１，ｔ２は、算出したｒ１，ｒ２を上記数１５の第２項に代入することで算出できる。
残るパラメータｔ３，ｗ，ｋ＝［ｋ１，ｋ２，ｋ３］を線形最小二乗法により算出する方法は、上記非特許文献７に詳しく記載されているため、ここでは記述を省略する。
このように、すべてのカメラパラメータが推定できたが、これらは上記数１０～数１２で表される部分問題の解である。そこで、これらを初期値としてニュートン法などの非線形最適化法を実行して、上記数４または数９を解けば、最終的に最適なカメラパラメータが推定できる。
ここで、上述の具体例を図３に示したステップＳ１１～Ｓ１４に沿って以下に説明する。
［ステップＳ１１］
まず、ステップＳ１１において、分割処理部２１は、上記数９で表される誤差関数を分割し、上記数１０、数１１、及び数１２の部分問題とする。
［ステップＳ１２］
次に、ステップＳ１２において、部分パラメータ算出部２２は、上記数１２で表される部分問題を、線形最小二乗問題として解いて、その解を出力する。部分パラメータ算出部２２は、まず上記数１６に基づいて３次元点が平面か非平面かの判定を行い、判定結果に基づいて平面か非平面かいずれかの解を計算しても良いし、判定を行わない、又は判定結果に関わらず、平面と非平面の両方について複数の解を計算してもよい。また、部分パラメータ算出部２２は、解が複数の場合に、すべての解を出力しても良いし、上記数１２を最小化する唯一の解を出力してもよい。
［ステップＳ１３］
次に、ステップＳ１３において、残存パラメータ算出部２３は、部分パラメータ算出部２２により計算された外部パラメータの一部を入力として、残存パラメータ、つまり、残りの外部パラメータと内部パラメータとを線形最小二乗法により解いて、出力する。ステップＳ１２で算出された解が複数であり、複数の解が入力された場合は、残存パラメータ算出部２３は、入力された解毎に、残存パラメータを算出する。このとき、部分パラメータ算出部２２は、すべてのパラメータを出力してもよいし、上記数４および部分問題のいずれかを最小化する唯一の解を出力してもよい。
［ステップＳ１４］
最後に、ステップＳ１４において、パラメータ最適化部３０は、部分パラメータ算出部２２と残存パラメータ算出部２３とで出力された１組以上のカメラパラメータを初期値として、上記数４の非線形最適化により誤差が最小となる最適解を出力する。ステップＳ１４の実行後、カメラパラメータ推定装置１０は動作を停止する。
また、ステップＳ１４において、パラメータ最適化部３０は、最適化前に入力パラメータを用いて、上記数４の誤差を計算し、誤差が閾値以下であれば最適化をせずに、入力パラメータをそのまま出力してもよい。また、ステップＳ１４においては、誤差関数として、上記数４以外にも、物理的な距離を表すいわゆる再投影誤差関数が採用されていてもよい。この場合、上記数４を最適化する場合よりも計算時間はかかるが、より高精度な解が得られる可能性がある。
［実施の形態における効果］
以上のとおり、本実施形態によれば、レンズ歪み係数を含む内部パラメータと外部パラメータとを、数値的に安定に、且つ、高精度に推定可能である。その理由は、以下のとおりである。
上記数４および数９は、２つのベクトルｍｉとＫ（ＲＸｉ＋ｔ）との方向が近しくなることを表している。換言すれば、上記数４および数９の最適解が求まったとすると、［ｍｉ］ｘＫ（ＲＸｉ＋ｔ）で表されるベクトルの各成分はそれぞれ０に近くなると期待される。そこで、上記数９を直接に解く代わりに、上記数９の各項を独立の制約付き最適化問題（上記数１０～数１２の部分問題）とみなして解いたとしても、それぞれの最適解は原問題である上記数９の最適解に非常に近いと考えられる。
更に、通常であれば計算量は未知数の数に依存して増加するが、本実施の形態では、外部パラメータのみに依存する部分問題（上記数１２）から順に解くことで、従来方法と同等以下の低演算量で実行が可能である。
また、従来方法は、上記数４または数１２を、拘束条件を含む多項式問題として定式化し、固有値分解に基づく数値的に不安定な方法を用いて解を算出しているが、本実施の形態では一部の拘束条件を無視した線形最小二乗問題として定式化し、数値的に安定な線形最小二乗法により解を求めた後で、拘束条件を満たすように補正する。以上の理由により、本実施の形態によれば、内部パラメータと外部パラメータとが、数値的に安定に、且つ、高精度に推定される。
［変形例］
本実施の形態は、上述した例に限定されるものではない。本実施の形態は、上述した例に対して、いわゆる当業者が理解し得る多様な変更を適用することが可能である。例えば、本実施の形態は、以下の変形例に示す形態によっても実施可能である。
変形例１
本実施の形態において、未知とする内部パラメータは、特に限定されるものではない。例えば、非特許文献４に開示されているように、推定すべき内部パラメータが焦点距離のみの場合、上記数１７においてｋ＝０として解けばよい。また、非特許文献３に開示されているように、すべての内部パラメータが既知である場合、上記数１７において、ｋ＝０、ｗ＝１／ｆとして解けばよい。変形例１においても、未知となる内部パラメータをどのように設定しようとも、上記数１２を用いて外部パラメータの一部を算出してから、残存パラメータを推定するという方法は変更する必要はない。
変形例２
本実施の形態では、カメラパラメータ推定装置１０は、一部の外部パラメータを外部より受け付けるために、外部パラメータ入力部を備えていてもよい。外部パラメータ入力部は、例えば、外部の加速度センサから、重力方向における加速度を特定するセンサ値を受け付けたり、ジャイロセンサから、角速度を特定するセンサ値を受け付けたりしてもよい。後者の場合、外部パラメータ入力部は、回転行列を算出し、得られた回転行列をパラメータ算出部２０に入力する。変形例２においては、推定すべき外部パラメータ数は減少する。但し、変形例２においても、上記数１２で未知の外部パラメータを算出してから、残りの内部パラメータを推定するという方法は変更する必要はない。
変形例３
本実施の形態において、外部パラメータはカメラの位置姿勢に限定されない。例えば、静止したカメラの前で動く対象物を撮影した場合、外部パラメータとして、対象物の位置姿勢を表すパラメータが推定される。カメラが動くか、対象物が動くかは相対的に等価なため、解くべき問題は変わらない。
変形例４
本実施の形態において、撮影に用いられるカメラは単眼カメラに限定されない。上記数３では、Ｘｉ とｍ ｉとの変換関係を、いわゆる単眼のピンホールカメラモデルを用いて説明しているが、本実施の形態では、複数の焦点を持つ一般カメラモデルが用いられていてもよい。焦点位置が事前に既知であれば、上述した例と同様の方法でカメラパラメータを推定可能である。
［プログラム］
本実施の形態におけるプログラムは、コンピュータに、図３に示すステップＳ１１～Ｓ１４を実行させるプログラムであれば良い。このプログラムをコンピュータにインストールし、実行することによって、本実施の形態におけるカメラパラメータ推定装置１０とカメラパラメータ推定方法とを実現することができる。この場合、コンピュータのプロセッサは、パラメータ算出部２０及びパラメータ最適化部３０として機能し、処理を行なう。
また、本実施の形態におけるプログラムは、複数のコンピュータによって構築されたコンピュータシステムによって実行されても良い。この場合は、例えば、各コンピュータが、それぞれ、パラメータ算出部２０、及びパラメータ最適化部３０のいずれかとして機能しても良い。
［物理構成］
ここで、本実施の形態におけるプログラムを実行することによって、カメラパラメータ推定装置１０を実現するコンピュータについて図４を用いて説明する。図４は、本発明の実施の形態におけるカメラパラメータ推定装置１０を実現するコンピュータの一例を示すブロック図である。
図４に示すように、コンピュータ１１０は、ＣＰＵ（CentralProcessing Unit）１１１と、メインメモリ１１２と、記憶装置１１３と、入力インターフェイス１１４と、表示コントローラ１１５と、データリーダ／ライタ１１６と、通信インターフェイス１１７とを備える。これらの各部は、バス１２１を介して、互いにデータ通信可能に接続される。なお、コンピュータ１１０は、ＣＰＵ１１１に加えて、又はＣＰＵ１１１に代えて、ＧＰＵ（GraphicsProcessing Unit）、又はＦＰＧＡ（Field-Programmable Gate Array）を備えていても良い。
ＣＰＵ１１１は、記憶装置１１３に格納された、本実施の形態におけるプログラム（コード）をメインメモリ１１２に展開し、これらを所定順序で実行することにより、各種の演算を実施する。メインメモリ１１２は、典型的には、ＤＲＡＭ（Dynamic Random Access Memory）等の揮発性の記憶装置である。また、本実施の形態におけるプログラムは、コンピュータ読み取り可能な記録媒体１２０に格納された状態で提供される。なお、本実施の形態におけるプログラムは、通信インターフェイス１１７を介して接続されたインターネット上で流通するものであっても良い。
また、記憶装置１１３の具体例としては、ハードディスクドライブの他、フラッシュメモリ等の半導体記憶装置が挙げられる。入力インターフェイス１１４は、ＣＰＵ１１１と、キーボード及びマウスといった入力機器１１８との間のデータ伝送を仲介する。表示コントローラ１１５は、ディスプレイ装置１１９と接続され、ディスプレイ装置１１９での表示を制御する。
データリーダ／ライタ１１６は、ＣＰＵ１１１と記録媒体１２０との間のデータ伝送を仲介し、記録媒体１２０からのプログラムの読み出し、及びコンピュータ１１０における処理結果の記録媒体１２０への書き込みを実行する。通信インターフェイス１１７は、ＣＰＵ１１１と、他のコンピュータとの間のデータ伝送を仲介する。
また、記録媒体１２０の具体例としては、ＣＦ（CompactFlash（登録商標））及びＳＤ（Secure Digital）等の汎用的な半導体記憶デバイス、フレキシブルディスク（Flexible Disk）等の磁気記憶媒体、又はＣＤ－ＲＯＭ（Compact DiskRead Only Memory）などの光学記憶媒体が挙げられる。
なお、本実施の形態におけるカメラパラメータ推定装置１０は、プログラムがインストールされたコンピュータではなく、各部に対応したハードウェアを用いることによっても実現可能である。更に、カメラパラメータ推定装置１０は、一部がプログラムで実現され、残りの部分がハードウェアで実現されていてもよい。また、カメラパラメータ推定装置１０は、各部が、それぞれ別々のユニットとして実現されていてもよい。
上述した実施の形態の一部又は全部は、以下に記載する（付記１）～（付記１２）によって表現することができるが、以下の記載に限定されるものではない。
（付記１）
カメラによる対象物の撮影時におけるカメラパラメータを推定するための装置であって、
前記対象物に関する６組以上の３次元点と、前記３次元点それぞれに対応する画像上の２次元点との変換を表す誤差関数を、前記カメラパラメータの依存関係と前記３次元点の分布の平面性とに基づいて、複数の部分問題に分割し、前記複数の部分問題それぞれにおいて最適解候補を算出する、パラメータ算出部と、
前記最適解候補を初期値として、前記誤差関数によって得られる誤差を最小とする最適解を求め、求めた前記最適解を最適なカメラパラメータとして出力する、パラメータ最適化部と、
を備えている、ことを特徴とするカメラパラメータ推定装置。
（付記２）
前記カメラパラメータが、前記カメラの外部に関する外部パラメータを含み、
前記パラメータ算出部が、
前記複数の部分問題のうち、前記外部パラメータの一部のみに依存する１つの部分問題について、前記最適解候補を算出し、
更に、算出した前記最適解候補と、前記３次元点及び前記２次元点とに基づいて、前記最適解候補が未だ算出されていない前記部分問題について、前記最適解候補を算出する、
付記１に記載のカメラパラメータ推定装置。
（付記３）
前記パラメータ最適化部が、前記パラメータ算出部によって算出された前記最適解候補と、前記３次元点及び前記２次元点とに基づいて、前記誤差関数を最適化することによって、前記誤差を最小とする最適解を求める、
付記１または２に記載のカメラパラメータ推定装置。
（付記４）
前記外部パラメータが、前記カメラの３次元区間における位置及び姿勢のいずれか一方または両方を含み、
前記カメラパラメータが、前記外部パラメータに加えて、前記カメラの内部に関する内部パラメータを含み、前記内部パラメータは、前記カメラにおける、レンズの焦点距離、光学中心、アスペクト比、せん断係数、及びレンズ歪み係数のうち、少なくとも１つを含む、
付記２に記載のカメラパラメータ推定装置。
（付記５）
カメラによる対象物の撮影時におけるカメラパラメータを推定するための方法であって、
（ａ）前記対象物に関する６組以上の３次元点と、前記３次元点それぞれに対応する画像上の２次元点との変換を表す誤差関数を、前記カメラパラメータの依存関係と前記３次元点の分布の平面性とに基づいて、複数の部分問題に分割し、前記複数の部分問題それぞれにおいて最適解候補を算出する、ステップと、
（ｂ）前記最適解候補を初期値として、前記誤差関数によって得られる誤差を最小とする最適解を求め、求めた前記最適解を最適なカメラパラメータとして出力する、ステップと、
を有する、ことを特徴とするカメラパラメータ推定方法。
（付記６）
前記カメラパラメータが、前記カメラの外部に関する外部パラメータを含み、
前記（ａ）のステップにおいて、
前記複数の部分問題のうち、前記外部パラメータの一部のみに依存する１つの部分問題について、前記最適解候補を算出し、
更に、算出した前記最適解候補と、前記３次元点及び前記２次元点とに基づいて、前記最適解候補が未だ算出されていない前記部分問題について、前記最適解候補を算出する、
付記５に記載のカメラパラメータ推定方法。
（付記７）
前記（ｂ）のステップにおいて、前記（ａ）のステップによって算出された前記最適解候補と、前記３次元点及び前記２次元点とに基づいて、前記誤差関数を最適化することによって、前記誤差を最小とする最適解を求める、
付記５または６に記載のカメラパラメータ推定方法。
（付記８）
前記外部パラメータが、前記カメラの３次元区間における位置及び姿勢のいずれか一方または両方を含み、
前記カメラパラメータが、前記外部パラメータに加えて、前記カメラの内部に関する内部パラメータを含み、前記内部パラメータは、前記カメラにおける、レンズの焦点距離、光学中心、アスペクト比、せん断係数、及びレンズ歪み係数のうち、少なくとも１つを含む、
付記６に記載のカメラパラメータ推定方法。
（付記９）
コンピュータによって、カメラによる対象物の撮影時におけるカメラパラメータを推定するためのプログラムであって、
前記コンピュータに、
（ａ）前記対象物に関する６組以上の３次元点と、前記３次元点それぞれに対応する画像上の２次元点との変換を表す誤差関数を、前記カメラパラメータの依存関係と前記３次元点の分布の平面性とに基づいて、複数の部分問題に分割し、前記複数の部分問題それぞれにおいて最適解候補を算出する、ステップと、
（ｂ）前記最適解候補を初期値として、前記誤差関数によって得られる誤差を最小とする最適解を求め、求めた前記最適解を最適なカメラパラメータとして出力する、ステップと、
を実行させる、プログラム。
（付記１０）
前記カメラパラメータが、前記カメラの外部に関する外部パラメータを含み、
前記（ａ）のステップにおいて、
前記複数の部分問題のうち、前記外部パラメータの一部のみに依存する１つの部分問題について、前記最適解候補を算出し、
更に、算出した前記最適解候補と、前記３次元点及び前記２次元点とに基づいて、前記最適解候補が未だ算出されていない前記部分問題について、前記最適解候補を算出する、
付記９に記載のプログラム。
（付記１１）
前記（ｂ）のステップにおいて、前記（ａ）のステップによって算出された前記最適解候補と、前記３次元点及び前記２次元点とに基づいて、前記誤差関数を最適化することによって、前記誤差を最小とする最適解を求める、
付記９または１０に記載のプログラム。
（付記１２）
前記外部パラメータが、前記カメラの３次元区間における位置及び姿勢のいずれか一方または両方を含み、
前記カメラパラメータが、前記外部パラメータに加えて、前記カメラの内部に関する内部パラメータを含み、前記内部パラメータは、前記カメラにおける、レンズの焦点距離、光学中心、アスペクト比、せん断係数、及びレンズ歪み係数のうち、少なくとも１つを含む、
付記１０に記載のプログラム。
以上、実施の形態を参照して本願発明を説明したが、本願発明は上記実施の形態に限定されるものではない。本願発明の構成や詳細には、本願発明のスコープ内で当業者が理解し得る様々な変更をすることができる。
この出願は、２０１７年３月１４日に出願された日本出願特願２０１７－０４８９１５を基礎とする優先権を主張し、その開示の全てをここに取り込む。
以上のように、本発明によれば、レンズ歪み係数を含む内部パラメータと外部パラメータとを高精度に推定することができる。本発明は、カメラパラメータの推定が求められる分野、例えば、ロボットの制御、画像処理等に有効である。
１０ カメラパラメータ推定装置
２０ パラメータ算出部
２１ 分割処理部
２２ 部分パラメータ算出部
２３ 残存パラメータ算出部
３０ パラメータ最適化部
１１０ コンピュータ
１１１ ＣＰＵ
１１２ メインメモリ
１１３ 記憶装置
１１４ 入力インターフェイス
１１５ 表示コントローラ
１１６ データリーダ／ライタ
１１７ 通信インターフェイス
１１８ 入力機器
１１９ ディスプレイ装置
１２０ 記録媒体
１２１ バス

Claims (6)

  1. カメラによる対象物の撮影時におけるカメラパラメータを推定するための装置であって、
    前記対象物に関する６組以上の３次元点と、前記３次元点それぞれに対応する画像上の２次元点との変換を表す誤差関数を、前記カメラパラメータの依存関係と前記３次元点の分布の平面性とに基づいて、複数の部分問題に分割し、前記複数の部分問題それぞれにおいて最適解候補を算出する、パラメータ算出部と、
    前記最適解候補を初期値として、前記誤差関数によって得られる誤差を最小とする最適解を求め、求めた前記最適解を最適なカメラパラメータとして出力する、パラメータ最適化部と、
    を備えている、ことを特徴とするカメラパラメータ推定装置。
  2. 前記カメラパラメータが、前記カメラの外部に関する外部パラメータを含み、
    前記パラメータ算出部が、
    前記複数の部分問題のうち、前記外部パラメータの一部のみに依存する１つの部分問題について、前記最適解候補を算出し、
    更に、算出した前記最適解候補と、前記３次元点及び前記２次元点とに基づいて、前記最適解候補が未だ算出されていない前記部分問題について、前記最適解候補を算出する、
    請求項１に記載のカメラパラメータ推定装置。
  3. 前記パラメータ最適化部が、前記パラメータ算出部によって算出された前記最適解候補と、前記３次元点及び前記２次元点とに基づいて、前記誤差関数を最適化することによって、前記誤差を最小とする最適解を求める、
    請求項１または２に記載のカメラパラメータ推定装置。
  4. 前記外部パラメータが、前記カメラの３次元区間における位置及び姿勢のいずれか一方または両方を含み、
    前記カメラパラメータが、前記外部パラメータに加えて、前記カメラの内部に関する内部パラメータを含み、前記内部パラメータは、前記カメラにおける、レンズの焦点距離、光学中心、アスペクト比、せん断係数、及びレンズ歪み係数のうち、少なくとも１つを含む、
    請求項２に記載のカメラパラメータ推定装置。
  5. カメラによる対象物の撮影時におけるカメラパラメータを推定するための方法であって、
    （ａ）前記対象物に関する６組以上の３次元点と、前記３次元点それぞれに対応する画像上の２次元点との変換を表す誤差関数を、前記カメラパラメータの依存関係と前記３次元点の分布の平面性とに基づいて、複数の部分問題に分割し、前記複数の部分問題それぞれにおいて最適解候補を算出する、ステップと、
    （ｂ）前記最適解候補を初期値として、前記誤差関数によって得られる誤差を最小とする最適解を求め、求めた前記最適解を最適なカメラパラメータとして出力する、ステップと、
    を有する、ことを特徴とするカメラパラメータ推定方法。
  6. コンピュータによって、カメラによる対象物の撮影時におけるカメラパラメータを推定するためのプログラムであって、
    前記コンピュータに、
    （ａ）前記対象物に関する６組以上の３次元点と、前記３次元点それぞれに対応する画像上の２次元点との変換を表す誤差関数を、前記カメラパラメータの依存関係と前記３次元点の分布の平面性とに基づいて、複数の部分問題に分割し、前記複数の部分問題それぞれにおいて最適解候補を算出する、ステップと、
    （ｂ）前記最適解候補を初期値として、前記誤差関数によって得られる誤差を最小とする最適解を求め、求めた前記最適解を最適なカメラパラメータとして出力する、ステップと、
    を実行させる、プログラム。
JP2019505764A 2017-03-14 2018-02-05 カメラパラメータ推定装置、カメラパラメータ推定方法、及びプログラム Active JP7052788B2 (ja)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
JP2017048915 2017-03-14
JP2017048915 2017-03-14
PCT/JP2018/003814 WO2018168255A1 (ja) 2017-03-14 2018-02-05 カメラパラメータ推定装置、カメラパラメータ推定方法、及びコンピュータ読み取り可能な記録媒体

Publications (2)

Publication Number Publication Date
JPWO2018168255A1 JPWO2018168255A1 (ja) 2020-01-23
JP7052788B2 true JP7052788B2 (ja) 2022-04-12

Family

ID=63523908

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2019505764A Active JP7052788B2 (ja) 2017-03-14 2018-02-05 カメラパラメータ推定装置、カメラパラメータ推定方法、及びプログラム

Country Status (3)

Country Link
US (1) US10719956B2 (ja)
JP (1) JP7052788B2 (ja)
WO (1) WO2018168255A1 (ja)

Families Citing this family (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US11727597B2 (en) * 2018-12-21 2023-08-15 Sony Group Corporation Calibrating volumetric rig with structured light
JP7111183B2 (ja) * 2019-01-09 2022-08-02 日本電気株式会社 カメラパラメータ推定装置、カメラパラメータ推定方法、及びプログラム
JP7151879B2 (ja) 2019-04-08 2022-10-12 日本電気株式会社 カメラ校正装置、カメラ校正方法、及びプログラム
WO2020246261A1 (ja) * 2019-06-05 2020-12-10 ソニー株式会社 移動体、位置推定方法、およびプログラム
CN113298868B (zh) * 2021-03-17 2024-04-05 阿里巴巴创新公司 模型建立方法、装置、电子设备、介质及程序产品
CN113516719B (zh) * 2021-04-21 2023-08-15 深圳臻像科技有限公司 一种基于多单应性矩阵的相机标定方法、系统及存储介质
CN113256540B (zh) * 2021-07-14 2021-11-19 智道网联科技(北京)有限公司 图像去畸变方法和装置、电子设备和计算机可读存储介质

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2010250452A (ja) 2009-04-14 2010-11-04 Tokyo Univ Of Science 任意視点画像合成装置
JP2015106287A (ja) 2013-11-29 2015-06-08 キヤノン株式会社 キャリブレーション装置及び方法

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2006067272A (ja) * 2004-08-27 2006-03-09 Matsushita Electric Ind Co Ltd カメラキャリブレーション装置及びカメラキャリブレーション方法
WO2010109730A1 (ja) * 2009-03-26 2010-09-30 アイシン精機株式会社 カメラ校正装置
JP2017163386A (ja) * 2016-03-10 2017-09-14 日本電気株式会社 カメラパラメータ推定装置、カメラパラメータ推定方法、及びプログラム

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2010250452A (ja) 2009-04-14 2010-11-04 Tokyo Univ Of Science 任意視点画像合成装置
JP2015106287A (ja) 2013-11-29 2015-06-08 キヤノン株式会社 キャリブレーション装置及び方法

Also Published As

Publication number Publication date
JPWO2018168255A1 (ja) 2020-01-23
US20200051279A1 (en) 2020-02-13
US10719956B2 (en) 2020-07-21
WO2018168255A1 (ja) 2018-09-20

Similar Documents

Publication Publication Date Title
JP7052788B2 (ja) カメラパラメータ推定装置、カメラパラメータ推定方法、及びプログラム
Usenko et al. The double sphere camera model
Geneva et al. Openvins: A research platform for visual-inertial estimation
Urban et al. Mlpnp-a real-time maximum likelihood solution to the perspective-n-point problem
JP5839971B2 (ja) 情報処理装置、情報処理方法及びプログラム
WO2013111229A1 (ja) カメラ校正装置、カメラ校正方法およびカメラ校正用プログラム
JP6547754B2 (ja) 三角測量装置、三角測量方法およびそのプログラム
US20140316698A1 (en) Observability-constrained vision-aided inertial navigation
JP7099627B2 (ja) カメラパラメータ推定装置、カメラパラメータ推定方法、及びプログラム
Albl et al. Degeneracies in rolling shutter sfm
US9865061B2 (en) Constructing a 3D structure
Kukelova et al. Hand-eye calibration without hand orientation measurement using minimal solution
Barath et al. Relative pose from sift features
JP7111183B2 (ja) カメラパラメータ推定装置、カメラパラメータ推定方法、及びプログラム
CN114049401A (zh) 双目相机标定方法、装置、设备及介质
JP2017163386A (ja) カメラパラメータ推定装置、カメラパラメータ推定方法、及びプログラム
Yu et al. An efficient and globally optimal method for camera pose estimation using line features
JP4970118B2 (ja) カメラ校正方法、そのプログラム、記録媒体、装置
JP2020041950A (ja) 測量装置、測量方法、及びプログラム
JP6260533B2 (ja) 位置姿勢推定装置、位置姿勢推定方法および位置姿勢推定プログラム
Micusik et al. Renormalization for initialization of rolling shutter visual-inertial odometry
JP2005063012A (ja) 全方位カメラ運動と３次元情報の復元方法とその装置及びプログラム並びにこれを記録した記録媒体
JP5215615B2 (ja) ３次元位置情報復元装置およびその方法
JP7546748B2 (ja) カメラ校正方法及びカメラ校正システム
WO2023166618A1 (ja) カメラパラメータ推定装置、カメラパラメータ推定方法、及びコンピュータ読み取り可能な記録媒体

Legal Events

Date Code Title Description
A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20190911

A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20210106

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20220301

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20220314

R151 Written notification of patent or utility model registration

Ref document number: 7052788

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R151