OpenCV 鳥瞰図の座標がわかりません

Question

OpenCVの鳥瞰図用のホモグラフィを計算するために渡す点の座標がわかりません。
参考にしたのは「詳解OpenCV3」という本のP672の例19-1のプログラムです。

##問題点
本に書いてあることをほぼそのまま写して実行しました。
元の画像は次の画像です。(解像度を少し下げました)

鳥瞰図の結果は次のようになりました。

少し見づらいですが、左上に元の画像の鳥瞰画像が出ています。
コードは次に示します。

#include <iostream>
#include <fstream>
#include <opencv2/opencv.hpp>
#include <stdlib.h>

const int WIDTH = 600;
const int HEIGHT = 600;

//入力は[board_w] = 6, [board_h] = 9, [intinsics.xml] = intrinsics.xml, [checker_image] = image.png
void help(char *argv[]){
    std::cout << "Call: " /Users/shimizumiwa/coyomi/apps/bird/main.cpp<< argv[0] << "[board_w] [board_h] [intrinsics.xml] [checker_image]" << std::endl;
    //std::cout << "Demonstrates Pyramid Lucas-Kanade optical flow." << std::endl;
}

//引数：[board_w], [board_h], [intrinsics.xml], [checker_image]
//[board_w]：チェスボードの横のマス数
//[board_h]：チェスボードの盾のマス数
//[intrinsics.xml]：カメラ定数とかが入ってるファイル
//[checker_image]：入力画像
int main(int argc, char *argv[]){
    
    //argcはプログラムの引数の数
    //./bird image.jpgで実行したとしたらargv[0] = bird,argv[1] = image.jpgでargcは2
    if(argc != 5){
        std::cout << "\nERROR: too few parameters\n";
		//help(argv);
		return -1;
	}

	//atoiは文字列を数値に変換
	int board_w = atoi(argv[1]);
	int board_h = atoi(argv[2]);
	int board_n = board_w * board_h;
    cv::Size board_sz = cv::Size2i(board_w, board_h);
	//cv::Size board_sz(board_w, board_h);
	//読み込みでファイルを開く
	cv::FileStorage fs(argv[3], cv::FileStorage::READ);
	cv::Mat intrinsic, distortion;

	fs["camera_matrix"] >> intrinsic;
	fs["distortion_coefficients"] >> distortion;
    if(!fs.isOpened() || intrinsic.empty() || distortion.empty()){
        std::cout << "Error: Couldn't load intrinsic parameters from "
             << argv[3] << std::endl;
        return -1;
    }
	fs.release();

	cv::Mat gray_image, image, image0 = cv::imread(argv[4], 1);
	if(image0.empty()){
		std::cout << "Error: Couldn't load image" << argv[4] << std::endl;
		return -1;
	}

	//画像の歪みを補正する
	cv::undistort(image0, image, intrinsic, distortion, intrinsic);
	//入力画像、出力画像、グレーにする定数
	cv::cvtColor(image0, gray_image, cv::COLOR_BGR2GRAY);

	//その平面上のチェスボードを得る
	std::vector<cv::Point2f> corners;
	//チェスボードの内側コーナー位置を求める
	bool found = cv::findChessboardCorners(
		image,
        //内側コーナーの個数
		board_sz,
        //コーナーの出力配列
		corners/*,
        //白黒画像に変換する際に固定閾値の代わりに適応的閾値を利用
		//輪郭抽出の際の誤った四角を除外するための基準追加(あると逆に上手くいかない)
		cv::CALIB_CB_ADAPTIVE_THRESH | cv::CALIB_CB_FILTER_QUADS*/);
	if(!found){
		std::cout << "Couldn't acquire checkerboard on " << argv[4]
		<< ", only found " << corners.size() << " of " << board_n
		<< " corners\n";
		return -1;
	}
	//これらのコーナーをサブピクセルの精度で得る
    std::vector<std::vector<cv::Point2f>> img_points;
    cv::Mat src_gray = cv::Mat(image.size(), CV_8UC1);
    cv::cvtColor(image, src_gray, cv::COLOR_BGR2GRAY);
    cv::find4QuadCornerSubpix(src_gray, corners, cv::Size(3,3));
    cv::drawChessboardCorners(image, board_sz, corners, found);
    img_points.emplace_back(corners);
    

	//画像と物体の点を得る
	//物体の点(r,c)
	//(0, 0), (board_w-1, 0), (0, board_h-1), (board_w-1, board_h-1)
	//これはコーナーがcorners[r*board_w + c]であることを意味するらしい

	cv::Point2f objPts[4], imgPts[4];
	objPts[0].x = 0;                objPts[0].y = 0;
	objPts[1].x = board_w - 1;        objPts[1].y = 0;
	objPts[2].x = 0;                objPts[2].y = board_h - 1;
	objPts[3].x = board_w - 1;        objPts[3].y = board_h - 1;
	imgPts[0] = corners[0];
	imgPts[1] = corners[board_w - 1];
	imgPts[2] = corners[(board_h - 1) * board_w];
	imgPts[3] = corners[(board_h - 1) * board_w + board_w - 1];
    
    
	//点をB,G,R,YELLOWの順で描画する
	cv::circle(image, imgPts[0], 9, cv::Scalar(255, 0, 0), 3);
	cv::circle(image, imgPts[1], 9, cv::Scalar(0, 255, 0), 3);
	cv::circle(image, imgPts[2], 9, cv::Scalar(0, 0, 255), 3);
	cv::circle(image, imgPts[3], 9, cv::Scalar(0, 255, 255), 3);

	//発見したチェスボードを描画する
	cv::drawChessboardCorners(image, board_sz, corners, found);
	cv::imshow("Checkers", image);
    cv::imwrite("./apps/bird/Checkers.png", image);

	//ホモグラフィを求める

	cv::Mat H = cv::getPerspectiveTransform(objPts, imgPts);

	//このビューの高さzをユーザーが調整できるようにする
	double Z = 1;
	cv::Mat birds_image = cv::Mat::zeros(HEIGHT + 100, WIDTH + 100, CV_8UC3);
    //Escキーで止まる
	for(;;){
		H.at<double>(2, 2) = Z;
		//ホモグラフィを用いて再マッピングする

		cv::warpPerspective(
			image,
			birds_image,
            //変換配列
			H,
            //出力画像の配列
			birds_image.size(),
			cv::WARP_INVERSE_MAP | cv::INTER_LINEAR,
			cv::BORDER_CONSTANT,
            //境界を黒で埋める
			cv::Scalar::all(0)
		);

		cv::imshow("Birds_Eye", birds_image);
		int key = cv::waitKey() & 255;
		if(key == 'u') Z += 0.5;
		if(key == 'd') Z -= 0.5;
		if(key == 27) break;
	}

    cv::imwrite("./apps/bird/Birds_Eye.png", birds_image);
   // cv::moveWindow("./apps/bird/Birds_Eye.png", 3000, 3000);
    
	//回転と平行移動ベクトルを表示
	std::vector<cv::Point2f> image_points;
	std::vector<cv::Point3f> object_points;
	for(int i = 0; i < 4; ++i){
		image_points.push_back(imgPts[i]);
		object_points.push_back(
			cv::Point3f(objPts[i].x, objPts[i].y, 0)
		);
	}

	cv::Mat rvec, tvec, rmat;
	cv::solvePnP(
        //オブジェクト座標系の3次元座標系
		object_points,
        //画像座標系の2次元座標点
		image_points,
        //カメラ行列
		intrinsic,
        //最初に歪みを補正したので歪み係数はゼロになる
		cv::Mat(),
        //回転ベクトルを出力
		rvec,
        //平行移動ベクトルを出力
		tvec
	);

	cv::Rodrigues(rvec, rmat);

	//出力と終了
	std::cout << "rotation matrix: "  << rmat  << std::endl;
	std::cout << "translation vector: "  << tvec  << std::endl;
	std::cout << "homography matrix: "  << H   << std::endl;
	std::cout << "inverted homography matrix:  "  << H.inv()  << std::endl;

	return 1;

}

自分のやりたいのは、この元の画像を回転せずにそのままの向きで鳥瞰図にすることです。
また、本に載っていた鳥瞰図取得の結果では、ウィンドウの下の真ん中あたりからから上に画像が広がっていくような鳥瞰図になっていました。

##特に知りたいこと
ホモグラフィを求めるための座標(objPts,imgPts)が間違っているのではないかと考えています。ですが、本のサンプルコードに書いてあるようにその座標がどうしてそう求められるのか理解していません。どのような考え方でこのobjPtsとimgPtsがこの求め方で出るのか教えていただきたいです。

##お陰様で出来ました

Answer 1

cv::Mat H = cv::getPerspectiveTransform(objPts, imgPts);

これで何を求めているのかを把握していますか？

• 引数の順序はこれで正しいですか？

(→後段で WARP_INVERSE_MAP を使っているから，これは正しい模様）

• objPts と imgPts の４点の対応関係は正しいですか？

→像が回転する要因になっていないか？

• objPts の値は実行時のコマンドライン引数由来であるようですが，一体どんな値を与えているのですか？

→コマンドライン引数では「何を」与える仕様としているのですか？
コマンドライン引数の使われ方を見るとちょっと謎です（コーナーの点の個数として使われているが，objPtsの座標値にも用いられている．これは正しいの？）

---

参考までに，最も単純なコードを示す．

• 入力画像はサイズがでかかったので適当に縮小した
• チェスボードパターンの四隅座標は，当該縮小した画像上で手作業（目視）で求めたものをハードコーディングしてある
• warpPerspective で WARP_INVERSE_MAP を使ってない

int main()
{
	//※"Board.png"は質問に貼られていた画像を縦横半分に縮小したものを使用した
	cv::Mat SrcImg = cv::imread( "Board.png" );
	if( SrcImg.empty() )return 0;
	
	//SrcImg 上のチェスボードパターンの四隅の座標
	cv::Point2f ImgPts[4] = {
		{ 218, 267 }, { 385, 267 }, { 196, 375 }, { 406, 373 }
	};

	//↑の４点の変換結果目標座標．
	//320*240サイズの結果画像の四隅から少し(10[pixel])だけ内側．
	cv::Mat ResultImg = cv::Mat::zeros( 240, 320, CV_8UC3 );
	cv::Point2f ObjPts[4] = {
		{ 10,10 }, { 310, 10, }, { 10, 230 }, { 310, 230 }
	};

	cv::Mat H = cv::getPerspectiveTransform( ImgPts, ObjPts );
	cv::warpPerspective(
		SrcImg, ResultImg,
		H, ResultImg.size(),
		cv::INTER_LINEAR, cv::BORDER_CONSTANT,
		cv::Scalar::all(0)
	);
	cv::imshow( "Result", ResultImg );
	cv::waitKey();
	return 0;
}

結果：

Comments

[fana]

簡単なサンプルを示した．

質問コードでは，objPts　の座標値に「パターン個数」の値をそのまま使っているから１桁の値になっており，結果として，変換結果画像でチェスボードの部分が数pixelな大きさになっている．

この点は後段での「Z の値がどうの」で倍率調整する想定なのだと見えるが，その調整方法で大丈夫なのか？　っていう．
（Zを 0.01 とかそういう値にできる必要があると思うのだが，今の質問コードだとそれができないのでは．
→Zの初期値を 100 とか 0.01 とか，いろいろ変えて何が起こるのかを見てみると良いのではないかと．
）

また，変換結果の像の位置（オフセット）はまともに考えられていない様子に見える．
（常にチェスパターンの４隅のうちのどこかが結果画像の(0,0)に来てしまうであろう）

[fana]

getPerspectiveTransform()でホモグラフィ行列を求めるとき
imgPts は入力画像上での位置(pixel)であり，
objPts はそれが出力画像上のどこに行くべきかという位置(pixel)である．

（ホモグラフィ行列を求む，というそのものズバリな名前の findHomography() も存在する）

[fana]

本のサンプルがどんなものかは知らないけど，
ここに示したように，最もシンプルな状態（その関数を使ってみるのに必要最低限な要素しか無い状態）を最初はやってみると良いんじゃないかな．

（何故か「サンプル」って　中途半端にアプリケーションっぽく　なってたりするけど，そういうの邪魔だよね）

[koomint]

ご回答ありがとうございます。

実は私も、objPtsでチェッカーボードのマスの数がそのまま用いられているのに疑問を持ちました。このobjPtsとimgPtsの座標の代入値は本に書いてあったものをそのまま写したので、自分もなぜこの計算方法で出すのかわかりません。間違えているのでは？と思いました。

Zの値については0.01とか小さい値になる程画像が左上に向かって小さくなって、0では全面真っ黒になることを確認しています。100のように大きくなるほど、画像が左上から右下まで拡大されていくようになります。

私が知りたいのは、なぜこのobjPtsとimgPtsの座標計算方法で本では上手くいってるのか、計算方法の理屈をもしご存知であればお伺いしたかったのです。
また、一つの画像に対してだけ鳥瞰図を得たいわけではなく、いろんな画像が入ってきたときに全てに対してプログラム一つで鳥瞰図をそれぞれ出したいので、objPtsやimgPtsの座標情報をcornersなどを使って固定値ではなくて変数として出したいです。

[koomint]

一応、cornersなど難しい要素を入れないで座標の値を自分で入れてやったものはうまくいったので、やっぱり座標の求め方が問題なんだなというのは把握しました。

[fana]

objPts の座標として「コーナー点の個数-1」を用いるというのは，
要は，チェスボードパターンの１個の四角形のサイズを 1x1 とする世界を考えていることに相当します．
（単位を「ボードの四角形の辺の長さ」という形での正規化を考えた世界，という感じ）
おそらくそういう意図があるのでしょう．
（その本に，そういう意味合いに関して何も説明が書いてないんですかね？）

↓

しかし，OpenCVの関数 warpPerspective は出力側の単位として[pixel]が用いられたものを期待しているので，
そのような正規化が成された世界に関して計算された行列 H をそのまま使っても，座標の単位に齟齬があるため，既述したように，めちゃくちゃ小さな結果画像が得られるわけです．
（「1x1」の「1」が，1[pixel]として扱われるから）

で，結果画像のズーム具合を Z なる値（H[2][2]の値）で後から調整しているわけですね．
（同次座標で計算した第三の要素の値を増減させることによって最終結果のスケールを変えている……のだと思うけど，もしもそういう意味なのであれば本来は H の第３行目の成分全てを同じ倍率でスケーリングするべきだと思う）

[fana]

> いろんな画像が入ってきたときに全てに対してプログラム一つで鳥瞰図をそれぞれ出したい

結局のところ，あなたが「最終的な結果をどういう絵にしたいか？」に依存するでしょう．

* チェスボードが左上らへんに写っている入力画像
* チェスボードが右下らへんに写っている入力画像
* …

などが存在し得る中で，チェスボードのコーナー点を用いて演算を行えば，objPtsの側が同じであれば（その画像でのチェスボードのコーナー点の位置がいつも結果画像でobjPtsの場所にくるので）「結果画像内で，入力画像の絵が占める部分」が異なる形になるでしょう．

[koomint]

試しに、変換前の点の座標にはcornersを用いた本に書いてある通りの座標情報を、変換後の点の座標にはfanaさんのように適当に表示したい場所の座標を入れてみました。
また、fanaさんの送ってくださったコードを見て、勘違いしていたところ(cv::warpPerspectiveのcv::WARP_INVERSE_MAPがいらなかった)が気づけたのでそこも変えてみたところ、下のようにうまくいきました。
(すいません画像が載せられないので質問の一番最後のとこに載せます)
全部が表示されているわけではないのでそこは調整が必要ですが、重要なとこは解決しました。

[koomint]

詳しく教えていただいてありがとうございました。
今回はZの変更についてはしないつもりではありましたが、あっても良いと思うので、スケーリングについてはこの後考えてみたいと思います。Hの三行目の成分を全て定数倍するなどはまだHの中身が理解しきれていないので調べていきたいと思っています。
本当にありがとうございました。助かりました。

[fana]

OpenCVのリファレンスの getPerspectiveTransform のところに H に関する数式が書いてあるハズです．
入力画像の座標(x,y) に対する，Hを用いた出力画像の座標が，こんな感じの式で書いてあるハズ（一部，記号が違うかもだけど）：

(k*x' k*y' k)^t = H * (x y 1)^t　　（※ここで ^t は転置の意味）

左辺の第三成分で他の成分を割ることで，結果画像上の座標(x', y')が得られる仕組み．（割るのが「同次座標」の約束事．詳しくは「同次座標」で調べてください）
で，ここでもしもHの三行目の成分を全部２倍に改ざんしたら，左辺の第三成分だけが2倍になる→なのでそれで割った結果の座標は元の半分になる．

[fana]

> 全部が表示されているわけではない

結果を見れば想像が付くでしょうが，「全部」を表示するには相当に巨大な結果画像になり得ます．
（そして像がめちゃくちゃに伸びます．）
実用上は「ある程度の範囲」で良いのではないかと．
（例えば，入力画像の四隅の座標値とかを↑の式に突っ込んで計算してみれば結果画像のどこにくるのかがわかるでしょう．）

[koomint]

そんな感じの式、見た気がします！そういうことだったんですね。
わざわざありがとうございます。
そうですね、色々試してみます。