2. 程式人生 > >OpenCV實現基於Zernike矩的亞畫素邊緣檢測

OpenCV實現基於Zernike矩的亞畫素邊緣檢測

阿新 • • 發佈：2019-02-04

    在做物體檢測時，由於成本和應用場合的限制，不能夠一味地增加相機的解析度，或者已經用了解析度很高的相機，但是視野範圍很大，仍然無法實現很高的精度，這時就要考慮亞畫素技術，亞畫素技術就是在兩個畫素點之間進行進一步的細分，從而得到亞畫素級別的邊緣點的座標（也就是float型別的座標），一般來說，現有的技術可以做到2細分、4細分，甚至很牛的能做到更高，通過亞畫素邊緣檢測技術的使用，可以節約成本，提高識別精度。

    常見的亞畫素技術包括灰度矩、Hu矩、空間矩、Zernike矩等，通過查閱多篇論文，綜合比較精度、時間和實現難度，選擇Zernike矩作為專案中的亞畫素邊緣檢測演算法。基於Zernike矩的邊緣檢測原理，下一篇文章詳細再總結一下，這裡大體把步驟說一下，並使用OpenCV實現。

    大體步驟就是，首先確定使用的模板大小，然後根據Zernike矩的公式求出各個模板係數，例如我選擇的是7*7模板（模板越大精度越高，但是運算時間越長），要計算M00，M11R，M11I，M20，M31R，M31I，M40七個Zernike模板。第二步是對待檢測影象進行預處理，包括濾波二值化等，也可以在進行一次Canny邊緣檢測。第三步是將預處理的影象跟7個Zernike矩模板分別進行卷積，得到7個Zernike矩。第四步是把上一步的矩乘上角度校正係數（這是因為利用Zernike的旋轉不變性，Zernike模型把邊緣都簡化成了x=n的直線，這裡要調整回來）。第五步是計算距離引數l和灰度差引數k，根據k和l判斷該點是否為邊緣點。以下是基於OpenCV的實現。

#include <math.h>
#include <opencv2/opencv.hpp>

using namespace std;
using namespace cv;

const double PI = 3.14159265358979323846;
const int g_N = 7;

Mat M00 = (Mat_<double>(7, 7) <<
	0, 0.0287, 0.0686, 0.0807, 0.0686, 0.0287, 0,
	0.0287, 0.0815, 0.0816, 0.0816, 0.0816, 0.0815, 0.0287,
	0.0686, 0.0816, 0.0816, 0.0816, 0.0816, 0.0816, 0.0686,
	0.0807, 0.0816, 0.0816, 0.0816, 0.0816, 0.0816, 0.0807,
	0.0686, 0.0816, 0.0816, 0.0816, 0.0816, 0.0816, 0.0686,
	0.0287, 0.0815, 0.0816, 0.0816, 0.0816, 0.0815, 0.0287,
	0, 0.0287, 0.0686, 0.0807, 0.0686, 0.0287, 0);

Mat M11R = (Mat_<double>(7, 7) <<
	0, -0.015, -0.019, 0, 0.019, 0.015, 0,
	-0.0224, -0.0466, -0.0233, 0, 0.0233, 0.0466, 0.0224,
	-0.0573, -0.0466, -0.0233, 0, 0.0233, 0.0466, 0.0573,
	-0.069, -0.0466, -0.0233, 0, 0.0233, 0.0466, 0.069,
	-0.0573, -0.0466, -0.0233, 0, 0.0233, 0.0466, 0.0573,
	-0.0224, -0.0466, -0.0233, 0, 0.0233, 0.0466, 0.0224,
	0, -0.015, -0.019, 0, 0.019, 0.015, 0);

Mat M11I = (Mat_<double>(7, 7) <<
	0, -0.0224, -0.0573, -0.069, -0.0573, -0.0224, 0,
	-0.015, -0.0466, -0.0466, -0.0466, -0.0466, -0.0466, -0.015,
	-0.019, -0.0233, -0.0233, -0.0233, -0.0233, -0.0233, -0.019,
	0, 0, 0, 0, 0, 0, 0,
	0.019, 0.0233, 0.0233, 0.0233, 0.0233, 0.0233, 0.019,
	0.015, 0.0466, 0.0466, 0.0466, 0.0466, 0.0466, 0.015,
	0, 0.0224, 0.0573, 0.069, 0.0573, 0.0224, 0);

Mat M20 = (Mat_<double>(7, 7) <<
	0, 0.0225, 0.0394, 0.0396, 0.0394, 0.0225, 0,
	0.0225, 0.0271, -0.0128, -0.0261, -0.0128, 0.0271, 0.0225,
	0.0394, -0.0128, -0.0528, -0.0661, -0.0528, -0.0128, 0.0394,
	0.0396, -0.0261, -0.0661, -0.0794, -0.0661, -0.0261, 0.0396,
	0.0394, -0.0128, -0.0528, -0.0661, -0.0528, -0.0128, 0.0394,
	0.0225, 0.0271, -0.0128, -0.0261, -0.0128, 0.0271, 0.0225,
	0, 0.0225, 0.0394, 0.0396, 0.0394, 0.0225, 0);

Mat M31R = (Mat_<double>(7, 7) <<
	0, -0.0103, -0.0073, 0, 0.0073, 0.0103, 0,
	-0.0153, -0.0018, 0.0162, 0, -0.0162, 0.0018, 0.0153,
	-0.0223, 0.0324, 0.0333, 0, -0.0333, -0.0324, 0.0223,
	-0.0190, 0.0438, 0.0390, 0, -0.0390, -0.0438, 0.0190,
	-0.0223, 0.0324, 0.0333, 0, -0.0333, -0.0324, 0.0223,
	-0.0153, -0.0018, 0.0162, 0, -0.0162, 0.0018, 0.0153,
	0, -0.0103, -0.0073, 0, 0.0073, 0.0103, 0);

Mat M31I = (Mat_<double>(7, 7) <<
	0, -0.0153, -0.0223, -0.019, -0.0223, -0.0153, 0,
	-0.0103, -0.0018, 0.0324, 0.0438, 0.0324, -0.0018, -0.0103,
	-0.0073, 0.0162, 0.0333, 0.039, 0.0333, 0.0162, -0.0073,
	0, 0, 0, 0, 0, 0, 0,
	0.0073, -0.0162, -0.0333, -0.039, -0.0333, -0.0162, 0.0073,
	0.0103, 0.0018, -0.0324, -0.0438, -0.0324, 0.0018, 0.0103,
	0, 0.0153, 0.0223, 0.0190, 0.0223, 0.0153, 0);

Mat M40 = (Mat_<double>(7, 7) <<
	0, 0.013, 0.0056, -0.0018, 0.0056, 0.013, 0,
	0.0130, -0.0186, -0.0323, -0.0239, -0.0323, -0.0186, 0.0130,
	0.0056, -0.0323, 0.0125, 0.0406, 0.0125, -0.0323, 0.0056,
	-0.0018, -0.0239, 0.0406, 0.0751, 0.0406, -0.0239, -0.0018,
	0.0056, -0.0323, 0.0125, 0.0406, 0.0125, -0.0323, 0.0056,
	0.0130, -0.0186, -0.0323, -0.0239, -0.0323, -0.0186, 0.0130,
	0, 0.013, 0.0056, -0.0018, 0.0056, 0.013, 0);



int main()
{

	Mat OriginalImage = imread("original.png", 0);
	Mat SubImage = OriginalImage;

	Mat NewSmoothImage;
	medianBlur(OriginalImage, NewSmoothImage, 13);

	Mat NewAdaThresImage;
	adaptiveThreshold(NewSmoothImage, NewAdaThresImage, 255, ADAPTIVE_THRESH_GAUSSIAN_C, THRESH_BINARY_INV, 7, 4);

	vector<Point2d> SubEdgePoints;

	Mat ZerImageM00;
	filter2D(NewAdaThresImage, ZerImageM00, CV_64F, M00, Point(-1, -1), 0, BORDER_DEFAULT);
	//////////filter2D( cannyImage, zerImageM00, CV_64F, M00, Point(-1,-1), 0, BORDER_DEFAULT);

	Mat ZerImageM11R;
	filter2D(NewAdaThresImage, ZerImageM11R, CV_64F, M11R, Point(-1, -1), 0, BORDER_DEFAULT);
	//////////filter2D( cannyImage, zerImageM11R, CV_64F, M11R, Point(-1, -1), 0, BORDER_DEFAULT);

	Mat ZerImageM11I;
	filter2D(NewAdaThresImage, ZerImageM11I, CV_64F, M11I, Point(-1, -1), 0, BORDER_DEFAULT);
	//////////filter2D( cannyImage, zerImageM11I, CV_64F, M11I, Point(-1, -1), 0, BORDER_DEFAULT);

	Mat ZerImageM20;
	filter2D(NewAdaThresImage, ZerImageM20, CV_64F, M20, Point(-1, -1), 0, BORDER_DEFAULT);
	//////////filter2D( cannyImage, zerImageM20, CV_64F, M20, Point(-1, -1), 0, BORDER_DEFAULT);

	Mat ZerImageM31R;
	filter2D(NewAdaThresImage, ZerImageM31R, CV_64F, M31R, Point(-1, -1), 0, BORDER_DEFAULT);
	//////////filter2D(cannyImage, zerImageM31R, CV_64F, M31R, Point(-1, -1), 0, BORDER_DEFAULT);

	Mat ZerImageM31I;
	filter2D(NewAdaThresImage, ZerImageM31I, CV_64F, M31I, Point(-1, -1), 0, BORDER_DEFAULT);
	//////////filter2D(cannyImage, zerImageM31I, CV_64F, M31I, Point(-1, -1), 0, BORDER_DEFAULT);

	Mat ZerImageM40;
	filter2D(NewAdaThresImage, ZerImageM40, CV_64F, M40, Point(-1, -1), 0, BORDER_DEFAULT);
	//////////filter2D(cannyImage, zerImageM40, CV_64F, M40, Point(-1, -1), 0, BORDER_DEFAULT);



	int row_number = NewAdaThresImage.rows;
	int col_number = NewAdaThresImage.cols;
	//使用7個的7*7的Zernike模板（其本質是個矩陣）M00、M11R、M11I、M20、M31R、M31I、M40，分別與影象進行卷積，得到每個畫素點的七個Zernike矩Z00、Z11R、Z11I、Z20、Z31R、Z31I、Z40
	//對於每個點，根據它的七個Zernike矩，求得距離引數l和灰度差引數k，當l和k都滿足設定的條件時，則判讀該點為邊緣點，並進一步利用上述7個Zernike矩求出該點的亞畫素級座標
	//如果l或k不滿足設定的條件，則該點不是邊緣點，轉到下一個點求解距離引數l和灰度差引數k
	for (int i = 0; i < row_number; i++)
	{
		for (int j = 0; j <col_number; j++)
		{
			if (ZerImageM00.at<double>(i, j) == 0)
			{
				continue;
			}

			//compute theta
			//vector<vector<double> > theta2(0);
			double theta_temporary = atan2(ZerImageM31I.at<double>(i, j), ZerImageM31R.at<double>(i, j));
			//theta2[i].push_back(thetaTem);

			//compute Z11'/Z31'
			double rotated_z11 = 0.0;
			rotated_z11 = sin(theta_temporary)*(ZerImageM11I.at<double>(i, j)) + cos(theta_temporary)*(ZerImageM11R.at<double>(i, j));
			double rotated_z31 = 0.0;
			rotated_z31 = sin(theta_temporary)*(ZerImageM31I.at<double>(i, j)) + cos(theta_temporary)*(ZerImageM31R.at<double>(i, j));

			//compute l
			double l_method1 = sqrt((5 * ZerImageM40.at<double>(i, j) + 3 * ZerImageM20.at<double>(i, j)) / (8 * ZerImageM20.at<double>(i, j)));
			double l_method2 = sqrt((5 * rotated_z31 + rotated_z11) / (6 * rotated_z11));
			double l = (l_method1 + l_method2) / 2;
			//compute k/h
			double k, h;
		
			k = 3 * rotated_z11 / 2 / pow((1 - l_method2*l_method2), 1.5);
			h = (ZerImageM00.at<double>(i, j) - k*PI / 2 + k*asin(l_method2) + k*l_method2*sqrt(1 - l_method2*l_method2)) / PI;

			//judge the edge
			double k_value = 20.0;
		
			double l_value = sqrt(2) / g_N;
		
			double absl = abs(l_method2 - l_method1);
			if (k >= k_value && absl <= l_value)
			{
				Point2d point_temporary;
				point_temporary.x = j + g_N*l*cos(theta_temporary) / 2;
				point_temporary.y = i + g_N*l*sin(theta_temporary) / 2;
				SubEdgePoints.push_back(point_temporary);
			}
			else
			{
				continue;
			}
		}
	}
	//顯示所檢測到的亞畫素邊緣
	for (size_t i = 0; i < SubEdgePoints.size(); i++)
	{
	    Point center_forshow(cvRound(SubEdgePoints[i].x), cvRound(SubEdgePoints[i].y));
		circle(OriginalImage, center_forshow, 1, Scalar(0, 97, 255), 1, 8, 0);
	}
	imshow("亞畫素邊緣", OriginalImage);
	
	waitKey(0);
	return 0;
}


---------------------------------------7.28更新分割線-------------------------------------------------

1.因為留言要原始碼的朋友比較多，所以我把上面的程式碼補完整了，新的程式碼能夠實現基本的過程，但是針對不同的應用，肯定不是最優的效果，可能需要根據不同的應用場合進行調整；

2.評論裡提到的雙層輪廓，我能想到的是因為對線進行邊緣檢測操作，因為線的兩側相當於都是前景和背景的邊緣，所以兩側各會有一層邊緣，不管是Canny（事實上不少Zernike論文都建議先用canny再用Zernike）還是二值化，只要送入Zernike的是線而不是區域，就會有雙層輪廓，其他的原因我暫時還沒想到，歡迎大家賜教；

3.雙層輪廓並不一定是壞事，比如做圓檢測，最後是要用這些邊緣點來擬合圓，那麼兩層輪廓擬合出來的圓放大了說應該是處於兩層輪廓中間，這樣的結果其實更符合實際情況，但是如果不是規則輪廓，可能會對邊緣點的表達和擬合有些問題；

相關推薦

OpenCV實現基於Zernike矩的亞畫素邊緣檢測

    在做物體檢測時，由於成本和應用場合的限制，不能夠一味地增加相機的解析度，或者已經用了解析度很高的相機，但是視野範圍很大，仍然無法實現很高的精度，這時就要考慮亞畫素技術，亞畫素技術就是在兩個畫素點之間進行進一步的細分，從而得到亞畫素級別的邊緣點的座標（也就是float

OpenCv-C++-亞畫素級別角點檢測(檢測子畫素中的corner的位置)

使用亞畫素級別角點檢測，返回角點的浮點數值，它的精度比整數畫素更準確。可以用cornerSubPix()函式將角點定位到子畫素，從而取得亞畫素級別的角點檢測效果。 使用函式： void cv::cornerSubPix ( InputArray image, InputOutputArray

OpenCV亞畫素角點cornerSubPixel()原始碼分析

　　上一篇部落格中講到了goodFeatureToTrack（）這個API函式能夠獲取影象中的強角點。但是獲取的角點座標是整數，但是通常情況下，角點的真實位置並不一定在整數畫素位置，因此為了獲取更為精確的角點位置座標，需要角點座標達到亞畫素（subPixel）精度。 1. 求取亞畫素精度的原理 　　找到一篇講

OpenCV：訪問影象中的畫素

OpenCV影象處理運算元都是一個函式。 作用： 接受一個輸入或多個輸入，產生輸出影象。 格式： g(x)=f(h(x))//單個輸入， 或 g(x)=f(h0(x),h1(x),...,hn(x))//多個輸入 影象處理變換中典型的操作：點操作 點操作： 影象亮度和對

論文學習 + 論文寫作 | 最前沿的亞畫素運動補償 + 視訊超分辨：Detail-revealing Deep Video Super-resolution

目錄 O. Abstract I. Introduction II. Relative Work III. Sub-pixel Motion Compensation (SPMC) 摘抄本 0. Abstract 1. Introduction 2. Re

cv::cornerSubPix()亞畫素角點檢測

轉自https://blog.csdn.net/guduruyu/article/details/69537083 cv::goodFeaturesToTrack()提取到的角點只能達到畫素級別，在很多情況下並不能滿足實際的需求，這時，我們則需要使用cv::cornerSubPix()對檢測到的角點

萌新學習手冊：亞畫素級別角點檢測

在現實應用中我們的角點往往並不是整數，所以為了提高我們尋找角點的精確程度，我們需要進行處理，一般有三種方法 1.插值方法 2.基於影象矩計算 3.圖線擬合（高斯曲面，多項式，橢圓曲面）其中高斯曲面最為常用 具體的使用過程就是我們先用Shi-Tomasi角點檢測 把角點

cornerSubPixel亞畫素角點檢測原始碼分析（最小二乘法）

OpenCV的中有cornerSubPixel（）這個API函式用來針對初始的整數角點座標進行亞畫素精度的優化，該函式原型如下： C++: void cornerSubPix(InputArray image, InputOutputArray corners, Siz

實現手機端的1px畫素 border的實現

實現1px畫素需要知道min-device-pixel-ratio裝置畫素比。如不清除的情移步瞭解下張鑫旭大神的裝置畫素比devicePixelRatio簡單介紹 @media (-webkit-device-pixel-ratio: 1.5),(device-pixel-ratio:

亞畫素與halcon

朋友發來兩個小專案，要求亞畫素精度。突然想問幾個問題：1、何為亞畫素？2、何為亞畫素精度？3、使用亞畫素測量，系統應注意什麼？ 1、何謂亞畫素？ 面陣攝像機的成像面以畫素為最小單位。例如某CMOS攝像晶片，其畫素間距為5.2微米。攝像機拍攝時，將物理世界中連續的影象進行了離

opencv學習（五）之畫素遍歷三種方式耗時分析

前面用兩篇介紹了畫素的顏色空間縮減、查詢表、遍歷畫素的三種方式、程式計時等，也說了一下每種方法的優缺點，現在用一個綜合型的程式進行對比。方式是用三種方式對lena影象（220x220）進行處理，使其顏色種類從256中變成64種。在顏色空間縮減方法中講過這種方式，

opencv讀取彩色/灰度圖片畫素值並存儲在本地檔案中c++程式碼例項及執行結果

c++程式碼彩色圖片#include<opencv2/opencv.hpp> #include<fstream> using namespace std; using namespace cv; int main(int argc, char* ar

opencv點選滑鼠讀取圖片畫素c++程式碼

程式碼#include<opencv2/opencv.hpp> using namespace std; void onMouse(int event, int x, int y, int flags, void *param) { cv::Mat *im =

OpenCV取影象中的每個畫素值

import cv2import numpy as pydef salt(img): IMAGE_FILE = img img = cv2.imread(IMAGE_FILE) print(img.shape[1]) for i in

Opencv Mat 類詳解以及畫素點基本讀取方法

class CV_EXPORTS Mat { public: //! default constructor Mat(); //! constructs 2D matrix of the specified size and type // (_type is CV_8UC1,

Opencv獲取並改變影象的畫素值

#include <iostream> using namespace std; #include "cv.h" #include "highgui.h" #pragma  comment(lib,"cxcore.lib") #pragma  comment(l

OpenCV從入門到放棄(五):畫素！

一.概念 1.影象本質上面是由數值組成的矩陣。矩陣中的一個元素對應一個畫素。 2.對於灰度影象（黑白影象），畫素是8位無符號數（CV_8U）。0表示黑色，255表示白色。對於彩色影象，是用三原色資料合成彩色。3個8位（CV_8UC3）的數值組成矩陣的一個元素

OpenCV提取某個特定點的畫素值

OpenCV 預設的顏色空間是 BGR，類似於RGB，但不是RGB。 即某個畫素的第一個值是Blue的值， 第二個是Green的值，第三個是R的值。 程式的核心是 用 cvGet2D(圖片， y座標， x座標)    獲取 CvScalar 物件。 注意y,x順序，不

opencv中遍歷每一個畫素點進行處理

轉載：http://blog.csdn.net/wwwsssZheRen/article/details/743155911.用動態地址操作畫素：Mat srcImage(100, 100, CV_8UC3, Scalar(200,20,100));  imshow("顯示影

opencv——操作影象中每一個畫素

以下均針對opencv1.0 第一種方法：使用cvGet2D及cvSet2D 對於單通道影象： IplImage* img = cvCreateImage(cvSize(640, 480), IPL_DEPTH_8U, 1); for