Harrisコーナー検出器 (Harris Corner Detector)

記事を共有する:

1. Harrisコーナー検出器 (Harris Corner Detector)の概要

Harrisコーナー検出器 (Harris Corner Detector)は,古典的なコーナー検出の代表手法である[Harris and Stephens, 1988].近傍パッチとのSSDを計算することでコーナー点を推定する『Moravecコーナー検出(2節)』を,頑健性と計算コストと両面で改善したコーナー検出器がHarrisコーナー検出器である.

この記事では,Harrisコーナー検出器(3節)と,その改善元の手法であるMoraveコーナー検出器(2節)を紹介する.

コンピュータビジョンでは「近傍画素の全てに対して,画素値が大きく変化している点」をコーナー点(corner point)と定義する.また,この定義に沿ったアルゴリズムを設計することで,コーナー検出も行う.

検出したコーナー点は,局所特徴(Local Features)の特徴記述子(feature descriptor)ベクトルを計算するために,前準備として検出しておくのが,関心点(interest point)であり,コーナー点も注目点の1種類である.

例えば,画像同士の位置合わせ(Registration)やモザイク合成(Image Stitching)などの応用において,コーナー点や注目点同士の対応付け(マッチング)が行われる [Sur, 2011].また,動画処理などでも,物体を追跡する際に,各画像から検出しておいたコーナー点やその他の注目点は,対応付け目的で重用される.

1.1 記事の構成

2節以降,記事の構成は以下の通りである:

  • 2節:Moravecコーナー検出器
    • 2.1「平坦/エッジ/コーナー」の判定.
    • 2.2 エネルギー関数:小移動前後のSSD
    • 2.3 コーナーネスを用いた検出
    • 2.4 Moravecコーナー検出器の課題
  • 3節:Harrisコーナー検出器
    • 3.1 改善された連続的局所エネルギー関数
    • 3.2 行列Mの解釈:「コーナー/エッジ/平坦」の識別
    • 3.3 コーナーネス計算の効率化
    • 3.4 Harrisコーナー検出のアルゴリズム
  • 4節:まとめ

Harrisコーナー検出器は,Moravec検出器の発展版手法である.従って,先に2節でMoravec検出器を紹介し,その後に3節でHarrisコーナー検出器について紹介する.

2. Moravecコーナー検出器

Moravecコーナー検出器 [Moravec, 1980], [Moravec, 1981] は,参照点$(x,y)$を中心とする局所パッチを,近傍8方向に$(u,v)$だけ平行移動させた時のSSD(Sum of Squared Differences)誤差を合計した「局所変化量$E(u,v)$」について極小点(最小値)を求めることにより,画像上のコーナー点を検出する手法である.

2.1 コーナー検出の基本:「平坦/エッジ/コーナー」の判定.

図1. Harrisコーナー検出の基本アイデア:SSD誤差量による3種分類

1節冒頭でも述べたように,移動後のパッチと移動前のパッチのSSD誤差が大きい点ほどコーナーに相当している.ここで「平坦/エッジ/コーナー」の3種類に分けて,変化量の合計を考えてみたい (図1).すると,「小移動前後でのパッチSSD誤差の量」が,これらの3種類を識別する値として使用できることがわかる.

小移動前後のパッチ間でSSD誤差量が大きいほど,「コーナー領域[図1-(c)]」に相当しており,パッチ同士も画素の様子が異なっているので特徴マッチングに使いやすいパッチである.逆に,小移動前後のパッチSSD誤差量が小さいほど,「平坦領域[図1-(a)]」や「エッジ領域[図1-(b)]」となり,特徴マッチングに向かない局所パッチであるとわかる.

まとめると,「近傍への小さな平行移動の前後で,パッチ間のSSD誤差が大きく出る点がコーナー点」という発想が,コーナー検出器の基本的アイデアである.また,このアイデアは「コーナー点同士は,点周辺のパッチ特徴を用いるとマッチング処理をおこないやすい」という特徴マッチング(Feature Matching)の話にも繋がる.

2.2 エネルギー関数:小移動前後のSSD

Harrisコーナー検出:局所窓の局所移動
図2. 変化量$E(u,v)$の計算.注目点$(x,y)$周囲の局所窓$\mathcal{W}$について,移動前パッチと$(u,v)$だけ平行移動後のパッチ間で,SSDを計算した値としてE(u,v)を算出.

入力画像$I$の座標$(x,y)$における輝度値を$I(x,y)$とする.$(x,y)$の周辺$6 \times 6$のパッチ窓$\mathcal{W}$を,ベクトル$(u,v)$により少しだけ平行移動することを考える(図2).

ここで,$(x,y)$がコーナーかどうかを判定するためのエネルギー関数$E(u,v)$は,$(u,v)$だけ平行移動した前後におけるSSD(パッチ間での輝度値の合計変化量)として,以下のように定義できる:

\begin{equation}
E(u,v) = \sum_{(x, y) \in \mathcal{W}} w(x, y) [I(x+u, y+v)-I(x,y)]^2 \tag{2.1}
\end{equation}

ここで$w(x,y)$は$(x,y)$近傍のパッチ窓$\mathcal{W}$内の画素だけを用いるための重み係数(窓関数)である.窓関数としては,1 or 0の矩形関数ではっきりパッチ内の画素だけ残すか,ガウス関数などで滑らかにパッチ内だけ残したりもする.

2.3 コーナーネスを用いた検出

コーナー点判定をしたい注目点$(x,y)$の,近傍全ての方向に$(u,v) \in \mathcal{S}$だけ平行移動させるとする(※ 45度ずつ近傍8画素方向がデフォルト).

この時,2.2節のとおり,平行移動した時の変化量$E(u,v)$の合計である指標値$c_{M}$を,最小化する点をコーナー点として識別できる:

\begin{equation}
c_{M} = \min \{ E(u,v) \mid (u,v) \in \mathcal{S} \} \tag{2.2}
\end{equation}

この「画素$(x,y)$のコーナー点度合い」を表す指標$c_{M}$を,コーナーネス(cornerness)と呼ぶ.

入力画像$I$の各座標$(x,y)$において式(2.2)からコーナーネス値$c_{M}$を計算する.そして,画像全体から,コーナーネス値の高い座標$(x,y)$を閾値処理などで選択することで,コーナー点を検出できる.

以上がMoravecコーナー検出のアルゴリズムである.

2.4 Moravecコーナー検出器の課題

$E(u,v)$の最小値だけを保持して,かつ離散的な8方向の平行移動のみしか考慮できていない$c_{M}$は,$(x,y)$の近傍局所における変化を正確に表現できた指標ではない.なぜなら,8方向しか考慮していないので,図1-(c)の「全方向でSSDが大きい」が満たせているかが判定できていないからである.よって,コーナー検出結果がさほど安定しない課題があった.

また,各画素で個別に,近傍8方向で変化量(SSD)をそれぞれ計算するゆえ,計算コストも高くリアルタイム性に欠けている課題もあった.

そこでMoravecコーナー検出器の改善版として,3節のHarrisコーナー検出器が提案される.Harrisコーナー検出器では,安定性と高速性が達成され,当時の標準的コーナー検出手法として定着することとなる.

3. Harrisコーナー検出器

Harrisコーナー検出器 [Harris and Stephens, 1988]は,2節のMoravecコーナー検出の問題点を改善したコーナー検出器である.具体的には「コーナー検出器とエッジ検出器を統合した検出器」として提案された.

平行移動量が小さいという仮定を置き,8方向に平行移動したパッチの代わりに画像の(偏)微分値を直接使用する.そして,エネルギー関数$E(u,v)$が,画像の2次モーメント行列$M$により局所的に近似されている.この近似表現のおかげで,コーナーの識別判定[図1-(c)]が,$M$の固有値を通じて安定して行える(3.2節).さらにその行列$M$の固有値も,特異値行列分解無しで計算できるようにしているおかげで,計算効率的にも向上している(3.3節).

3.1~3.3はアルゴリズムを構成する個々の部品の導出を行っている.よって,最後の3.4節にて全体を通したアルゴリズムを紹介し,3.1~3.3節を復習して全体的に再度俯瞰しやすくなる構成にした.

3.1 改善された連続的局所エネルギー関数

Harrisコーナー検出器は,局所エネルギー関数$E(u,v)$ (SSD誤差)を近似し,その近似式の一部分を成す行列$\bm{M}$の固有値の範囲次第で,コーナーネスを用いたコーナー検出が可能となる仕組みに改良された.

Moravecコーナー検出器と同様に,局所パッチの平行移動前後の画素値変化量(SSD)を用いて,コーナーおよびエッジを検出したい.しかし,Harris検出器では,近傍8方向移動を通じたパッチ計算は行わず,2次偏微分画像から直接$E(u,v)$の最大化が可能なアルゴリズムになっている.

まず,Harrisコーナー検出器の局所エネルギー関数$E(u,v)$を以下のように定義する:

\begin{align}
E(u,v) = \sum_{(x,y) \in \mathcal{W}} w(x,y)[I(x+u,y+v)-I(x,y)]^2 \tag{3.1}
\end{align}

この式(3.1)は式(2.1)と同じ式を意味しているが,Harrisコーナー検出器では平行移動量(u,v)が非常に小さいと仮定しており,$(u,v)$に対して拘束条件 $u^2+v^2 = 1$ を加えている.

この「$(u,v)$の移動量がとても小さい」という仮定から,式(3.1)の$I(x+u,y+v)$部分を,テイラー級数による関数近似により以下のとおり変形できる:

\begin{equation}
I(x+u,y+v) \approx I(x,y) + I_x(x,y)u + I_y(x,y)v \tag{3.2}
\end{equation}

ここで,式(3.2)を用いると$E(u,v)$全体を以下のように近似できる:

\begin{equation}
E(u,v) \approx \sum_{(x,y) \in \mathcal{W}} w(x,y) \left[ I_x(x,y) u + I_y(x,y) v \right]^2 \tag{3.3}
\end{equation}

この1次近似である式(3.3)は,行列$\bm{M}$を用いて以下のように書き直すことができる:

\begin{equation}
E(x,y)=
\begin{bmatrix}
u & v
\end{bmatrix}
\bm{M}
\begin{bmatrix}
u \\
v
\end{bmatrix}\tag{3.4}
\end{equation}

ここで,式(3.4)の$\bm{M}$は,以下のような成分からなる行列である:

\begin{align}
\bm{M} &= \sum_{(x,y) \in \mathcal{W}}
\begin{bmatrix}
I_x^2(x,y) & I_x(x,y) I_y(x,y) \\
I_x(x,y) I_y(x,y) & I_y^2(x,y)
\end{bmatrix} \\
&= \begin{bmatrix}
\sum_{(x,y) \in \mathcal{W}} I_x^2(x,y) & \sum_{(x,y) \in \mathcal{W}} I_x (x,y) I_y(x,y) \\
\sum_{(x,y) \in \mathcal{W}} I_x(x,y) I_y(x,y) & \sum_{(x,y) \in \mathcal{W}} I_y^2(x,y)
\end{bmatrix}\tag{3.5}
\end{align}

以上の導出により,$E(x,y)$の近似として行列$\bm{M}$を式(3.5)で計算し,コーナー検出の指標に使用すれば良いとわかる.そして,式(3.5)は,1次微分画像と2次微分画像の画素値から計算できることがわかる.

つまり,画素ごとに周辺の値から計算を行う手間が生じていたMoravecコーナー検出と違い,移動させる局所近傍を周辺1画素に絞ったHarrisコーナー検出では,順に画像全体のフィルタリングを行っていくだけで,楽に各処理を進めることができるようになった(3.4節で,全体の処理手順をまとめる).

3.2 行列Mの解釈:「コーナー/エッジ/平坦」の識別

この節では,行列$\bm{M}$の幾何的解釈を行い,$\bm{M}$の固有値がどういった値の場合に「平坦/エッジ/コーナー」に相当するかを見ていきたい.

$\bm{M}$は画像$\nabla I$における$(x,y)$近傍の曲率,そのパッチ内の和となっており,これは入力画像中の局所的なパッチ内の2次モーメントの和に相当している(式(3.5)).つまり式(3.4)を見ると,$E(u,v)$の2次関数を水平断面で見ると楕円形を成している.よって$\bm{M}$を特異値分解により対角化すると,固有値$\lambda_1, \lambda_2$が楕円のそれぞれ第1第2主軸方向となり,$\lambda_1, \lambda_2$の値次第で,コーナーなのかエッジなのかを判定できる.

図3は,$\bm{M}$の固有値の値次第で,どの領域がそれぞれ「コーナー/エッジ/平坦」領域の画像パッチになっているかを表した図である(元論文のFigure 5に相当).

Harrisコーナ検出:コーナー/エッジのMの固有値による識別
図3. Harrisコーナ検出:コーナー/エッジのMの固有値による識別:
水色の円・楕円は,それぞれの領域に対応する行列$\bm{M}$が示すパッチ内の幾何的形状.

図3のように,$\lambda_1, \lambda_2$の値の大小関係から「コーナー/エッジ/平坦」領域の判別が可能となる.従って$\bm{M}$の固有値分解を行うとその固有値からコーナーもエッジも検出できる.

$\bm{M}$の2つの固有値$\lambda_1, \lambda_2$が大きな値を持つときは,点$(x,y)$がコーナーに相当する.一方で,$\lambda_1 \gg \lambda_2$もしくは$\lambda_1 \ll \lambda_2$ であれば,エネルギー関数の形はそれぞれ縦長,横長になるのでエッジ領域に相当する.

しかし,このまま全画素で$\bm{M}$の特異値行列分解を実行して固有値を求めようとしてしまうと,計算コストが高くなってしまう.そこで,Harrisコーナー検出器では,代わりに「コーナネス計算の効率化」を行うことで固有値計算を目的とする特異値行列分解を回避する.

3.3 コーナーネス計算の効率化

Harrisコーナー検出器では,パラメータ$k >0$の範囲の場合は$\lambda_1$と$\lambda_2$の値は計算しないまま,以下のコーナーネス$c_{H}$を代わりに計算する:

$$ c_{H} = \mathrm{det}(\bm{M}) – k \mathrm{tr}^2(\bm{M}) \tag{3.6}$$

(※ $k$の値には経験的に決められた0.04-0.06の範囲の値を用いる)

この(3.6)の右辺の値は,固有値の一般的な性質から$\lambda_1, \lambda_2$を用いて変形すると以下のようにも表せる:

$$c_{H} = \lambda_1\lambda_2 – k \cdot (\lambda_1+\lambda_2) \tag{3.7}$$

しかし,上記の式(3.6)の $\mathrm{det}$ と $\mathrm{tr}$ を用いた$c_{H}$の計算を行なった方が,行列の固有値(特異値)分解を行わないで済む計算上の利点がある,従って,Harrisコーナー検出器では式(3.6)の方をコーナーネス$c_{H}$の計算に採用している.

また,3種のパッチ領域は,コーナーネスの値を用いると以下のように識別できる:

  • $c_{H}$が大きい非負の値 → コーナー
  • $c_{H}$が小さい非負の値 → 平坦
  • $c_{H} < 0$ → エッジ

よって,$c_{H}$を全画素で求めたのちに,局所的な山のそれぞれにおける極大値点を求めることで,安定したコーナー点検出を達成することができる.

3.4 Harrisコーナー検出器のアルゴリズム(おさらい)

最後に3.1~3.3節を踏まえて,Harrisコーナー検出器のアルゴリズムを処理順に短く示すことで,全体の処理を振り返っておきたい.

Harrisコーナー検出では,グレースケール画像$I$に対して,以下の3つの処理(3.4.1~3.4.3)を順に行うことにより,コーナー点の検出を行う:

3.4.1 Mの計算

各画素$(x,y)$におけるコーナーネス値を計算するために,まずは各画素の$\bm{M}$を計算したい.そこで,まず以下の3つの手順で前処理を行う(3.1節):

  1. 微分画像$I_x,I_y$を微分フィルタで計算する.
  2. $\bm{M}$の構成要素である「微分値の積」$I_x^2,I_y^2,I_{xy}$の画像を,$I_x,I_y$を用いて計算する.
  3. 式(3.5)の計算:ガウシアンフィルタ$g$を畳み込み,各画素の『「微分値の積」の和』に相当する $g *(I_x^2), g*(I_y^2), g*(I_{xy})$をそれぞれ計算する(ガウシアンフィルタの窓が,パッチ和を取る局所窓$\mathcal{W}$のサイズに相当する).

これら普通の画像フィルタリングだけで計算が済む1~3の前処理を行なったのち,(手順2で計算した微分画像値ではなく)手順3で計算した『「微分の積」の和』の値を用いて $\bm{M}$ を計算する:

\begin{equation}
\bm{M}(x,y) = \begin{bmatrix}
I_x^2(x,y) & I_x(x,y) I_y(x,y) \\
I_x(x,y) I_y(x,y) & I_y^2(x,y)
\end{bmatrix}\tag{3.8}
\end{equation}

この式(3.8)の計算は,式(3.5)と等価であり,これは「全平行移動方向の値を(連続的に)加える」ことに相当している.Moravecコーナー検出器(2節)では,SSDは8方向の$(u,v)$のみSSD和に加えていたのが,連続的に全シフト方向でSSDの和を算出しており,どの局所移動方向へのコーナーも検出できることがHarrisコーナー検出器のポイントである(ただし,式(3.1)で $u^2+v^2 = 1$の拘束を課していることには注意).

3.4.2 コーナーネス値を全画素で計算 (3.3節)

$\bm{M}$の$\mathrm{det}(M)$と$\mathrm{tr}(M)$から,コーナーネス$c_{H}$を式(3.6)により各画素で計算する.

これにより,各座標の画素値にコーナーネスの応答値が収納された画像が手に入る.

3.4.3 非極大値抑制により,局所極大値をコーナー点として出力.

コーナーネスの応答値画像では複数の山が出ており,手順(2)の時点では,各山の点のうちどの1点がコーナーであるかがまだ未確定である.

そこで,以下の二段処理で,応答値画像から各局所のコーナーネス極大点のみをコーナー点として残す処理を行う:

NMS処理の手順

  1. 事前に決めておいた閾値$t$以上の$c_{H} > t$を示す座標のみにc_{H}の値を残した,「閾値処理された応答値画像」を計算.
  2. 「閾値処理された応答値画像」に,NMS(Non-maximal suppression, 非極大値抑制) 処理を行い,局所の各山におけるコーナーネスの極大点のみを残す.残った極大点の集合が,最終的なコーナ点として出力される.

これで,Harrisコーナー検出の全ての処理は完了である.

4. まとめ

この記事では,前身であるMoravecコーナー検出器について紹介したのち(2節),その改善版であるHarrisコーナー検出器を紹介した.

3.4節で全体の手順をおさらいした通り,Moravecコーナー検出では離散的8方向のみコーナー判定をしていたのを,「連続化」して計算コストも下げた手法である.SSD和 $E(u,v)$ の計算において,超局所しか加味しない近似計算のおかげで,微分画像群から各画素のコーナーネス値を低計算コストで推定できる.そのコーナーネス値を用いることで,閾値処理とNMSだけでフィルタリング処理によりコーナー点推定できる.こうして計算効率も良い手法となり,その後Harrisコーナー検出は,標準的なコーナー検出手法として定着した.

参考書籍

References

  • [Farinella et al., 2013] Farinella, Giovanni Maria, Sebastiano Battiato, and Roberto Cipolla, eds. Advanced Topics in Computer Vision. Springer London, 2013.
  • [Harris and Stephens, 1988] Harris, Chris, and Mike Stephens. “A combined corner and edge detector.” Alvey vision conference. Vol. 15. No. 50. 1988.
  • [Moravec, 1980] Moravec Hans Peter. Obstacle avoidance and navigation in the real world by a seeing robot rover. Diss. Stanford University, 1980.
  • [Moravec, 1981] Moravec Hans Peter “Rover Visual Obstacle Avoidance.” IJCAI. Vol. 81. 1981.
  • [Sur, 2011] CVPR2011 tutorial – Tools and Methods for Image Registration – Frederic Sur “Similarity and affine invariant point detectors and descriptors”.

外部参照リンク