勉強しないとな~blog

ちゃんと勉強せねば…な電気設計エンジニアです。

OpenCVやってみる-22. エピポーラ幾何

今回は、OpenCVチュートリアルの以下のページを参考に進めてみたいと思います。

OpenCV: Epipolar Geometry

エピポーラ幾何の理論はチュートリアルサイトの説明を参照で。

今回のテーマ画像

今回は、春のパン祭り台紙の2視点分の画像でマッチング、Epipolar lineの描画をしてみたいと思います。
似たような特徴点が複数あるので難しそうと思いましたが、RANSAC等適用すればうまくマッチングできるのではと。

import cv2
img1 = cv2.imread('harupan_210402_1.jpg')
img2 = cv2.imread('harupan_210402_2.jpg')
img1= cv2.resize(img1, None, fx=0.2, fy=0.2, interpolation=cv2.INTER_AREA)
img2= cv2.resize(img2, None, fx=0.2, fy=0.2, interpolation=cv2.INTER_AREA)

img12 = np.hstack((img1, img2))
cv2.imshow('Image', img12)
cv2.waitKey(0)

f:id:nokixa:20211016205646p:plain

画像サイズはいずれも806x605になっています。

>>> img1.shape
(806, 605, 3)
>>> img2.shape
(806, 605, 3)

特徴点検出、特徴量算出

SIFTを使いました。

sift = cv2.SIFT_create()
kp1, des1 = sift.detectAndCompute(img1, None)
kp2, des2 = sift.detectAndCompute(img2, None)

img1_kp = cv2.drawKeypoints(img1, kp1, None, flags=cv2.DRAW_MATCHES_FLAGS_DRAW_RICH_KEYPOINTS)
img2_kp = cv2.drawKeypoints(img2, kp2, None, flags=cv2.DRAW_MATCHES_FLAGS_DRAW_RICH_KEYPOINTS)

img12_kp = np.hstack((img1_kp, img2_kp))
cv2.imshow('Keypoints', img12_kp)
cv2.waitKey(0)

f:id:nokixa:20211017010320p:plain

マッチング

ここでは総当たりマッチングを使います。

bf = cv2.BFMatcher()
matches = bf.match(des1, des2)
matches = sorted(matches, key = lambda x:x.distance)

img_matches = cv2.drawMatches(img1, kp1, img2, kp2, matches[:200], None, flags=cv2.DrawMatchesFlags_NOT_DRAW_SINGLE_POINTS)
cv2.imshow('Matches', img_matches)
cv2.waitKey(0)

f:id:nokixa:20211016011731p:plain

いまいちな結果なので、ratio testも入れてみます。

matches_rt = bf.knnMatch(des1, des2, k=2)
good_matches = []
for m,n in matches_rt:
    if m.distance < 0.75*n.distance:
        good_matches.append(m)

good_matches = sorted(good_matches, key = lambda x:x.distance)
img_good_matches = cv2.drawMatches(img1, kp1, img2, kp2, good_matches[:200], None, flags=cv2.DrawMatchesFlags_NOT_DRAW_SINGLE_POINTS)
cv2.imshow('Good matches', img_good_matches)
cv2.waitKey(0)

f:id:nokixa:20211016014617p:plain

これでもまだまだですが、この結果を使っていきます。

ちなみに、good_matchesには1169個のマッチが含まれていました。

>>> len(good_matches)
1169

Fundamental matrix計算

Fundamental matrixの計算をやってみます。
この過程でマッチング結果のinlier、outlierを判定してくれるので、 それも見てみます。

pts1 = []
pts2 = []

for i,m in enumerate(good_matches):
    pts1.append(kp1[m.queryIdx].pt)
    pts2.append(kp2[m.trainIdx].pt)

import numpy as np
pts1 = np.int32(pts1)
pts2 = np.int32(pts2)
F, mask = cv2.findFundamentalMat(pts1, pts2, cv2.FM_LMEDS)

pts1_inliers = pts1[mask.ravel() == 1]
pts2_inliers = pts2[mask.ravel() == 1]

good_matches_inliers = [m for i,m in enumerate(good_matches[:200]) if mask[i,0] == 1]

img_good_matches_inliers = cv2.drawMatches(img1,kp1,img2,kp2,good_matches_inliers[:200],None,flags=cv2.DrawMatchesFlags_NOT_DRAW_SINGLE_POINTS)
cv2.imshow('Good matches(inliers)', img_good_matches_inliers)
cv2.waitKey(0)

f:id:nokixa:20211017025759p:plain

しっかり合っているものだけ残っているようです。 すごい!

ちなみに、マッチング点が1169点あったのに対して、inlierは866点になっていました。

>>> len(pts1_inliers)
866

Epipolar line描画

では計算したFundamental matrixを元にEpipolar lineを描画するとどうなるかというと、

def drawlines(img1,img2,lines,pts1,pts2):
    ''' img1 - image on which we draw the epilines for the points in img2
        lines - corresponding epilines '''
    r,c,ch = img1.shape
    for r,pt1,pt2 in zip(lines,pts1,pts2):
        color = tuple(np.random.randint(0,255,3).tolist())
        x0,y0 = map(int, [0, -r[2]/r[1] ])
        x1,y1 = map(int, [c, -(r[2]+r[0]*c)/r[1] ])
        img1 = cv2.line(img1, (x0,y0), (x1,y1), color,1)
        img1 = cv2.circle(img1,tuple(pt1),5,color,-1)
        img2 = cv2.circle(img2,tuple(pt2),5,color,-1)
    return img1,img2

pts1_inliers_200 = pts1_inliers[0:200,:]
pts2_inliers_200 = pts2_inliers[0:200,:]

lines1 = cv2.computeCorrespondEpilines(pts2_inliers_200.reshape(-1,1,2),2,F)
lines1 = lines1.reshape(-1,3)
img5,img6 = drawlines(img1.copy(),img2.copy(),lines1,pts1,pts2)
lines2 = cv2.computeCorrespondEpilines(pts1_inliers_200.reshape(-1,1,2),1,F)
lines2 = lines2.reshape(-1,3)
img3,img4 = drawlines(img2.copy(),img1.copy(),lines2,pts2,pts1)
img35 = np.hstack((img3,img5))
cv2.imshow('Epipolar lines', img35)
cv2.waitKey(0)

f:id:nokixa:20211021013954p:plain

こんな感じになっちゃいました。思っていたのと違う…

試しにFundamental matrixの計算で使うマッチング点数を減らしてみたりしましたが、

F, mask = cv2.findFundamentalMat(pts1[0:200,:], pts2[0:200,:], cv2.FM_LMEDS)

pts1_inliers = pts1[mask.ravel() == 1]
pts2_inliers = pts2[mask.ravel() == 1]

pts1_inliers_200 = pts1_inliers[0:200,:]
pts2_inliers_200 = pts2_inliers[0:200,:]

lines1 = cv2.computeCorrespondEpilines(pts2_inliers_200.reshape(-1,1,2),2,F)
lines1 = lines1.reshape(-1,3)
img5,img6 = drawlines(img1.copy(),img2.copy(),lines1,pts1,pts2)
lines2 = cv2.computeCorrespondEpilines(pts1_inliers_200.reshape(-1,1,2),1,F)
lines2 = lines2.reshape(-1,3)
img3,img4 = drawlines(img2.copy(),img1.copy(),lines2,pts2,pts1)
img35 = np.hstack((img3,img5))
cv2.imshow('Epipolar lines', img35)
cv2.waitKey(0)

f:id:nokixa:20211021015009p:plain

これまた変な結果に。

チュートリアルの最後のほうをよく見てみると、特徴点が1つの平面上にあると正しくFundamental matrixを計算できないという記載がありました。

Fundamental matrixの計算について

考察してみます。
以下のサイトのpdfを参考にしています。ケンブリッジ大学出版の教科書です。

Computer Vision Models

Amazonで購入も可能なようです。

Amazon | Computer Vision: Models, Learning, and Inference | Prince, Simon J. D. | Machine Vision

Fundamental matrixについて

Fundamental matrixがどういうものかというと、途中経過は端折りますが、

任意の実空間座標\boldsymbol{w}と、2視点からの画像上でこれに対応する点\boldsymbol{x_1}, \boldsymbol{x_2} (同次座標で)について、以下の式が成り立つ

\boldsymbol{x_1} \boldsymbol{F} \boldsymbol{x_2} = 0 \tag{1}

というものです。

\boldsymbol{x_1} = \lbrack x_1 y_1 1 \rbrack ^T , \boldsymbol{x_2} = \lbrack x_2 y_2 1 \rbrack ^T , \boldsymbol{F} = \begin{bmatrix} f_{11} & f_{12} & f_{13} \\
f_{21} & f_{22} & f_{23} \\
f_{31} & f_{32} & f_{33} \end{bmatrix}

というように要素を書き下すと、式(1)は

x_2 x_1 f_{11} + x_2 y_1 f_{12} + x_2 f_{13} + y_2 x_1 f_{21} + y_2 y_1 f_{f22} + y_2 f_{23} + x_1 f_{31} + y_1 f_{32} + f_{33} = 0 \tag{2} \lbrack x_2 x_1 \: x_2 y_1 \: x_2 \: y_2 x_1 \: y_2 y_1 \: y_2 \: x_1 \: y_1 \: 1 \rbrack \boldsymbol{f} = 0 \tag{3} (\boldsymbol{f} = \lbrack f_{11} \: f_{12} \: f_{13} \: f_{21} \: f_{22} \: f_{23} \: f_{31} \: f_{32} \: f_{33} \rbrack ^T)

と書けます。
\boldsymbol{F}はスケールの任意性があるので、自由度が8となっていて、8個以上のマッチングペアがあれば一意に決まります。

ここで、マッチングペアの点に対する実空間座標が全て同じ平面上にあった場合を考えます。

実空間座標を\boldsymbol{w}=\lbrack u \; v \; w \rbrack ^T とすると、画像上の点(同次座標)\boldsymbol{x_1}, \boldsymbol{x_2}u, v, wを使って以下のように表されます。

\lambda_1 \boldsymbol{x_1} = \boldsymbol{\Omega_1} \boldsymbol{w} + \boldsymbol{\tau_1}
\lambda_2 \boldsymbol{x_2} = \boldsymbol{\Omega_2} \boldsymbol{w} + \boldsymbol{\tau_2}

要素ごとに書くと、

\lambda_1 x_1 = \Omega_{1_{11}} u + \Omega_{1_{12}} v + \Omega_{1_{13}} w + \tau_{1_1}
\lambda_1 y_1 = \Omega_{1_{21}} u + \Omega_{1_{22}} v + \Omega_{1_{23}} w + \tau_{1_2}
\lambda_1 = \Omega_{1_{31}} u + \Omega_{1_{32}} v + \Omega_{1_{33}} w + \tau_{1_3}
\lambda_2 x_2 = \Omega_{2_{11}} u + \Omega_{2_{12}} v + \Omega_{2_{13}} w + \tau_{2_1}
\lambda_2 y_2 = \Omega_{2_{21}} u + \Omega_{2_{22}} v + \Omega_{2_{23}} w + \tau_{2_2}
\lambda_2 = \Omega_{2_{31}} u + \Omega_{2_{32}} v + \Omega_{2_{33}} w + \tau_{2_3}

式(2)に\lambda_1 \lambda_2を掛けると、全ての項が\lambda_1 x_1, \lambda_2 x_2, \lambda_1, \lambda_2の積で表されるので、上の6式を代入してみると、式(2)における係数がいずれも u, v, w, uv, uw, vw, u^2, v^2, w^2, 1 の線型結合で表されることがわかります。 なので、式(3)について8個の線形独立な係数ベクトルを用意することができ、\boldsymbol{f}を一意に(定数倍の自由度を除いて)決めることができます。

ただし、u, v, wの間に相関がある、例えばすべての実空間座標が1つの平面上にあると、 係数ベクトルがより少ない要素の線形結合になり、8個の線形独立なものを用意できなくなります。そうすると、一意に\boldsymbol{f}を定めることができなくなると。

なんとなくうまくまとまらないけどそんな感じ。

ここまで

今回の内容はここまでにします。
本当はちゃんとFundamental matrixが求められれば3次元座標の推定もできるのですが…

参考

tex書き方参考。特にはてなブログでの書き方が特殊だったので。

はてなブログでTex記法を使って行列を書く時の注意点 - 医療職からデータサイエンティストへ

hatena blogでtexによる下付き文字がうまく表示されない - MEMOcho-

LaTeXコマンド - 空白 - 水平方向のスペース - quad, hspace