OpenCV中對的旋轉一些思考
目錄
1.問題描述
數字影象處理是一門很有意思的學問,在現實生活中往往一個很簡單的問題在數字影象中有時會非常複雜,旋轉便是一類非常有意思的問題。如何在離散影象中高精度、快速求解影象的旋轉角度,這個問題我思考了很長時間,下面會使用三種不同的演算法逐一計算輪廓的旋轉。
如圖1-1所示,這是一幅鳥的輪廓:
如果這副影象發生了旋轉,會出現什麼樣的情況呢?如圖1-2所示,綠色的線表示旋轉後的輪廓,在這裡我設定了旋轉角度為50°,旋轉中心為輪廓的形心。
這也就是說,給定白色和綠色的兩條輪廓,我們需要解出它們的旋轉角度(50°);
2.旋轉的三種解法
2.1 應用迭代法進行求解
演算法的基本引數:
1.給定輪廓相似程度的度量方式,這裡我採用了ShapeContextDistance;
2.給定旋轉的步長,也就是每次匹配時綠色輪廓旋轉的變化量,這裡為了實驗方便步長為10°;
3.給定迭代的終止條件,這裡為了快速迭代,迭代次數為10;
如何迭代:
每次旋轉後計算ShapeContextDistance,達到閾值或者達到迭代上限即跳出迴圈。
code:
#include<opencv2/opencv.hpp> #include<iostream> using namespace std; using namespace cv; //計算影象上所有的輪廓,並計算它的質心和對角線長度 void FindBlobs(Mat img, vector<vector<Point>> &contours, vector<Point2f> &MassCentre, vector<float>&DiagonalLength) { //判斷影象是否為8位單通道 if (img.empty() || img.depth() != CV_8UC1) { cout << "Invalid Input Image!"; exit(-1); } //計算輪廓 vector<Vec4i> hierarchy; findContours(img, contours, hierarchy, RETR_EXTERNAL, CHAIN_APPROX_SIMPLE); //計算輪廓的最小外接矩形 vector<Rect> boundRect(contours.size()); vector<Moments> mu(contours.size()); for (size_t i = 0; i < contours.size(); i++) { //計算外界矩形 boundRect[i] = boundingRect(Mat(contours[i])); //計算輪廓矩 mu[i] = moments(contours[i], false); //當binaryImage=true時,所有的非零值都視為1 //計算質心 MassCentre.push_back(Point2f(static_cast<float>(mu[i].m10 / mu[i].m00), static_cast<float>(mu[i].m01 / mu[i].m00))); //計算外接矩形的對角距離 DiagonalLength.push_back((float)norm(boundRect[i].tl() - boundRect[i].br())); } } //對輪廓進行平移變化 void TransformContour(vector<vector<Point>> contours, vector<vector<Point>> &contours_Trans, Vec2f Translate) { //判斷 if (contours.size() == 0) { cout << "Invalid input!"; exit(-1); } //平移變換 contours_Trans = contours; for (size_t i = 0; i < contours.size(); i++) { for (int idx = 0; idx < contours[i].size(); idx++) { contours_Trans[i][idx] = Point(contours[i][idx].x + Translate(0), contours[i][idx].y + Translate(1)); } } } //旋轉輪廓 //1.基於仿射變化的2階方陣M計算旋轉位置 Point RotatePoint(const Mat &M, const Point &p) { Point2f rp; rp.x = (float)(M.at<double>(0, 0)*p.x + M.at<double>(0, 1)*p.y + M.at<double>(0, 2)); rp.y = (float)(M.at<double>(1, 0)*p.x + M.at<double>(1, 1)*p.y + M.at<double>(1, 2)); return rp; } //2.計算選擇輪廓 void RotateContour(vector<vector<Point>> contours, vector<vector<Point>> &contours_Rotated, double Angle, Point2f Centre) { //判斷 if (contours.size() == 0) { cout << "Invalid input!"; exit(-1); } //計算仿射矩陣 Mat M= getRotationMatrix2D(Centre, Angle, 1.0); //計算旋轉輪廓 contours_Rotated = contours; for (size_t i = 0; i < contours.size(); i++) { for (int idx = 0; idx < contours[i].size(); idx++) { contours_Rotated[i][idx] = RotatePoint(M, contours[i][idx]); } } } //對輪廓進行簡單的取樣使得輪廓點數量為300 方便計算ShapeContextDistance static vector<Point> simpleContour(vector<vector<Point>> _contoursQuery, int n = 300) { //當前資料點 vector <Point> contoursQuery; for (size_t border = 0; border<_contoursQuery.size(); border++) { for (size_t p = 0; p<_contoursQuery[border].size(); p++) { contoursQuery.push_back(_contoursQuery[border][p]); } } //增補資料點至n=300 int dummy = 0; for (int add = (int)contoursQuery.size() - 1; add<n; add++) { contoursQuery.push_back(contoursQuery[dummy++]); } //將資料點均勻化 random_shuffle(contoursQuery.begin(), contoursQuery.end()); vector<Point> cont; for (int i = 0; i<n; i++) { cont.push_back(contoursQuery[i]); } return cont; } int main() { //1 以GRAY的形式讀入,並二值化 Mat src = imread("rotate_consider.png",0); if (src.empty()) { cout << "Invalid Input Image!"; exit(-1); } threshold(src,src,50,255,CV_THRESH_BINARY); //2 計算輪廓及其質心、對角距離 vector<vector<Point>> contours; vector<Point2f> MassCentre; vector<float>DiagonalLength; FindBlobs(src, contours, MassCentre, DiagonalLength); //3 對輪廓進行平移變換,構造影象空間mContourSpace, 並將輪廓的質心座標移動到新座標系下的ptCCentre位置 vector<vector<Point>> contours_Trans(contours.size()); Mat mContourSpace(Size(DiagonalLength[0], DiagonalLength[0]), CV_8UC3, Scalar(0)); Point2f ptCCentre(DiagonalLength[0] / 2, DiagonalLength[0] / 2); Vec2f Translation(ptCCentre.x - MassCentre[0].x, ptCCentre.y - MassCentre[0].y); TransformContour(contours, contours_Trans, Translation); drawContours(mContourSpace, contours_Trans, 0, Scalar(255, 255, 255), 2, 8); //4 對輪廓進行旋轉變換 vector<vector<Point>> contours_Rotated(contours.size()); double Angle = -50.0; RotateContour(contours_Trans, contours_Rotated, Angle, ptCCentre); Mat _mContourSpace = mContourSpace.clone(); drawContours(_mContourSpace, contours_Rotated, 0, Scalar(0, 255, 0), 2, 8); //5 基於迭代思想來計算角度 //5.1 建立shapecontextdistance物件指標 Ptr <ShapeContextDistanceExtractor> mysc = cv::createShapeContextDistanceExtractor(); float bestMatchAngle = 0; float bestDis = FLT_MAX; float test_angle = 0; vector<Point> QueryContour= simpleContour(contours_Trans); //10次迭代,每次運算后角度增加10° for (int i = 0; i < 10; i++) { Mat mContourMatch = mContourSpace.clone(); test_angle += 10; vector<vector<Point>> tem_contours; RotateContour(contours_Rotated, tem_contours, test_angle, ptCCentre); vector<Point> TestContour = simpleContour(tem_contours); float dis = mysc->computeDistance(QueryContour, TestContour); drawContours(mContourMatch, tem_contours, (int)0, Scalar(0, 255, 0), 2, 8); putText(mContourMatch, format("Distance: %f ", dis), Point2f(ptCCentre.x + 30, 20), CV_FONT_HERSHEY_COMPLEX, 0.5, Scalar(0, 0, 255), 0.5); if (dis<bestDis) { bestMatchAngle = test_angle; bestDis = dis; } } //繪製計算結果 Mat dst= mContourSpace.clone(); vector<vector<Point>> result_rotate; RotateContour(contours_Rotated, result_rotate, bestMatchAngle, ptCCentre); drawContours(dst, result_rotate, (int)0, Scalar(0, 0, 255),2,8); return 0; }
演算法的運算結果:
圖 初始影象
圖 演算法的第1-3次迭代,計算的shapecontext距離值分別為:25.7、14.8、1.23
圖 演算法的第4-6次迭代,計算的shapecontext距離值分別為0.095、0.068、0.076
圖 演算法的第7-10次迭代,計算的shapecontext距離值分別為1.45、7.9、28.8、30.2
圖 演算法返回的最終解,angle=50°
2.2 應用特徵橢圓進行求解
特徵橢圓計算旋轉角度方法的數學理論推導詳見 Peter Corke著作的《Robotics Vision and Control Page:351-353》;在此只給出簡單公式圖片:
演算法的實現:
#include<opencv2/opencv.hpp>
#include<iostream>
using namespace std;
using namespace cv;
//尋找最大輪廓
vector<Point> FindBigestContour(Mat &src) {
int imax = 0;
int imaxcontour = -1;
std::vector<std::vector<Point> >contours;
findContours(src, contours, CV_RETR_LIST, CV_CHAIN_APPROX_SIMPLE);
for (int i = 0; i<contours.size(); i++) {
int itmp = contourArea(contours[i]);
if (imaxcontour < itmp) {
imax = i;
imaxcontour = itmp;
}
}
return contours[imax];
}
//輪廓的特徵橢圓的特徵矩陣
void computeEllipse(const vector<Point> &contour,Point ¢er,Size &EllipseSize,
float &angle)
{
//計算矩
Moments mc = moments(contour,false);
//計算形心
center = Point(mc.m10/mc.m00,mc.m01/mc.m00);
//計算特徵矩陣
Mat J(2, 2, CV_32F);
J.at<float>(0, 0) = mc.mu20;
J.at<float>(0, 1) = mc.mu11;
J.at<float>(1, 0) = mc.mu11;
J.at<float>(1, 1) = mc.mu02;
//計算J的特徵向量和特徵值
Mat eigenValues;
Mat eigenVectors;
eigen(J, eigenValues, eigenVectors);
//計算橢圓的尺寸
Mat size;
sqrt(4*eigenValues / mc.m00,size);
int a = size.at<float>(0);
int b = size.at<float>(1);
EllipseSize = Size(a, b);
//計算橢圓長軸的角度
angle = 180*(atan2(eigenVectors.at<float>(0,1), eigenVectors.at<float>(0,
0)))/3.14;
}
int main()
{
//1.準備實驗的影象素材
Mat src=imread("rotate_consider_src.png",0);
Mat match = imread("rotate_consider_match.png", 0);
threshold(src,src,50,255,CV_THRESH_BINARY);
threshold(match, match, 50, 255, CV_THRESH_BINARY);
//2.儲存輪廓
vector<Point> contour_src = FindBigestContour(src);
vector<Point> contour_match = FindBigestContour(match);
//3.計算特徵橢圓
Point center_src, center_match;
Size size_src,size_match;
float angle_src,angle_match;
computeEllipse(contour_src,center_src,size_src,angle_src);
computeEllipse(contour_match, center_match, size_match, angle_match);
//4.繪製
vector<vector<Point>> contours;
contours.push_back(contour_src);
contours.push_back(contour_match);
Mat drawMat(src.size(),CV_8UC3,Scalar(0));
drawContours(drawMat,contours,0,Scalar(255,0,0),2,8);
drawContours(drawMat, contours, 1, Scalar(255, 255, 255), 2, 8);
//src
ellipse(drawMat,center_src,size_src,angle_src,0,360,Scalar(0,0,255),2,8);
circle(drawMat,center_src,5,Scalar(0,0,255),-1);
//match
ellipse(drawMat, center_match, size_match, angle_match, 0, 360,
Scalar(0,255,0), 2, 8);
//單獨繪製src
Mat img_01(src.size(), CV_8UC3, Scalar(0));
drawContours(img_01, contours, 0, Scalar(255, 0, 0), 2, 8);
ellipse(img_01, center_src, size_src, angle_src, 0, 360, Scalar(0, 0, 255),
2, 8);
circle(img_01, center_src, 5, Scalar(0, 0, 255), -1);
//單獨繪製match
Mat img_02(src.size(), CV_8UC3, Scalar(0));
drawContours(img_02, contours, 1, Scalar(255, 0, 0), 2, 8);
ellipse(img_02, center_match, size_match, angle_match, 0, 360,
Scalar(0, 255, 0), 2, 8);
circle(img_02, center_src, 5, Scalar(0, 255,0), -1);
//5.角度計算結果
putText(drawMat, format("angleDiff: %f ", fabs(angle_src-angle_match))
, Point2f(600, 40),
CV_FONT_HERSHEY_COMPLEX, 0.5, Scalar(0, 255, 255), 0.5);
return 0;
}
演算法的計算結果:
圖 計算輪廓的特徵橢圓
圖 演算法返回的最終解,angle=50.027°
2.3 應用PCA主成分分析的方法求解
PCA是機器學習裡面進行資料降維的常用方法之一,不懂的小夥伴建議去讀相關文獻。
具體思路是計算兩個輪廓的主向量,並計算主向量之間的夾角。
實現code:
#include<opencv2/opencv.hpp>
#include<iostream>
using namespace std;
using namespace cv;
//對輪廓點進行PCA分析
void contourPCA(vector<Point> &contour,Mat &img, float &angle)
{
//1 重新佈置資料點
Mat data_pts = Mat(contour.size(), 2, CV_64FC1);
for (int i = 0; i < data_pts.rows; ++i)
{
data_pts.at<double>(i, 0) = contour[i].x;
data_pts.at<double>(i, 1) = contour[i].y;
}
//2 進行PCA分析
PCA pca_analysis(data_pts, Mat(), CV_PCA_DATA_AS_ROW);
//3 計算輪廓中心點
Point2f center = Point2f(pca_analysis.mean.at<double>(0, 0),pca_analysis.mean.at<double>(0, 1));
//4 生成 特徵向量矩陣 和 特徵值矩陣
Mat eigen_vecsMat = pca_analysis.eigenvectors;
Mat eigen_valMat = pca_analysis.eigenvalues;
//5 儲存特徵值和特徵向量
vector<Point2f> eigen_vecs(2); //儲存PCA分析結果,其中0組為主方向,1組為垂直方向
vector<float> eigen_val(2);
for (int i = 0; i < 2; ++i)
{
eigen_vecs[i] = Point2d(eigen_vecsMat.at<double>(i, 0), eigen_vecsMat.at<double>(i, 1));
eigen_val[i] = eigen_valMat.at<double>(i, 0);
}
//6 計算外界矩形的交點
//6.1 計算輪廓的最大外接矩形約束交點計算的範圍
Rect boundRect = boundingRect(contour);
//6.2 計算主方向與矩形輪廓的交點
//6.2.1 計算方向1與矩形的交點
float k1 = eigen_vecs[0].y / eigen_vecs[0].x;//斜率k1
angle = (atan2(eigen_vecs[0].y, eigen_vecs[0].x))*180/3.14;
Point2f pt1 = Point2f(boundRect.x, k1*(boundRect.x - center.x) + center.y);
Point2f pt2 = Point2f((boundRect.x + boundRect.width), k1*((boundRect.x + boundRect.width) - center.x) + center.y);
//6.2.2 計算方向2與矩形的交點
float k2 = eigen_vecs[1].y / eigen_vecs[1].x;//斜率k1
Point2f pt3 = Point2f(boundRect.x, k2*(boundRect.x - center.x) + center.y);
Point2f pt4 = Point2f((boundRect.x + boundRect.width), k2*((boundRect.x + boundRect.width) - center.x) + center.y);
//7 繪圖
line(img, pt1, pt2, Scalar(0, 255, 255),2,8);
line(img, pt3, pt4, Scalar(0, 255, 0),2,8);
circle(img, center, 4, Scalar(0, 0, 255), -1);
}
//尋找最大輪廓
vector<Point> FindBigestContour(Mat &src) {
int imax = 0;
int imaxcontour = -1;
std::vector<std::vector<Point> >contours;
findContours(src, contours, CV_RETR_LIST, CV_CHAIN_APPROX_SIMPLE);
for (int i = 0; i<contours.size(); i++) {
int itmp = contourArea(contours[i]);
if (imaxcontour < itmp) {
imax = i;
imaxcontour = itmp;
}
}
return contours[imax];
}
int main()
{
//1.準備實驗的影象素材
Mat src = imread("rotate_consider_src.png", 0);
Mat match = imread("rotate_consider_match.png", 0);
threshold(src, src, 50, 255, CV_THRESH_BINARY);
threshold(match, match, 50, 255, CV_THRESH_BINARY);
//2.儲存輪廓
vector<Point> contour_src = FindBigestContour(src);
vector<Point> contour_match = FindBigestContour(match);
vector<vector<Point>> contours(2);
contours[0]= contour_src;
contours[1] = contour_match;
//3.對輪廓進行主成分分析
Mat drawMat_src(src.size(), CV_8UC3, Scalar(0));
Mat drawMat_match(src.size(), CV_8UC3, Scalar(0));
drawContours(drawMat_src,contours,0,Scalar(255,255,255),2,8);
drawContours(drawMat_match, contours, 1, Scalar(255, 0, 0), 2, 8);
float angle_src,angle_match;
contourPCA(contour_src,drawMat_src,angle_src);
contourPCA(contour_match, drawMat_match, angle_match);
//5.角度計算結果
putText(drawMat_src, format("angleDiff: %f ", angle_src), Point2f(600, 40),
CV_FONT_HERSHEY_COMPLEX, 0.5, Scalar(255, 255, 255), 0.5);
putText(drawMat_match, format("angleDiff: %f ", angle_match), Point2f(600, 40),
CV_FONT_HERSHEY_COMPLEX, 0.5, Scalar(255, 0, 0), 0.5);
return 0;
}
演算法的計算結果(angle=48.52°):
圖 輪廓的主向量用綠色和黃色畫出,圖1的主向量角度angle1=31.04°,圖2的主向量角度angle2=79.56°
3. 結果分析
結合演算法的效率和精度,選擇方法2計算旋轉角度最佳,方法1和方法3可在特定的條件下應用。
4.相關資料:
1.answeOpenCV論壇:
http://answers.opencv.org/question/113492/orientation-of-two-contours/
http://answers.opencv.org/question/28489/how-to-compare-two-contours-translated-from-one-another/
http://answers.opencv.org/question/168357/way-to-filter-out-false-positives-in-template-matching/
http://answers.opencv.org/question/51486/template-matching-is-wrong-with-specific-reference-image/
2.GitHub:
https://github.com/Smorodov/LogPolarFFTTemplateMatcher/blob/master/fftm.cpp
3.部落格資料:
3.1 影象矩
https://blog.csdn.net/kuweicai/article/details/79027388
3.2 影象的平移、映象和旋轉
https://blog.csdn.net/qq_20823641/article/details/51925091
3.3 字元識別與區域定位
https://blog.csdn.net/u012556077/article/details/47126311
4.參考文獻
1. Book:Robotics Vision and Control Author:Peter Corke Page:351-353