(57)【特許請求の範囲】
【請求項１】
（ｉ）外部からの光の入射を避けるために不透明な材料で作られたボックス（１）と、レンズスコープへのアクセスを可能にするための１つまたは複数の開口部（２）であって、各開口部は、外部からの光の入射を最小化するゴム膜（３）によって囲まれており、前記ボックスは、透明または半透明の材料からなるセパレータ（６）によって上部隔室と下部隔室の２つの隔室に分割され、前記上部隔室は、目視検査を実施するための平均値をより良く再現するための流体で満たされることが可能であり、
（ｉｉ）前記上部隔室内であって前記セパレータの上に配置された透明または半透明の材料からなる既知のパターンを有する平面グリッド（７）と、
（ｉｉｉ）前記ボックスの前記下部隔室（５）に取り付けられた、前記平面グリッドの背面照明拡散のための光源と、
（ｉｖ）逐次処理機能および並列処理機能をもつ演算装置と、
（ｖ）前記演算装置に接続されたディスプレイ装置と、
（ｖｉ）前記演算装置に接続された、ユーザーコマンドを入力するためのインターフェース装置と、
を備える装置。
【請求項２】
カメラによって観察される平面チェッカーボードパターンであって、この平面チェッカーボードパターンは、２つの連続する斜め方向の白い四角の中心にある識別可能な特徴を最終的に含むことができ、前記２つの四角は、それらの間にある前記グリッド角を特定する、請求項１に記載の装置。
【請求項３】
前記装置であって、前記平面パターンを観察するための複数の可能な入射点を含むことができ、これらの入射点は、レンズプローブのタイプによって、また前記平面パターンを異なる視点から見るために、複数の形状およびサイズ有することができる、請求項１に記載の装置。
【請求項４】
自動カメラキャリブレーションと、可能なレンズプローブ回転に応じた投影モデル更新と、リアルタイム画像歪み補正とのための方法であって、リアルタイム画像歪み補正は、
（ｉ）キャリブレーションボックスを用いて、チェッカーボードパターンからなる１つまたは複数の画像を取得し、
（ｉｉ）取得されたキャリブレーション画像におけるグリッド角を検出して、画像−平面間対応を確立し、
（ｉｉｉ）少なくとも１つのキャリブレーション画像における前記画像−平面間対応を用いて、焦点距離、主点、スキュー、アスペクト比、およびレンズ歪みのカメラ内部パラメータを決定し、
（ｉｖ）画像／動画を連続的に取得し、
（ｖ）各フレーム時刻において、円形領域とフレーム領域との間の境界輪郭と、三角マークの画像座標とを検出し、
（ｖｉ）（ｖ）の結果を用いて、カメラヘッドに対する前記レンズプローブの角度変位をロバストに決定し、
（ｖｉｉ）（ｉｉｉ）の初期カメラキャリブレーションと、（ｖｉ）の相対回転情報とを所与として、新しい投影モデルを利用して内部キャリブレーション情報を更新し、
（ｖｉｉｉ）（ｖｉｉ）で更新されたカメラキャリブレーションを用いて放射歪曲を補正するリアルタイム画像修正を実行する、
操作を実行する、方法。
【請求項５】
請求項４に記載の方法であって、前記操作（ｉｉ）は、前記チェッカーボードパターンを有する前記平面における角点と、前記キャリブレーション画像における角点との間の対応を自動的に確立する操作からなり、前記ステップ（ｉｉ）は、
ａ．複数の画像処理技術を用いて、中央画像領域内の角またはシード点を検出し、
ｂ．グリッド面内の前記角またはシード点の座標を推測し、
ｃ．平面グリッド座標内の複数の点を、前記カメラ内部パラメータと、レンズ歪みと、前記グリッド面と前記カメラとの間の相対的位置姿勢とを用いて、画像座標内の複数の点に写像する関数を決定し、
ｄ．（ｃ）で決定された前記関数を用いて前記平面グリッド角を前記画像内に投影し、
ｅ．投影された平面グリッド角の周囲の局所的探索を通して前記画像内の新しい角の位置を特定し、
ｆ．確立された角対応の数が十分であると考えられるまで、ステップ（ｃ）から（ｅ）を繰り返すステップをさらに含む方法。
【請求項６】
請求項４に記載の方法であって、前記操作（ｉｉｉ）は、請求項５で確立された角対応を用いて１つのキャリブレーション画像から前記カメラ内部パラメータおよびレンズ歪みを推定するステップからなり、前記操作（ｉｉｉ）は、
ａ．前記平面グリッド内の持ち上げ点座標を前記画像平面内の持ち上げ点座標内へ写像する関数を計算するステップと、
ｂ．焦点距離、主点、アスペクト比、スキュー、歪みパラメータ、およびグリッド面とカメラとの間の相対位置姿勢、回転、および並進を決定するために、（ａ）の前記写像関数を分解するステップと、
を含む方法。
【請求項７】
請求項６に記載の方法であって、さらに反復非線形最適化ステップであって、
ａ．前記推定パラメータの精緻化と、
ｂ．異なるキャリブレーション画像によって提供された前記推定値の統合と、
ｃ．請求項５および６に記載された方法の計算で想定されたモデル関数とは異なるモデル関数を用いた放射歪曲の推定と、
の目的のうち１つまたは複数を達成する反復非線形最適化ステップによって補完される方法。
【請求項８】
請求項４に記載の方法であって、前記操作（ｖ）は、円形画像（ＣＩ）の前記境界輪郭を教示なくリアルタイムで追跡するステップと、前記三角マークを検出するステップとからなり、前記操作（ｖ）は、
ａ．前記ＣＩ境界の最後の推定値を部分的または全体的に含む前記元の画像におけるリング領域を定義するステップと、
ｂ．前記前回の境界推定値が中央の垂直線内に写像された新しい縞状画像をレンダリングするステップと、
ｃ．画像処理技術を用いて前記縞状画像内の境界点を検出するステップと、
ｄ．前記検出された点を前記元の画像へ写像して、ロバスト推定を用いて円錐曲線のフィッティングを行うステップと、
ｅ．前記境界輪郭推定値が収束するまでステップ（ａ）から（ｄ）を繰り返すステップと、
ｆ．前記縞状画像内でわずかに右にずれている補助垂直線を走査することによって、前記三角マークをロバスト検出して、最大強度を有する画素位置を選択し、その画素位置を前記元の画像へ写像するステップと、
を含む方法。
【請求項９】
請求項４に記載の方法であって、前記レンズプローブが回転運動する場合のレンズ位置に関する角度変位の推定における前記操作（ｖｉ）は、
ａ．前記境界輪郭および前記三角マークを検出するステップと、
ｂ．現在の角度変位を推測する確率フィルタへ（ａ）の情報を供給するステップであって、相対回転は、画像情報だけを用いて決定される、ステップと、
を含む方法。
【請求項１０】
請求項４に記載の方法であって、前記操作（ｖｉｉ）は、前記レンズプローブが回転運動するとき、初期基準位置におけるカメラ内部キャリブレーションと、初期基準位置に対する前記レンズプローブの現在の角度変位とを所与として、焦点距離、主点、アスペクト比、スキュー、歪みパラメータを含む現在のフレームの内部キャリブレーションおよび投影モデルを推測するステップからなる、方法。
【請求項１１】
請求項４に記載の方法であって、前記操作（ｖｉｉｉ）は、並列処理に基づく画像ワーピングを通じて前記放射歪曲をリアルタイム補正するステップと、従来の演算処理装置システムを利用するステップとからなり、前記方法は、
ａ．ホスト演算装置から並列処理装置へデータを移動するステップと、
ｂ．前記並列処理装置の画像を処理するステップと、
ｃ．最終処理画像を前記ホスト演算装置に戻すステップと、
ｄ．結果を可視化するステップと、
を含む方法。
【請求項１２】
コンピュータ媒体に実装され、プロセッサによって実行される請求項４〜１１のいずれかに記載の方法ステップを実行するコンピュータプログラム。
【請求項１３】
請求項８，１１のいずれかに記載の方法であって、前記並列処理装置の並列演算機能を用いて加速される、方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、キャリブレーションを行い、レンズスコープの対称軸の周りの回転を決定し、それに応じて投影モデルを更新し、最良の画像品質のために並列処理を用いてリアルタイムで画像の放射歪曲を補正する、高性能な方法、モデル、および装置に関する。
【背景技術】
【０００２】
光学機器は、通常、所与の対象物と幾何学的に一貫性を有する画像であって、３次元における各点が画像内の１点に投影された画像を生成する必要がある。画像は、通常、いくつかの所定のイメージングモデルに従って形成され、本発明の場合には、ピンホールモデルを想定した所定のイメージングモデルに従って形成されている。実用光学系が理想挙動から逸脱すると、結果として得られる画像に収差が導入されてしまう。ここで、収差とはレンズ歪みであり、それは、カラー撮像装置およびグレースケール撮像装置の両方に存在している。その性質は、主に画素の半径方向の幾何学的な変位であり、樽型効果または糸巻き型効果をもたらす。レンズ歪みは、多くの比較的単純な画像解析作業において悩みの種となっている。レンズ歪みは、実像で行われる測定の精度を広範囲にわたって低下させ、通常の遠近投影カメラにおける画素が半径方向に沿って大きく変位することがあり得る。さらに、奥行き表現を変化させることにより視覚的にも影響を及ぼす。
【０００３】
レンズ歪みを解決するための１つの既知の方法は、カメラシステムの完全なキャリブレーションに関連している。完全なカメラキャリブレーションとは、本質的には、カメラの内部パラメータおよび外部パラメータを取得することを意味する［８］。多くのアプリケーションでは、あるキャリブレーション対象に関する内部カメラパラメータおよび相対的配向を完全に補足する必要はない。非測定アプリケーションは、歪みが必ずしもカメラの内部パラメータおよび外部パラメータ全体を推定することなく決定されている点で異なっている。
【０００４】
１９６０年代に、写真測量の領域に従事していたDuane C. Brownは、直線は直線として撮像されなければならないという自明の理に基づいて、レンズ歪みを判断する方法を提案した。この技術は、（Brown, 1971）で発表され、（Fryer and Brown, 1986）で拡張され、「プラムライン（plumb line）」法として知られるようになった。包括的な歴史的見直しが、Clarke and Fryer（1998）においてなされた。この技術は、マシンビジョン業界で採用され、簡易版のプラムライン法が提示され、例えば、Prescott and McLean（1997b）（また特許としては、Prescott and McLean (1997a) XP000683415）、Haneishi et al. (1995a)（また特許としてはHaneishi et al. (1995b) XP000527216）、Poulo and Gorman (1999) US6002525、およびAsari et al. (1999)は、すべて、共線点からなる画像を用いて歪みを補正する類似した公理について記載している。これらの方法は、歪みを推定するだけなので、非計測キャリブレーションと呼ばれることもある。
【０００５】
「プラムライン」に基づくキャリブレーションに内在する問題は、最適化／探索が、（未知の）直線パラメータおよび（これもまた未知の）歪みのパラメータの両方に関して実行されなければならないことである。別の手法が、例えば、Devernay and Faugeras (2001)、Tzu-Hung (2003)、およびBing (1999) EP0895189において採用され、この手法では、歪んだ線分にぴたりと合う線を最小化するために、歪みパラメータを反復的に調整する。客観的誤差と歪みパラメータの間には、十分に検証された数学的関係が存在せず、したがって、分析導関数（analytical derivatives）が全く利用できない。これらは、収束の遅延という結果をもたらし、Swaminathan and Nayar (2000)のような特別な手段が取られない限り、歪み程度の上昇に対して不安定になる可能性がある。Swaminathan and Nayar (2000)は、この非計測手法において通常の歪みのない空間ではなく、歪みのある空間の目的関数を構築した。
【０００６】
別の手法が、Ahmed and Farag (2001)に示唆されており、ここでは、検出された線の曲率を使用して、歪み微分方程式のパラメータを推定する。しかし、シミュレーションの結果は、ノイズ存在下で酷いパフォーマンスを示し、また一方で実際の結果は、量的な評価を欠いていた。
【０００７】
より標準的な歪みをキャリブレーションする方法は、カメラの外部および内部パラメータの同時推定を伴う。Tsaiの方法（Tsai, 1987）は、反復的数値最適化スキームを介して、１つの歪みパラメータと焦点距離などの複数の内部パラメータとを、制御点の集合の３次元位置を所与として、同時に推定することを含んでいる。この手法の欠点は、既知の３次元制御点を必要とし、その見返りとしては、単純同然の歪みプロファイルに対して比較的低い精度しか提供されないことである。この原理のアルゴリズムの変形は、Weng et al.（1992）などの数名の著者によって提案されている。Wei and Ma (1994)は、レンズ歪みに対してより適切なモデルを使用している。これらの方法もまた、既知の３Ｄ制御点を必要とする。
【０００８】
歪み補正後の画像のレンダリングが、Heikkila and Silven (1997)で調査されている。Heikkila (2000)は、３Ｄ制御点か、または２Ｄ制御点からなる複数の画像集合を必要とする類似技術について記載している。また、複数の２Ｄ制御点の集合に基づく別の方法が、Zhang (1998, 2000)およびSturm and Maybank (1999)において提案されている。この技術は、交互線形最小二乗法を通して歪みに対処し、その後、歪みは全ての推定パラメータを含む数値最小化で反復的に調整される。これらの技術の相対的複雑性は、レンズ歪みが含まれることで増加する。
【０００９】
一方で、内部パラメータおよび外部パラメータを完全に補足するのではなく、歪みの除去だけが必要となる場合も多い。例えば、実像における多視点幾何の推定であって、特にレンズの歪みを調節するための技術が開発されている。Zhang (1996)は、歪みパラメータと基本行列とを同時に推定する可能性について調査している。その結果、これは、ノイズが低くかつ歪みが高ければ可能であると結論付けている。Fitzgibbon (2001)（特許としては、Fitzgibbon (2003) GB2380887）、Micusik and Pajdla (2003)、およびBarreto and Daniilidis (2004)が、歪みに対する別のモデルを使用しており、結果として多項式固有値問題（polynomial Eigen value problem）や、より信頼性の高い歪みおよび幾何の推定に至っている。Stein (1997)は、逆の手法を取り、歪みパラメータを推定するための客観的誤差として基本行列の推定における誤差を使用した。
【００１０】
別の歪みキャリブレーションの方法では、歪みのないシーンを選択するために制御点の対応は破棄される。その後、これらのシーンがカメラシステムによって結像され、その上に画像アライメント処理が実施され、歪みが補正される。Lucchese and Mitra (2003)は、歪んだ画像を基準画像に（強度に関して）一致するまでワープ処理する技術について記載している。粗いフィルタを使用して位置合わせを見出す同様の技術が、Tamaki (2002)（特許としては、Tamaki et al. (2002) US2002057345）に記載されている一方で、Sawhney and Kumar (1999)（特許としては、Kumar et al. (1998) WO9821690）は、歪みのない基準画像を必要としない位置合わせ方法について記載している。代わりに、複数の画像が、例えば、いくつかの異なるビューから得られるパノラマ画像などのモザイク画像の生成のために位置合わせされ、歪みが同時に推定される。これらの技術は、Tamaki (2002)において２０分と引用されたように、非常に高い演算オーバーヘッドを有する。
【００１１】
非計測キャリブレーション方法の最後の分類は、歪みによって誘発された周波数領域における高次相関に基づいている。Farid and Popescu (2001)が、ある技術について記載しているが、その性能は通常のカメラキャリブレーション手法と比較して悪く、また画像内容に若干依存しているようである。Yu (2004)は、さらに別の歪みモデルを用いてこの手法を発展させ、ソース画像が通常のキャリブレーション対象である場合に、通常のカメラキャリブレーションを用いて達成される精度接近について報告している。最後に、湾曲したＣＣＤアレイを製造する手段が、Gary (2003) US2003141433によって示唆されている。 レンズ歪みプロファイルが、アレイにおいて線分列毎にコピーされ、こうして結果として得られた画像は歪みがないとされている。
【先行技術文献】
【非特許文献】
【００１２】
【非特許文献１】［１］ D. Kirk and W.-m. Hwu, Programming Massively Parallel Processors: A Hands-on Approach. Morgan Kaufmann, 2010.
【非特許文献２】［２］ T. Yamaguchi, M. Nakamoto, Y. Sato, Y. Nakajima, K. Konishi, M. Hashizume, T. Nishii, N. Sugano, H. Yoshikawa, K. Yonenobu, and S. Tamura, “Camera Model and Calibration Procedure for Oblique-Viewing Endoscope,” in MICCAI, 2003, pp. 373-381.
【非特許文献３】［３］ C. Wu, B. Jaramaz, and S. G. Narasimhan, “A Full Geometric and Photometric Calibration Method for Oblique-viewing Endoscope,” International Journal of Computer Aided Surgery, vol. 15, pp. 19-31, 2010.
【非特許文献４】［４］ N. Fukuda, Y. Chen, M. Nakamoto, and T, “A scope cylinder rotation tracking method for oblique-viewing endoscopes without attached sensing device,” Software Engineering and Data Mining, no. 1, pp. 684-687, 2010.
【非特許文献５】［５］ J. Barreto, J. Santos, P. Menezes, and F. Fonseca, “Ray-based Calibration of Rigid Medical Endoscopes,” in OMNIVIS, 2008.
【非特許文献６】［６］ B. Chapman, G. Jost, and R. Van Der Pass, Using OpenMP: Portable Shared Memory Parallel Programming (Scientific Computation and Engineering Series).The MIT Press, 2008.
【非特許文献７】［７］ L. Velho, A. C. Frery, and J. Gomes, Image Processing for Computer Graphics and Vision.      Springer Publishing Company, Incorporated, 2008.
【非特許文献８】［８］ R. I. Hartley and A. Zisserman, Multiple View Geometry in Computer Vision, 2nd ed. Cambridge University Press, ISBN: 0521540518, 2004.
【非特許文献９】［９］ W. Smith, N. Vakil, and S. Maislin, “Correction of Distortion in Endoscope Images,” IEEE Transactions on Medical Imaging, vol. 11,no. 1, pp. 117-122, Jan 1992.
【非特許文献１０】［１０］ K. VijayanAsari, S. Kumar, and D. Radhakrishnan, “A new approach for nonlinear distortion correction in endoscopic images based on least squares estimation,” IEEE Transactions on Medical Imaging, vol. 18,no. 4, pp. 345-354, 1999.
【非特許文献１１】［１１］ J. Helferty, C. Zhang, G. McLennan, and W. Higgins, “Videoendoscopic distortion correction and its application to virtual guidance of endoscopy,” IEEE Transactions on Medical Imaging, vol. 20, no. 7, pp.605-617, 2001.
【非特許文献１２】［１２］ S. M. Bozic, Digital and Kalman Filtering. London: Edward Arnold, 1979.
【非特許文献１３】［１３］ J. P. Barreto, “A unifying geometric representation for central projection systems,” Comput. Vis. Image Underst., vol. 103, no. 3, pp. 208-217, 2006.
【非特許文献１４】［１４］ A. Fitzgibbon, “Simultaneous linear estimation of multiple view geometry and lens distortion,” in CVPR, vol. 1, 2001, pp.I-125-I-132 vol.1.
【非特許文献１５】［１５］ T. Stehle, M. Hennes, S. Gross, A. Behrens, J. Wulff, and T. Aach, “Dynamic Distortion Correction for Endoscopy Systems with Exchangeable Optics,” in Bildverarbeitungfur die Medizin 2009. Berlin: Springer, 2009, 142-146.
【非特許文献１６】［１６］ J.-Y. Bouguet. Camera Calibration Toolbox for Matlab. ［Online］. Available: http://www.vision.caltech.edu/bouguetj/calib doc/index.html#ref
【非特許文献１７】［１７］ Z. Zhang, “Flexible camera calibration by viewing a plane from unknown orientations,” in ICCV, 1999, pp. 666-673.
【非特許文献１８】［１８］ R. Shahidi, M. Bax, C. Maurer, J. Johnson, E. Wilkinson, B. Wang,J. West, M. Citardi, K. Manwaring, and R. Khadem, “Implementation, calibration and accuracy testing of an image-enhanced endoscopy system,”IEEE Transactions on Medical Imaging, vol. 21, no. 12, pp. 1524-1535, Jan 2002.
【非特許文献１９】［１９］ C. Wengert, M. Reeff, P. Cattin, and G. Sz’ekely, “Fully automatic endoscope calibration for intraoperative use,” in Bildverarbeitungfur die Medizin 2006, 2006, pp. 419-423.
【非特許文献２０】［２０］ J. Mallon and P. F. Whelan, “Which pattern? Biasing aspects of planar calibration patterns and detection methods,” Pattern Recognition Letters, vol. 28, no. 8, pp. 921-930, Jan 2007.
【非特許文献２１】［２１］ J. Barreto, J. Roquette, P. Sturm, and F. Fonseca, “Automatic Camera Calibration Applied to Medical Endoscopy,” in BMVC, 2009.
【非特許文献２２】［２２］ S. D. Buck, F. Maes, A. D’Hoore, and P. Suetens, “Evaluation of a novel calibration technique for optically tracked oblique laparoscopes,” Proceedings of the 10th international conference on Medical image computing and computer-assisted intervention-Volume Part I, pp. 467-474, 2007.
【非特許文献２３】［２３］ A. Fitzgibbon, M. Pilu, and R. Fisher, “Direct least square fitting of ellipses,” IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 21, no. 5, pp. 476-480, May 1999.
【非特許文献２４】［２４］ M. A. Fischler and R. C. Bolles, “Random Sample Consensus: A Paradigm for Model Fitting with Applications to Image Analysis and Automated Cartography,” Commun. ACM, vol. 24, no. 6, pp. 381-395, 1981.
【発明の概要】
【課題を解決するための手段】
【００１３】
本明細書において提供される解決手段は、概して医療および産業において通常使用されるカメラなどの光学素子の完全な幾何学的キャリブレーションと、その後のリアルタイムでの遠近補正画像のレンダリングに基づく。キャリブレーションは、当該入射光の３次元方向に各画素を割り当てる好適な写像関数のパラメータの決定の上に成り立つ。このような解決手段の実用的な実装は非常に簡単であり、容易に取得可能なキャリブレーション対象の画像表示を１つだけ取得するためのカメラであって、特別設計されたキャリブレーション装置の内部に組み立てることができるカメラと、演算プラットフォームの並列実行機能を使用してリアルタイムで動作する、コンピュータにより実装された処理パイプラインと、が必要となる。本発明は、コンピュータ媒体に実装され、プロセッサによって実行されるとすべての方法ステップを実行するコンピュータプログラムについても言及している。
【００１４】
以降の説明、添付の図面、および特許請求の範囲から本発明のより完全な理解が得られるだろう。
【図面の簡単な説明】
【００１５】
添付の図面は、本発明のさらなる理解を提供するために含まれ、本明細書に組み込まれて本出願の一部を構成し、本発明の実施の形態を説明し、明細書と共に本発明の原理を説明する働きをする。
【図１】ＣＣＤカメラ（ａ）にボアスコープを組み合わせた硬性医療用内視鏡である。本システムは、アクセスが困難または制限されている人体体腔の視覚化を可能にし、外科手術および診断の両方に広く使用され、（ｂ）は、３０°斜視関節鏡によって取得された膝関節の内部の画像である。
【図２】平面チェッカーボードパターンからなる１つの画像からの内視鏡カメラのキャリブレーションを示す図である。（ａ）のキャリブレーション画像は、（ｂ）に示された構成を使用して、制約のない方法（カメラと平面との間の相対的位置姿勢については何も仮定しない）で取得される。図示された構成は、照明条件を制御するアクリルボックスで構成されている。平面グリッドは背面照明され、ユーザーが操作しなくても画像四隅のロバストかつ正確な検出が可能になっている。
【図３】カメラをキャリブレーションして画像／動画の放射歪曲を補正するために提案された本システムの別のモジュールを示す模式図である。左側は内視鏡レンズをＣＣＤカメラに接続した後に一度だけ実行される初期化手順に関し、この手順では、図２の構成を用いて１つのキャリブレーション画像の取得が必要となる。右側は内視鏡ビデオストリーム内の放射歪曲を補正するために、各フレーム時刻ｉで実行される処理パイプラインを示している。
【図４】境界輪郭の追跡を示す図である。模式図は、左から右に、時刻ｒで取得された元のフレーム、アフィン変換Ｓを用いたワープ処理画像、デカルト座標から球面座標への変更によって得られた極画像に関する。破線の重なり(overlay)表示は前回の境界推定に関し、×は現在の輪郭における検出点である。Ｓは、Ωｉ−１を単位円内に写像する。網掛け領域は、極歪曲（polar distortion）を受ける画像領域であって、探索が実行される領域である。
【図５】（ａ）は、レンズが完全に回転する場合の主点のプロットｃ、三角マークｐ、境界中心ｗを示す図であり、（ｂ）および（ｃ）は、内視鏡カメラの内部行列を更新するために考慮されるモデルおよび仮定を示す図である。レンズは、平面Ｉ’上に仮想画像を投影し、それはＣＣＤカメラによってＩ内へ結像される。レンズの相対回転は、光軸に平行かつ平面ＩおよびＩ’に直交する軸ｌの周囲となる。ｃは、主点を表し、ｑは、誘発された画像回転の中心である。
【図６】三角マークが正確に検出された場合（ａ）とマークの情報が全くない場合（ｂ）との画像回転中心ｑｉの算出を示す図である。
【図７】角度変位αの関数としての再投影のＲＭＳ誤差を示す図である。青線は、回転補正が全くない場合の誤差を表している（値は、画像四隅のＲＭＳ変位である）。赤線は、我々のモデル（式１１）を用いた場合の再投影誤差を示している。α＝０のとき、プロットされた値はＳＩＣキャリブレーションの再投影誤差である。
【図８】レンズプローブ回転を伴う内視鏡動画像列の放射歪曲補正を示す図である。元のフレームとワープ処理後のフレームとが、それぞれ上列と下列に示され、（ａ）は、初期キャリブレーションが実行される基準位置（ａ＝０）を示している。（ｂ）は、レンズ回転運動の補正なし（左）と補正あり（右）の歪曲結果を比較している。（ｃ）と（ｄ）は、シーンを知覚する上でのＲＤ補正の利点を説明する例である。中央画像領域の寸法を維持するというオプションがあるため、ワープ処理後のフレームは大きくなっている。
【図９】異なる出力画像解像度に対するフレームごとの実行時間を示すグラフである。入力画像は、１６００×１２００画素の解像度を有し、出力画像は、可変矩形サイズを有する。右側にある表は、単純なＣＰＵ解決手段に対して最適化されたハイブリッドＣＰＵ＋ＧＰＵ手法を用いて達成された加速を示している。
【図１０】ハイブリッドＣＰＵ＋ＧＰＵ演算装置の具体的ケースに対する放射歪曲補正アルゴリズムステップを示す図である。緑色の点線ブロックは、割り当てられた装置メモリを表している。ＣＰＵで演算される紫色のブロック以外は、すべての処理がＧＰＵで行われる。
【図１１】ＳＩＣとBouguetキャリブレーションツールボックス（Bouguet calibration toolbox）との間のキャリブレーションパラメータ比較を示す表である。
【図１２】異なる画像解像度に対して最適化されたハイブリッドＣＰＵ＋ＧＰＵ解決手段を実行する全システムを用いて達成されたミリ秒単位での実行時間および１秒あたりのフレーム数（ｆｐｓ）を示す表である。
【図１３】１つの画像（ソフトウェアスクリーンショット）からのキャリブレーションを示す図である。
【図１４】キャリブレーションボックスの模式図である。（１）−外部からの光の入射を避けるために不透明な防水素材を使用して構築されたボックスの上部隔室と、（２）−レンズプローブへのアクセスを可能にするための開口部と、（３）−外部からの光の入射を最小化する不透明材料膜と、（４）−画像取得装置と、（５）拡散光源が取り付けられた、ボックスの下部隔室と、（６）−下部隔室と上部隔室との間のセパレータと、（７）−セパレータの上に置かれた透明または半透明の材料からなる既知のパターンを有する平面グリッドと、を示している。
【発明を実施するための形態】
【００１６】
キャリブレーションと、レンズ回転の場合の投影パラメータのオンライン更新と、リアルタイムでの放射歪曲補正と、のための完全な解決手段を提案する。上記解決手段は、図３の模式図に示されたモジュールとブロックおよび図１４の装置を含み、オペレータは、図２の構成を用いてチェッカーボードパターンからなる１つの画像を取得することによって解決手段を開始する。
【００１７】
本明細書において提示されたシステムは、標準的なＣＣＤカメラに搭載されたボロスコープレンズプローブを主として採用した医学的内視鏡検査のために元々作られた。だが、我々の解決手段は、正確な幾何学的キャリブレーションとリアルタイム画像補正の恩恵を享受することができる、著しい放射歪曲を有するカメラを使用する他の応用領域へ拡張することが可能である。これらの領域の例としては、多くの場合、魚眼レンズ、広角レンズ、小型レンズ、または低品質の光学素子を備えたカメラによって取得された画像／動画に基づく、監視、工業検査、自動車、ロボットナビゲーションなどが含まれる。
【００１８】
本特許出願に開示された装置は、既知のパターンからなる１つの画像を取得するために使用されるキャリブレーションボックス（図１４）によって構成される。取得された画像は、新規キャリブレーションアルゴリズム（図１３）への入力として使用され、それによって、リアルタイムかつレンズ回転に影響されない画像歪み補正を行うために必要なカメラパラメータが提供される。
【００１９】
キャリブレーション画像は、外部（図１４）からの光の入射を避けるために不透明な材料を使用して構築されたボックスを使用して取得される。ボックスは、レンズプローブ（図１４−３）へのアクセスを与えるための１つ以上の開口部を有し、各開口部は外側（図１４−２）からの光の入射を最小限に抑えるゴム膜によって囲まれている。ボックスは、透明または半透明の材料のセパレータ（図１４−６）によって２つの隔室に分割されている。透明または半透明の材料で作られた既知のパターン（図１４−７）を用いた平面グリッドが、上部隔室（図１４−１）内でセパレータの上面に配置される。平面グリッドの背面照明を拡散させる光源が、下部隔室（図１４−５）内に配置される。上部隔室には、目視検査を実施する平均値をより良く再現するために流体を充填することができる。本装置は、中央処理装置（ＣＰＵ）および画像取得用のカメラに接続されたグラフィックスプロセッシングユニット（ＧＰＵ）などの逐次処理および並列処理の機能をもつ演算装置を備えている。可視化およびユーザーコマンドの入力のためには、演算装置に接続されたディスプレイおよびインターフェース装置も必要である。
【００２０】
自動カメラキャリブレーション、レンズ回転の場合の投影パラメータのオンライン更新、および並列処理を使用したリアルタイム画像修正／補正の方法をまとめると以下のようになる。キャリブレーションフレーム内の四隅が検出され、カメラの内部パラメータＫ０と放射歪曲ξの両方が、ユーザーが操作なしで推定される。この簡単な初期化ステップの後、各取得画像に対して図３の右側の処理パイプラインが実行される。各フレーム時刻ｉにおいて、境界輪郭Ωｉだけでなく三角マークＰｉの位置も検出する。検出結果は推定フィルタの入力として使用され、推定フィルタは、キャリブレーションフレームにおける境界Ω０およびマーカー位置ｐ０を所与として、カメラヘッドに対するレンズプローブの可能性のある動きに起因する相対画像回転を推定する。この２次元回転は、新しい適応投影モデルに基づくカメラキャリブレーションの更新の入力として機能する角度αｉおよび固定点ｑｉによってパラメータ化される。最後に、現在の幾何学的キャリブレーションKi，ξを使用して入力フレームをワープ処理して、放射歪曲を補正する。この処理パイプラインは、画像ワーピングや境界検出などの計算集約的なタスクとともにリアルタイムで実行され、ＧＰＵの並列実行機能を使用して効率的に実行される。
【００２１】
Ａ．チェッカーボードパターンからなる１つまたは複数の画像の取得
ユーザーの操作は、図１４に詳しく説明された図２（ｂ）のキャリブレーションボックスを用いた１つまたは複数のキャリブレーション画像の取得に限られている。チェッカーボード面に対するカメラの位置姿勢は任意であり、図１４のボックスは、画像四隅の自動検出のための照明条件を制御するために作られている。取得された画像は、図２（ａ）と似ている。
【００２２】
Ｂ．キャリブレーション画像におけるグリッド角の検出と画像−平面対応関係の確立
直線に曲がりがあるために、強度の放射歪曲のある画像内で四隅を位置づけることには多くの問題が生じうることに注意しなければならない。
フレーム取得後の処理は以下の通りである。
１）歪みがあまり顕著でない中央画像領域における少なくとも１２の角の位置を特定する
この操作は、標準的画像処理技術を使用したヒューリスティックアルゴリズムによって実行される。画像の角ｘとグリッド点ｇとの間の矛盾のない対応関係は、チェッカーボードパターンの四角を数えることによって達成される。
２）画像−平面対応関係は、平面グリッド座標における点を画像座標における点へ写像する関数を推定するために使用され、これがＤＬＴのような手法（セクションＢの式７）を用いた推定値
0005989113.tif 000002
となる。
３）チェッカーボードパターンは、画像周縁において角の仮定を生成する射影変換（homography）を用いて画像平面上に投影される。これらの仮定は、標準的な画像角ファインダを適用することで確認され、精緻化される。
４）ステップ２)および３）は、確立された対応関係が満たされたと考えられるまで繰り返される。
【００２３】
Ｃ．少なくとも１つのキャリブレーション画像における画像−平面の対応関係を用いたカメラ内部パラメータの決定
図３の初期化手順の目的は、レンズプローブが基準位置Ω０，ｐ０にあるときの内部キャリブレーション行列Ｋ０と、放射歪曲ξとを決定することである。カメラキャリブレーションは、文献に何百回も言及されるようなよく研究されるテーマである。最も広く使用されているソフトウェアは、おそらく、平面グリッドの少なくとも３つの画像から汎用カメラをキャリブレーションするZhangの方法［１７］を実行するBouguetツールボックス［１６］である。残念ながら、Bouguetツールボックスは、本提案のアプリケーションの使いやすさ要件を満たしていない。実際には、正確な結果を達成するための入力画像の数は３をはるかに上回り、ＲＤを有する画像内のグリッド角の検出には、かなりのユーザー操作が必要となる。本方法は、少なくとも１画像で機能する。
【００２４】
数名の著者が医療用内視鏡の内部キャリブレーションの具体的な問題について言及している［９］〜［１１］、［１８］、［１９］。しかし、これらの方法では、ＯＲでの使用には非実用的であるか、あるいは、キャリブレーション点の自動検出を可能にするために円形の点パターンを用いるために結果の精度が損なわれる［２０］。
【００２５】
ボアスコープを搭載したカメラは、複雑な光学配置を有する複合型システムである。投影は中央にあり［５］、画像歪みはいわゆる分割モデルによってうまく記述される［２１］、［１４］、［１３］。
【００２６】
Barreto et al.は、上述の投影モデルの後続のカメラは、任意の位置で取得された１つの平面チェッカーボードパターンの画像からキャリブレーションされることが可能であることを示している［２１］。ｇを平面同次座標において表現されたチェッカーボードの１点とすると、ｘは画像平面における対応点であり、Ｈは、平面とカメラとの間の相対的位置姿勢の射影変換符号化である［８］。チェッカーボード平面上の点ｇは、次の関数を介して画像に写像される。
【数１】
0005989113.tif 000003
ここで、Γξは、画像の放射歪曲を考慮した非線形投射関数であり［21］、Ｋ０は、次のような構造を持つ内部パラメータの行列である。
【数２】
0005989113.tif 000004
ここで、f、ａ、ｓ、は、それぞれ焦点距離、アスペクト比、およびスキューを表し、また、c=（cx，cy）Tは、画像主点の非同次座標である。
【００２７】
0005989113.tif 000005
と
0005989113.tif 000006
を、二次ベロネーゼマップにしたがって、点ｇとｘとの持ち上げ表現に対応する次元６を有する同次ベクトルとする［１３］。以下のことが証明される。
【数３】
0005989113.tif 000007
ここで、
0005989113.tif 000008
は、６×６の行列である。各画像−平面対応関係が、持ち上げ射影変換に対して３つの線形制約を課すため、
0005989113.tif 000009
行列は直接線形変換（ＤＬＴ）のような手法を用いて、少なくとも１２点の対応から推定することができる。
0005989113.tif 000010
を所与とすると、カメラ固有値Ｋ０、歪みパラメータξ、および元の射影変換Ｈは、簡単に因数分解することができる。この初期カメラキャリブレーションは、さらに標準的な反復非線形最適化を用いて、精緻化することができる。
【００２８】
Ｄ．円形領域とフレーム領域との間の境界輪郭検出
円形領域とフレーム領域とを分離して画像の輪郭を見つけることは、有意義な映像内容の境界を定めるためだけ（図１（ｂ）参照）でなく、カメラヘッドに対するレンズプローブの回転を推測するためにも重要である。連続するフレームにまたがって境界輪郭と三角マークとを追跡する本提案の手法は、Fukuda et al.［４］およびStehle et al.［１５］による研究に関連している。前者は、従来の画像処理技術を用いて三角マークを抽出することにより、斜視型内視鏡のレンズの回転を推測する。この方法は、円形領域の位置が操作中に変化しないことを前提としているが、それは一般的にありえないことであり、したがって推測がリアルタイムで実行できるかどうかは不明である。Stehle et al.は、半径方向において検出されたエッジ点に円錐曲線を合わせることによって複数のフレームにまたがる境界輪郭を追跡することを提案している。我々のアルゴリズムとの主な違いは、我々は、極画像をレンダリングするためにハイブリッドシリアル＋パラレルの実装を使用し、水平線に沿って探索を行うということにある。この方式により、低い演算コストでロバスト性と正確性を両立させることができる。
【００２９】
レンズプローブがカメラに対して移動するので、輪郭位置は複数のフレームにまたがって変化し、そのため初期オフライン推定の使用は不可能となる。境界検出は各フレーム時刻で実行される必要があり、それは選択されたアルゴリズムの計算複雑性に制約を課すことになる。いくつかの問題があるために円形領域をセグメント化するための単純な手法の使用が不可能になっている。光は、多くの場合、フレーム領域へ広がる。円形画像は、撮像されたシーンおよび照明条件に応じて暗い領域を有することができる。多くの場合、ハイライト、鏡面度、および彩度があり、これらはセグメント化の能力に影響を与える。
【００３０】
決定対象の曲線は、常に、自由度（ＤＯＦ）５を有する楕円Ωであると仮定するのが妥当である［８］。したがって、我々は、ロバスト性と準決定的ランタイムを達成する前に、この形状を用いて複数のフレームにわたって境界を追跡することを提案する。図４（ａ）に示されるように、Ωｉ−１を、フレーム時刻ｉ−１での曲線推定値とする。現在のフレームｉの境界輪郭は、次のように更新される。
１）部分的または全体的に最後の境界推定値を含む元の画像（図４（ｂ））中のリング領域の画定。
２）前回の境界推定値が中央の垂直線内へ写像された新しい縞状画像（図４（ｃ））のレンダリング。
その中心が原点にある単位円内に任意の楕円を写像するアフィン変換が常に存在することは周知である［８］。このような変換Ｓは次式で与えられる。
【数４】
0005989113.tif 000011
ここで、rは、Ωｉ−１の短軸と長軸との間の比率であり、φは、長軸と水平方向との間の角度であり、（wi-1，x，wi-1，y）は、円錐中心Ｗｉ−１の非同次座標である。変換Ｓは、図４（ｂ）の中間結果を生成するために使用され、極画像（図４（ｃ））は、デカルト座標から球面座標への変更を適用することによって得られる。
３）縞状画像における境界点を複数の画像処理技術を使用して検出。
図４（ｃ）の縞状画像内のエッジ点は、（ラプラシアン・ガウシアンなどの）１次元エッジフィルタを適用し、遷移を見つけるために右から左に水平線を走査することによって検出される。
４）球面座標ｘ＝（ρ，θ）で表現されるエッジ点は、式１の関数ＦＳによって、元の像点に写像される。
【数５】
0005989113.tif 000012
現在の円錐境界Ωｉは、ロバスト円錐フィッティング［２３］を用いて、また起こり得る外れ値の悪影響を避けるためにＲＡＮＳＡＣ［２４］を用いて、最終的に推定される。
５）１）〜４）のステップは、境界輪郭推定値が収束するまで繰り返される。
円形画像の境界輪郭を推定するために提示されたステップは、演算装置の並列処理機能を用いて加速される。
【００３１】
Ｅ．レンズが回転運動をするときのレンズプローブのカメラに対する角度変位の推定
セクションＤに記載された方法が正確に収束したのち、境界輪郭は極画像（図４（ｃ））の中央の垂直線内に写像され、わずかに右にずれている補助垂直線を走査することによって三角マークのロバスト検出が可能になり、最大強度を有する画素位置を選択する。三角マークの位置は、各フレームでキャリブレーション位置に対する現在の角度変位を推測する確立フィルタへ供給される。
【００３２】
Ｆ．キャリブレーションとレンズプローブの角度変位に基づく現在のフレームの投影モデルの更新
レンズとカメラヘッドとの間の相対運動は、キャリブレーションパラメータの変化を引き起こし、それによって歪みを補正するための定数モデルの使用が困難になる［５］。この問題の解決策を提案する研究はいくつかあるものの［２］〜［４］、［１５］、［２２］、それらのほとんどは、レンズの回転を決定するために追加の器具類を必要とするという欠点を有する［２］、［３］、［２２］。画像情報にのみ基づいて相対運動を推定する少数の方法には、ロバスト性がないか［４］、あるいはカメラのパラメータ一式を更新することが不可能である［１５］。このセクションでは、実験によって、また光学配置の概念的理解によって、斜視型内視鏡用の、同時駆動される新しい内部カメラモデルを提案する。プローブは、画像平面に直交するが光軸と必ずしも一致しない軸を中心に回転するものとする。カメラモデルを更新するために必要なレンズ回転のパラメータは、境界輪郭と三角マークの画像情報を入力として受信するロバストなＥＫＦによって推定される。我々の動的キャリブレーション方式は、［１５］に対して、２つの重要な利点がある。つまり、（ｉ）投影モデル全体が、ＲＤ輪郭曲線だけではなく、レンズ回転の関数として更新されることと、（ｉｉ）レンズの回転は、三角マークがない場合でも推定することができることである（図１（ｂ）参照）。シーンの３次元点の再投影における説得力のある実験結果によって本手法の有効性が確認される。
【００３３】
カメラモデルにおけるレンズとカメラヘッドとの間の相対回転の影響を評価するために、我々は、レンズプローブを完全に１回転させながら、１０個のキャリブレーション画像を取得した。セクションＣの方法論を用いて角度位置毎にカメラキャリブレーションを推定し、境界Ωと三角マークの両方をセクションDおよびEに説明されたように配置した。図５（ａ）は、主点ｃ、境界中心ｗ、レンズマークｐの結果をプロットしたものである。これらの３つのパラメータは、ほぼ完全な同心軌跡を記述しているため、画像平面に直交する軸を中心とする回転を使用して、レンズ回転のカメラ固有値に対する影響をモデル化するのが合理的と思われる。この考えは、既にWu et al.［３］によって進められているが、そこでは、軸が常に主点を通過すると考えており、これは、我々の実験で示されたように、一般的には、適用できない前提である。
【００３４】
図５（ｂ）、（ｃ）の方式は、カメラ固有値におけるレンズ回転の影響を記述するための本提案のモデルの考えを提示することを目的としている。ボアスコープレンズが、遠端に配置された平面Ｉ’上に仮想画像を投影すると仮定する。Ｉ’は、カメラ接眼部を通して直接見ることによって観察できる画像として考えることができる。Kｃは、この仮想投影の内部行列であり、ｃ’は、光軸が平面と交わる点の位置である。次に、カメラヘッドは、ＣＣＤ面ＩがＩ’と完全に平行でありかつ光軸と直交するように接眼レンズに接続されていると仮定する。Ｉへの投影は、固有値Ｋｈを有し、主点ｃは、ｃ’の画像である。そのため、Ｋｈが内部行列推定値
【数６】
0005989113.tif 000013
であるとき、
【数７】
0005989113.tif 000014
と因数分解できる。
ここで、ｆｃは、ボアスコープレンズの焦点距離であり、また、ｆｈは、メートル単位を画素に変換するカメラヘッドの焦点距離である。
【００３５】
次に、レンズプローブが、軸ｌを中心に角度αだけ回転している場合を考える（図５（ｃ））。ｌは、仮想平面Ｉ’に直交すると仮定されるが、必ずしもレンズ軸と一致してはいない。この場合、点ｃ’は、偏角αを有する円の弧を描き、また、IとＩ’は平行なので、同じことがその画像ｃでも発生する。カメラヘッドおよび回転後のボアスコープによって形成される複合光学系の内部行列は、
【数８】
0005989113.tif 000015
となる。
ここで、Ｒα，ｑ’は、点ｑ’を中心とするαだけの平面回転であり、ここで軸ｌは、Ｉ’と交差している。
【数９】
0005989113.tif 000016
【００３６】
ｑ’の位置は、回転によって不変であることが明らかであり、同じことが画像ｑ〜Ｋｈｑ’にも言える。Ｋｈの具体的構造を考慮して、式８を以下のように書き直すことが可能である。
【数１０】
0005989113.tif 000017
【００３７】
我々は、レンズプローブの回転に対応しており、図８（ａ）の観察結果と一致する内視鏡カメラの投影モデルを導出した。初期化手順では、任意の基準位置（α＝０）でカメラキャリブレーションＫ０，ξを推定する。あるフレーム時刻ｉにおいて、内部パラメータ行列は、
【数１１】
0005989113.tif 000018
となる。
ここで、αｉは、レンズの相対的角度変位であり、ｑｉは、回転中に固定されたままの像点である。放射歪曲がレンズの特性であるため、パラメータξは、カメラヘッドに対する相対運動の影響を受けない。したがって、式２より、カメラ座標フレームにおいて表される一般の３次元点Ｘは、
【数１２】
0005989113.tif 000019
において結像されることになる。
【００３８】
各フレーム時刻における内部パラメータ行列の更新には、相対的角度変位αｉと画像回転中心ｑｉを知ることが必要となる。次に、これらのパラメータがどのようにして境界輪郭Ωの位置と三角マークｐとから推測できるかについて説明する。
【００３９】
ｗｉおよびｗ０を、それぞれ現在の基準フレームにおける境界輪郭ΩｉおよびΩ０の中心とする。同様に、ｐｉおよびｐ０は、その２つの画像における三角マークの位置である。ｗｉ，ｗ０とｐｉ，ｐ０の両方が、推定対象であるパラメータをもつ平面回転Ｒαｉ，ｑｉによって関連付けられると仮定する。この状況は図６（ａ）に図示されており、回転中心ｑｉは、ｗｉ，ｗ０およびｐｉ，ｐ０によって定義される線分の二等分線の交点でなければならないことがわかる。ｑｉがわかれば、回転角度αｉの推定は自明である。三角マークが未知であるとき（存在しないか、または検出できない場合）には、常に、ｑｉの推定には少なくとも３つの識別可能な境界輪郭（図６（ｂ））が必要となる。
【００４０】
制約のある状況を回避して、ｗおよびｐの測定の誤差に対するロバスト性を向上させるために、回転パラメータを推定する確率ＥＫＦ［１２］を使用する。状態の遷移は、運動および定常回転中心のための一定速度モデルを前提としている。
式は、状態変数
【数１３】
0005989113.tif 000020
において線形である。
Ｔは、フレーム取得間隔δｔに依存しており、
【数１４】
0005989113.tif 000021
となる。
αｉおよびｑｉにおいて、計測式は非線形
【数１５】
0005989113.tif 000022
となり、
三角マークの検出に失敗したときにはいつでも最後の２式は破棄される。
【００４１】
本提案のモデルは、複数の異なる角度αで取得されたチェッカーボードパターンの画像上にグリッドの角を再投影することによって有効性が証明された。ＳＩＣキャリブレーション［２１］を基準位置（α＝０）に対して実行して、行列Ｋ０と、ＲＤパラメータξと、グリッド点の３次元座標との決定が可能になった。次に、キャリブレーションパターンに対する相対的位置姿勢を維持するために、レンズプローブを動かさずにカメラヘッドを注意深く回転させた。画像回転中心ｑｉおよび角度変位αｉが、図６の幾何学的形状を用いてフレーム毎に推定された。最後に、式１１を用いて３次元のグリッド点をフレーム上に投影し、実際の画像の角位置との誤差距離を測定した。図７は、複数の異なる角度変位αの再投影誤差のＲＭＳをプロットしている。値は２画素から５画素の間で変化しているが、規則的な挙動(systematic behavior)は観察されない。［２］と［３］は内部および外部キャリブレーションの両方を考慮して追加の計測を採用しているが、我々は内部パラメータのみに焦点を当てた。直接比較はできないが、我々の再投影誤差が、［３］より小さく、また画像中心に近い点のみを考慮した［２］と同等であることは言及する価値がある。以上のことから、またすべての憶測に反して、実験結果は、本提案のモデルの有効性を明らかに証明した。
【００４２】
Ｇ．画像ワーピングを介した放射歪曲リアルタイム補正
このセクションでは、視覚化システムの最終的な出力である正確な遠近画像のレンダリングについて説明する。［７］で指摘されているように、特定の変換による画像の効率的なワーピングは、逆写像法（inverse mapping method）を用いて行われるべきである。
したがって、所望の歪みのない画像内の点yを、元の歪曲したフレーム内の点ｘへ写像する関数Ｆを導出しなければならない。式１１から、
【数１６】
0005989113.tif 000023
が得られる。
Ｋｙは、歪みのない画像の特定の特徴（例えば、中心、解像度）を特定し、Ｒ＿αｉ，ｑ”は、ワーピング結果を元の向きに戻すように回転させ、ｑ”は、回転中心ｑｉの逆投影
【数１７】
0005989113.tif 000024
である。
【００４３】
中央領域内のオブジェクトの縮尺を維持するために、画像周辺部を拡大して、歪みのない中央領域の大きさを保つ。これは、セクションＤの境界輪郭の半径からワープ処理された画像の大きさｕを計算することにより行われる。ｒｄを原点と点Ｋｉ−１ｐ０との間の距離（歪曲した半径）とする。所望の画像サイズｕは、ｕ＝ｆｒｕ
によって与えられ、ここでｆはカメラの焦点距離であり、ｒｕは歪みのない半径である。したがって、行列Ｋｙは、
【数１８】
0005989113.tif 000025
でなければならず、ワープ処理された画像は、画像回転中心ｑｉを写像した軌跡である。図８は、内視鏡動画像列のうちいくつかのフレームのＲＤ補正結果を示している。これらの例は、シーンの認知における改善と、画像形状の補正の間にレンズプローブ回転を考慮することの重要性（図８（ｂ））とをはっきりと示している。
【００４４】
補正された画像のレンダリングには、高性能計算リソースによるリアルタイムのデータ処理が必要となる。並列処理装置で（この場合はＧＰＵにおいて）ＲＤを補正するために、我々のアルゴリズムのいくつかの部分を並列化することを提案する。
４つの手法、すなわち、
（ｉ）ＣＰＵベースの単純で純粋な解決手段と、
（ｉｉ）ｘ８６ファミリの従来の汎用市販マルチコアＣＰＵなどの共有メモリアーキテクチャに使用することができるOpenMP2directives［６］を用いた仮想ＣＰＵバージョンと、
（ｉｉｉ）最適化されていないハイブリッドＣＰＵ＋ＧＰＵバージョンと、
（ｉｖ）最適化されたＣＰＵ＋ＧＰＵ複合バージョンと、を比較することによって、我々のハイブリッドシリアル＋パラレル解決手段の効果を評価する。
【００４５】
完全なランニングシステムの実験構成は、２．４０ＧＨｚのインテルCoreTM2 QuadVIII（登録商標）のＣＰＵと、４ギガバイトのＲＡＭと、GeForce 8800 GTXのＧＰＵとを搭載したワークステーションからなる。ＧＰＵは、各々が１．３５ＧＨｚで動作する１２８個のストリームプロセッサ（ＳＰ）と、合計で７６８メガバイトのビデオＲＡＭとを有する。ソースコードは、ＧＣＣバージョン４.３.４にコンパイルされ、ＣＵＤＡバージョン３．４が使用された。
【００４６】
図９は、上述の複数の解決手段によって達成された（１つのフレームを処理するのに必要となる）実行時間を比較している。ＧＰＵがＣＰＵよりもどれだけ速いかの測定である加速表(speedup table)には、単純なＣＰＵと最適化されたハイブリッドＣＰＵ＋ＧＰＵとの間の実行時間の対比が示されている。
図９に提示された比較は、ＣＰＵがＧＰＵほど効率的に大きなワーピング操作（アルゴリズムの中で最も時間のかかるタスク）を処理することができないことを示している。表は、３０００×３０００画素以上の出力サイズをもつ画像の場合、ＧＰＵが実行を加速するだけではなく、キャッシュミスペナルティの増加によりＣＰＵの性能が低下することを示している。
【００４７】
図１０は、放射歪曲補正アルゴリズムのシリアル＋パラレルハイブリッド実行ステップと、ＣＰＵおよびＧＰＵで実行するハイブリッドプログラムの開発において採用された実装方式の詳細を示している。ＧＰＵでのアルゴリズムの並列化は、主に３つのステップに分かれている。
（ｉ）画像変換：画像は、ＲＧＢチャンネルに分割され、グレースケール変換が実行される。
（ｉｉ）境界推定：グレースケール画像は、並列処理ユニットメモリにバインドされ、輪郭強調カーネルを伴う球面座標への写像（式１）が起動される。
（ｉｉｉ）放射歪曲補正：画像変換段階からの３つの画像チャンネルが、並列演算処理ユニットのメモリにバインドされ、ＲＤ補正カーネルが起動される。
【００４８】
ＣＰＵで計算される中間ステップ（図１０においては分離されている）がある。極画像の各ラインは、楕円がフィッティングされた輪郭点を求めて走査される。楕円のパラメータがセクションＦで説明されたＥＫＦに供給され、キャリブレーションパラメータが更新されて、歪み補正カーネルに引数として渡される。このステップは、関係するプロセスのシリアル特性のため、ＧＰＵではなくＣＰＵで実装されている。
【００４９】
最適化されたハイブリッドＣＰＵ＋ＧＰＵ解決手段は、調停（coalescence）として知られる、スレッドのグループ毎に１つのメモリアクセスの実行を可能にするデータプリアライメントに基づいている［１］。そうではなく、ＧＰＵによって処理されるデータがアライメントされていない場合には、調停されたアクセスが行われず、各スレッドのグループに対していくつかのメモリトランザクションが発生し、実行中のカーネルの性能を大幅に低下させる。
【００５０】
画像のアルファチャンネルは使用されていないが、完全に調停したアクセスを実行するためには、メモリレイアウト要件を満たすことが必要である。転送対象のデータ量の増大により、転送時間で１０．６％の損失が取り込まれるが、一方で達成された調停によりカーネルの実行時間は６６．１％低減する。要するに、調停の実装によって、最適化されていないハイブリッドＣＰＵ＋ＧＰＵの実装と比較して、演算時間が３８．７％節約される。
【００５１】
図１２の表は、１５０フレーム動画像列の処理中に、図１０の各ステップにかかる時間を示している（時間は、解像度毎に処理されるすべてのフレームの平均時間の値を表している）。画像変換（Ｉ．Ｃ．）時間には、ホストから装置のグーバルメモリへの入力画像転送時間も含まれる。境界推定（Ｂ．Ｅ．）時間には、ＣＰＵで実行される、輪郭点の抽出、楕円フィッティング、およびＥＫＦ推定の各処理も考慮されている。ここで留意すべきは、入力画像の球面空間と輪郭強調フィルタとへの写像は、両方ともＧＰＵで計算され、表ＩＩに示された境界推定演算時間のうちの１〜３ミリ秒を費やすことである。放射歪曲補正（Ｒ．Ｄ．Ｃ．）時間には、装置のグローバルメモリからホストへの転送時間も考慮されている。
【図２ａ】
0005989113.tif 000026
【図２ｂ】
0005989113.tif 000027
【図４】
0005989113.tif 000028
【図６（ａ）】
0005989113.tif 000029
【図６（ｂ）】
0005989113.tif 000030
【図１３】
0005989113.tif 000031
【図１４】
0005989113.tif 000032
【図１】
0005989113.tif 000033
【図３】
0005989113.tif 000034
【図５】
0005989113.tif 000035
【図７】
0005989113.tif 000036
【図８】
0005989113.tif 000037
【図９】
0005989113.tif 000038
【図１０】
0005989113.tif 000039
【図１１】
0005989113.tif 000040
【図１２】
0005989113.tif 000041