2014-10-05 16:17:46 gongem 阅读数 12147
  • OpenCV3.2 Java图像处理视频学习教程

    OpenCV3.2 Java图像处理视频培训课程:基于OpenCV新版本3.2.0详细讲述Java OpenCV图像处理部分内容,包括Mat对象使用、图像读写、 基于常用核心API讲述基本原理、使用方法、参数、代码演示、图像处理思路与流程讲授。主要内容包括opencv像素操作、滤波、边缘提取、直线与圆检测、形态学操作与分水岭、图像金子塔融合重建、多尺度模板匹配、opencv人脸检测、OpenCV跟Tomcat使用实现服务器端图像处理服务。

    4253 人正在学习 去看看 贾志刚

1. 边缘检测的概念

边缘检测是图像处理与计算机视觉中极为重要的一种分析图像的方法,至少在我做图像分析与识别时,边缘是我最喜欢的图像特征。边缘检测的目的就是找到图像中亮度变化剧烈的像素点构成的集合,表现出来往往是轮廓。在对现实世界的图像采集中,有下面4种情况会表现在图像中时形成一个边缘。

  1. 深度的不连续(物体处在不同的物平面上);
  2. 表面方向的不连续(如正方体的不同的两个面);
  3. 物体材料不同(这样会导致光的反射系数不同);
  4. 场景中光照不同(如被树萌投向的地面);


    image上面的图像是图像中水平方向7个像素点的灰度值显示效果,我们很容易地判断在第4和第5个像素之间有一个边缘,因为它俩之间发生了强烈的灰度跳变。在实际的边缘检测中,边缘远没有上图这样简单明显,我们需要取对应的阈值来区分出它们。

    原图 Roberts边缘 Prewitt边缘
    image image image
    Sobel边缘 Log边缘 Canny边缘
    image image image
     


    2. 边缘检测的基本方法

    2.1 一阶微分边缘算子

    一阶微分边缘算子也称为梯度边缘算子,它是利用图像在边缘处的阶跃性,即图像梯度在边缘取得极大值的特性进行边缘检测。梯度是

    一个矢量,它具有方向θ和模|ΔI|

    梯度的方向提供了边缘的趋势信息,因为梯度方向始终是垂直于边缘方向,梯度的模值大小提供了边缘的强度信息。

    在实际使用中,通常利用有限差分进行梯度近似。对于上面的公式,我们有如下的近似:

    2.2 Roberts边缘检测算子

    1963年,Roberts提出了这种寻找边缘的算子。Roberts边缘算子是一个2x2的模板,采用的是对角方向相邻的两个像素之差。从图像处理的实际效果来看,边缘定位较准,对噪声敏感。在Roberts检测算子中:

    可以导出Roberts在点(i+1/2,j+1/2)处的水平与竖直边缘检测卷积核为:

    2.3 Prewitt边缘检测算子

    Prewitt利用周围邻域8个点的灰度值来估计中心的梯度,它的梯度计算公式如下:

    所以,Prewitt的卷积核为:

    2.4 Sobel边缘检测算子

    比起Prewitt算子,Sobel也是用周围8个像素来估计中心像素的梯度,但是Sobel算子认为靠近中心像素的点应该给予更高的权重,所以Sobel算子把与中心像素4邻接的像素的权重设置为2或-2。

    Sobel边缘检测算子的卷积核为:

    Sobel进行边缘检测的实现可以参考我原来写的一篇博文:图像特征检测:sobel边缘检测,重要的是梯度图像计算后的阈值的确定与边缘的非极大值抑制算法,Roberts与Prewitt原理与sobel一致。

    3. 二阶微分算子

    学过微积分我们都知道,边缘即是图像的一阶导数局部最大值的地方,那么也意味着该点的二阶导数为零。二阶微分边缘检测算子就是利用图像在边缘处的阶跃性导致图像二阶微分在边缘处出现零值这一特性进行边缘检测的。

    对于图像的二阶微分可以用拉普拉斯算子来表示:

    我们在像素点(i,j)3×3的邻域内,可以有如下的近似:

    对应的二阶微分卷积核为:

    所以二阶微分检测边缘的方法就分两步:1)用上面的Laplace核与图像进行卷积;2)对卷积后的图像,取得那些卷积结果为0的点。

    虽然上述使用二阶微分检测边缘的方法简单,但它的缺点是对噪声十分敏感,同时也没有能够提供边缘的方向信息。为了实现对噪声的抑制,Marr等提出了LOG的方法。

    为了减少噪声对边缘的影响,首先图像要进行低通滤波,LOG采用了高斯函数作为低通滤波器。高斯函数为:

    上面的公式中σ决定了对图像的平滑程度。高斯函数生成的滤波模板尺寸一般设定为6σ+1(加1是会了使滤波器的尺寸为奇数)。使用高斯函数对图像进行滤波并对图像滤波结果进行二阶微分运算的过程,可以转换为先对高斯函数进行二阶微分,再利用高斯函数的二阶微分结果对图像进行卷积运算:

    关于LOG算子的计算在上一篇文章:图像特征提取:斑点检测中有实现的代码。

    4. Canny边缘检测

    canny边缘检测实际上是一种一阶微分算子检测算法,但为什么这里拿出来说呢,因为它几乎是边缘检测算子中最为常用的一种,也是个人认为现在最优秀的边缘检测算子。Canny提出了边缘检测算子优劣评判的三条标准:

    Canny边缘检测之所以优秀是因为它在一阶微分算子的基础上,增加了非最大值抑制和双阈值两项改进。利用非极大值抑制不仅可以有效地抑制多响应边缘,而且还可以提高边缘的定位精度;利用双阈值可以有效减少边缘的漏检率。

    Canny边缘检测主要分四步进行:

    1. 去噪声;
    2. 计算梯度与方向角;
    3. 非最大值抑制;
    4. 滞后阈值化;

    其中前两步很简单,先用一个高斯滤波器对图像进行滤波,然后用Sobel水平和竖直检测子与图像卷积,来计算梯度和方向角。

    非极大值抑制

    图像梯度幅值矩阵中的元素值越大,说明图像中该点的梯度值越大,但这不不能说明该点就是边缘(这仅仅是属于图像增强的过程)。在Canny算法中,非极大值抑制是进行边缘检测的重要步骤,通俗意义上是指寻找像素点局部最大值,将非极大值点所对应的灰度值置为0,这样可以剔除掉一大部分非边缘的点。

    image

    根据上图可知,要进行非极大值抑制,就首先要确定像素点C的灰度值在其8值邻域内是否为最大。图中蓝色的线条方向为C点的梯度方向,这样就可以确定其局部的最大值肯定分布在这条线上,也即出了C点外,梯度方向的交点dTmp1和dTmp2这两个点的值也可能会是局部最大值。因此,判断C点灰度与这两个点灰度大小即可判断C点是否为其邻域内的局部最大灰度点。如果经过判断,C点灰度值小于这两个点中的任一个,那就说明C点不是局部极大值,那么则可以排除C点为边缘。这就是非极大值抑制的工作原理。

    在理解的过程中需要注意以下两点:

    1. 中非最大抑制是回答这样一个问题:“当前的梯度值在梯度方向上是一个局部最大值吗?” 所以,要把当前位置的梯度值与梯度方向上两侧的梯度值进行比较;
    2. 梯度方向垂直于边缘方向。

    但实际上,我们只能得到C点邻域的8个点的值,而dTmp1和dTmp2并不在其中,要得到这两个值就需要对该两个点两端的已知灰度进行线性插值,也即根据图中的g1和g2对dTmp1进行插值,根据g3和g4对dTmp2进行插值,这要用到其梯度方向。

    滞后阈值化

    由于噪声的影响,经常会在本应该连续的边缘出现断裂的问题。滞后阈值化设定两个阈值:一个为高阈值Th,一个为低阈值Tl。如果任何像素边缘算子的影响超过高阈值,将这些像素标记为边缘;响应超过低阈值(高低阈值之间)的像素,如果与已经标记为边缘的像素4-邻接或8-邻接,则将这些像素也标记为边缘。所以不整个过程描述如下:

    1. 如果该像素的梯度值小于Tl,则该像素为非边缘像素;
    2. 如果该像素的梯度值大于Th,则该像素为边缘像素;
    3. 如果该像素的梯度值介于TlTh之间,需要进一步检测该像素的3×3邻域内的8个点,如果这8个点内有一个或以上的点梯度超过了Th,则该像素为边缘像素,否则不是边缘像素。
     


    Canny边缘检测的实现


    #include <iostream>


    #include "opencv2/core/core.hpp"
    #include "opencv2/highgui/highgui.hpp"
    #include "opencv2/imgproc/imgproc.hpp"




    using namespace std;
    using namespace cv;




    void getCannyEdge(const Mat& imgSrc, Mat& imgDst, double lowThresh = -1, double highThresh = -1, double sigma = 1);
    double getGaussianThresh(const Mat& inputArray, double percentage);


    int main(int argc, char** argv)
    {
        Mat image = imread(argv[1]);
        Mat imgCanny;
        getCannyEdge(image,imgCanny);
        t = ((double)getTickCount() - t)/getTickFrequency();
        imwrite("imgCanny.png", imgCanny);


        return 0;
    }


    void getCannyEdge(const Mat& imgSrc, Mat& imgDst, double lowThresh, double highThresh, double sigma)
    {
        Mat gray;
        if (imgSrc.channels() == 3)
        {
            cvtColor(imgSrc, gray, CV_BGR2GRAY);
        }
        else
        {
            gray = imgSrc.clone();
        }
        gray.convertTo(gray, CV_64F);
        gray = gray / 255;
        
        double gaussianDieOff = .0001;
        double percentOfPixelsNotEdges = .7; // Used for selecting thresholds
        double thresholdRatio = .4;   // Low thresh is this fraction of the high


        int possibleWidth = 30;
        double ssq = sigma * sigma;
        for (int i = 1; i <= possibleWidth; i++)
        {
            if (exp(-(i * i) / (2* ssq)) < gaussianDieOff)
            {
                possibleWidth = i - 1;
                break;
            }
        }


        if (possibleWidth == 30)
        {
            possibleWidth = 1; // the user entered a reallly small sigma
        }


        // get the 1D gaussian filter
        int winSz = 2 * possibleWidth + 1;
        Mat gaussKernel1D(1, winSz, CV_64F);
        double* kernelPtr = gaussKernel1D.ptr<double>(0);
        for (int i = 0; i < gaussKernel1D.cols; i++)
        {
            kernelPtr[i] = exp(-(i - possibleWidth) * (i - possibleWidth) / (2 * ssq)) / (2 * CV_PI * ssq);
        }


        
        // get the derectional derivatives of gaussian kernel
        Mat dGaussKernel(winSz, winSz, CV_64F);
        for (int i = 0; i < dGaussKernel.rows; i++)
        {
            double* linePtr = dGaussKernel.ptr<double>(i);
            for (int j = 0; j< dGaussKernel.cols; j++)
            {
                linePtr[j] = - (j - possibleWidth) * exp(-((i - possibleWidth) * (i - possibleWidth) + (j - possibleWidth) * (j - possibleWidth)) / (2 * ssq)) / (CV_PI * ssq);
            }
        }




        /* smooth the image out*/
        Mat imgSmooth;
        filter2D(gray, imgSmooth, -1, gaussKernel1D);
        filter2D(imgSmooth, imgSmooth, -1, gaussKernel1D.t());
        /*apply directional derivatives*/


        Mat imgX, imgY;
        filter2D(imgSmooth, imgX, -1, dGaussKernel);
        filter2D(imgSmooth, imgY, -1, dGaussKernel.t());


        Mat imgMag;
        sqrt(imgX.mul(imgX) + imgY.mul(imgY), imgMag);
        double magMax;
        minMaxLoc(imgMag, NULL, &magMax, NULL, NULL);


        if (magMax > 0 )
        {
            imgMag = imgMag / magMax;
        }


        
        if (lowThresh == -1 || highThresh == -1)
        {
            highThresh = getGaussianThresh(imgMag, percentOfPixelsNotEdges);
            lowThresh = thresholdRatio * highThresh;
        }








        Mat imgStrong = Mat::zeros(imgMag.size(), CV_8U);
        Mat imgWeak = Mat::zeros(imgMag.size(), CV_8U);
        
        
        for (int dir = 1; dir <= 4; dir++)
        {
            Mat gradMag1(imgMag.size(), imgMag.type());
            Mat gradMag2(imgMag.size(), imgMag.type());
            Mat idx = Mat::zeros(imgMag.size(), CV_8U);
            if (dir == 1)
            {
                Mat dCof = abs(imgY / imgX);
                idx = (imgY <= 0 & imgX > -imgY) | (imgY >= 0 & imgX < -imgY);
                idx.row(0).setTo(Scalar(0));
                idx.row(idx.rows - 1).setTo(Scalar(0));
                idx.col(0).setTo(Scalar(0));
                idx.col(idx.cols - 1).setTo(Scalar(0));
                for (int i = 1; i < imgMag.rows - 1; i++)
                {
                    for (int j = 1; j < imgMag.cols - 1; j++)
                    {
                        gradMag1.at<double>(i,j) = (1 - dCof.at<double>(i,j)) * imgMag.at<double>(i,j + 1) + dCof.at<double>(i,j) * imgMag.at<double>(i - 1,j + 1);
                        gradMag2.at<double>(i,j) = (1 - dCof.at<double>(i,j)) * imgMag.at<double>(i,j - 1) + dCof.at<double>(i,j) * imgMag.at<double>(i + 1,j - 1);
                    }
                }
            }
            else if(dir == 2)
            {
                Mat dCof = abs(imgX / imgY);
                idx = (imgX > 0 & -imgY >= imgX) | (imgX < 0 & -imgY <= imgX);
                for (int i = 1; i < imgMag.rows - 1; i++)
                {
                    for (int j = 1; j < imgMag.cols - 1; j++)
                    {
                        gradMag1.at<double>(i,j) = (1 - dCof.at<double>(i,j)) * imgMag.at<double>(i - 1,j) + dCof.at<double>(i,j) * imgMag.at<double>(i - 1,j + 1);
                        gradMag2.at<double>(i,j) = (1 - dCof.at<double>(i,j)) * imgMag.at<double>(i + 1,j) + dCof.at<double>(i,j) * imgMag.at<double>(i + 1,j - 1);
                    }
                }
            }
            else if(dir == 3)
            {
                Mat dCof = abs(imgX / imgY);
                idx = (imgX <= 0 & imgX > imgY) | (imgX >= 0 & imgX < imgY);
                for (int i = 1; i < imgMag.rows - 1; i++)
                {
                    for (int j = 1; j < imgMag.cols - 1; j++)
                    {
                        gradMag1.at<double>(i,j) = (1 - dCof.at<double>(i,j)) * imgMag.at<double>(i - 1,j) + dCof.at<double>(i,j) * imgMag.at<double>(i - 1,j - 1);
                        gradMag2.at<double>(i,j) = (1 - dCof.at<double>(i,j)) * imgMag.at<double>(i + 1,j) + dCof.at<double>(i,j) * imgMag.at<double>(i + 1,j + 1);
                    }
                }
            
            }
            else
            {
                Mat dCof = abs(imgY / imgX);
                idx = (imgY <0 & imgX <= imgY) | (imgY > 0 & imgX >= imgY);
                for (int i = 1; i < imgMag.rows - 1; i++)
                {
                    for (int j = 1; j < imgMag.cols - 1; j++)
                    {
                        gradMag1.at<double>(i,j) = (1 - dCof.at<double>(i,j)) * imgMag.at<double>(i,j - 1) + dCof.at<double>(i,j) * imgMag.at<double>(i - 1,j - 1);
                        gradMag2.at<double>(i,j) = (1 - dCof.at<double>(i,j)) * imgMag.at<double>(i,j + 1) + dCof.at<double>(i,j) * imgMag.at<double>(i + 1,j + 1);
                    }
                }
            }


            Mat idxLocalMax = idx & ((imgMag >= gradMag1) & (imgMag >= gradMag2));




            imgWeak = imgWeak | ((imgMag > lowThresh) & idxLocalMax);
            imgStrong= imgStrong| ((imgMag > highThresh) & imgWeak);


        }


        imgDst = Mat::zeros(imgWeak.size(),imgWeak.type());
        for (int i = 1; i < imgWeak.rows - 1; i++)
        {
            uchar* pWeak = imgWeak.ptr<uchar>(i);
            uchar* pDst = imgDst.ptr<uchar>(i);
            uchar* pStrPre = imgStrong.ptr<uchar>(i - 1);
            uchar* pStrMid = imgStrong.ptr<uchar>(i);
            uchar* pStrAft = imgStrong.ptr<uchar>(i + 1);
            for (int j = 1; j < imgWeak.cols - 1; j++)
            {
                if (!pWeak[j])
                {
                    continue;
                }
                if (pStrMid[j])
                {
                    pDst[j] = 255;
                }
                if (pStrMid[j-1] || pStrMid[j+1] || pStrPre[j-1] || pStrPre[j] || pStrPre[j+1] || pStrAft[j-1] || pStrAft[j] ||pStrAft[j+1])
                {
                    pDst[j] = 255;
                }
            }
        }
    }


    double getGaussianThresh(const Mat& inputArray, double percentage)
    {
        double thresh = -1.0;
        // compute the 64-hist of inputArray
        int nBins = 64;
        double minValue, maxValue;
        minMaxLoc(inputArray, &minValue, &maxValue, NULL, NULL);
        double step = (maxValue - minValue) / nBins;


        vector<unsigned> histBin(nBins,0);
        for (int i = 0; i < inputArray.rows; i++)
        {
            const double* pData = inputArray.ptr<double>(i);
            for(int j = 0; j < inputArray.cols; j++)
            {


                int index = (pData[j] - minValue) / step;
                histBin[index]++;
            }
        }
        unsigned cumSum = 0; 
        for (int i = 0; i < nBins; i++)
        {
            cumSum += histBin[i];


            if (cumSum > percentage * inputArray.rows * inputArray.cols)
            {
                thresh = (i + 1) / 64.0;
                break;
            }
        }
        return thresh;
        
    }
  5. 文章引子:

         

    http://www.cnblogs.com/ronny/p/4001910.html

原图 Roberts边缘 Prewitt边缘
image image image
Sobel边缘 Log边缘 Canny边缘
image image image
2014-11-16 16:50:47 jia20003 阅读数 60462
  • OpenCV3.2 Java图像处理视频学习教程

    OpenCV3.2 Java图像处理视频培训课程:基于OpenCV新版本3.2.0详细讲述Java OpenCV图像处理部分内容,包括Mat对象使用、图像读写、 基于常用核心API讲述基本原理、使用方法、参数、代码演示、图像处理思路与流程讲授。主要内容包括opencv像素操作、滤波、边缘提取、直线与圆检测、形态学操作与分水岭、图像金子塔融合重建、多尺度模板匹配、opencv人脸检测、OpenCV跟Tomcat使用实现服务器端图像处理服务。

    4253 人正在学习 去看看 贾志刚

图像处理之Canny 边缘检测

一:历史

Canny边缘检测算法是1986年有John F. Canny开发出来一种基于图像梯度计算的边缘

检测算法,同时Canny本人对计算图像边缘提取学科的发展也是做出了很多的贡献。尽

管至今已经许多年过去,但是该算法仍然是图像边缘检测方法经典算法之一。

二:Canny边缘检测算法

经典的Canny边缘检测算法通常都是从高斯模糊开始,到基于双阈值实现边缘连接结束

。但是在实际工程应用中,考虑到输入图像都是彩色图像,最终边缘连接之后的图像要

二值化输出显示,所以完整的Canny边缘检测算法实现步骤如下:

1.      彩色图像转换为灰度图像

2.      对图像进行高斯模糊

3.      计算图像梯度,根据梯度计算图像边缘幅值与角度

4.      非最大信号压制处理(边缘细化)

5.      双阈值边缘连接处理

6.      二值化图像输出结果

三:各步详解与代码实现

1.      彩色图像转灰度图像

根据彩色图像RGB转灰度公式:gray  =  R * 0.299 + G * 0.587 + B * 0.114

将彩色图像中每个RGB像素转为灰度值的代码如下:

int gray = (int) (0.299 * tr + 0.587 * tg + 0.114 * tb);

2.      对图像进行高斯模糊

图像高斯模糊时,首先要根据输入参数确定高斯方差与窗口大小,这里我设置默认方

差值窗口大小为16x16,根据这两个参数生成高斯卷积核算子的代码如下:

		float kernel[][] = new float[gaussianKernelWidth][gaussianKernelWidth];
		for(int x=0; x<gaussianKernelWidth; x++)
		{
			for(int y=0; y<gaussianKernelWidth; y++)
			{
				kernel[x][y] = gaussian(x, y, gaussianKernelRadius);
			}
		}

获取了高斯卷积算子之后,我们就可以对图像高斯卷积模糊,关于高斯图像模糊更详

细的解释可以参见这里:http://blog.csdn.net/jia20003/article/details/7234741实现

图像高斯卷积模糊的代码如下:

// 高斯模糊 -灰度图像
int krr = (int)gaussianKernelRadius;
for (int row = 0; row < height; row++) {
	for (int col = 0; col < width; col++) {
		index = row * width + col;
		double weightSum = 0.0;
		double redSum = 0;
		for(int subRow=-krr; subRow<=krr; subRow++)
		{
			int nrow = row + subRow;
			if(nrow >= height || nrow < 0)
			{
				nrow = 0;
			}
			for(int subCol=-krr; subCol<=krr; subCol++)
			{
				int ncol = col + subCol;
				if(ncol >= width || ncol <=0)
				{
					ncol = 0;
				}
				int index2 = nrow * width + ncol;
				int tr1 = (inPixels[index2] >> 16) & 0xff;
				redSum += tr1*kernel[subRow+krr][subCol+krr];
				weightSum += kernel[subRow+krr][subCol+krr];
			}
		}
		int gray = (int)(redSum / weightSum);
		outPixels[index] = gray;
	}
}

3.      计算图像X方向与Y方向梯度,根据梯度计算图像边缘幅值与角度大小

高斯模糊的目的主要为了整体降低图像噪声,目的是为了更准确计算图像梯度及边缘

幅值。计算图像梯度可以选择算子有Robot算子、Sobel算子、Prewitt算子等。关于

图像梯度计算更多的解释可以看这里:

http://blog.csdn.net/jia20003/article/details/7664777。

这里采用更加简单明了的2x2的算子,其数学表达如下:


// 计算梯度-gradient, X放与Y方向
data = new float[width * height];
magnitudes = new float[width * height];
for (int row = 0; row < height; row++) {
	for (int col = 0; col < width; col++) {
		index = row * width + col;
		// 计算X方向梯度
		float xg = (getPixel(outPixels, width, height, col, row+1) - 
				getPixel(outPixels, width, height, col, row) + 
				getPixel(outPixels, width, height, col+1, row+1) -
				getPixel(outPixels, width, height, col+1, row))/2.0f;
		float yg = (getPixel(outPixels, width, height, col, row)-
				getPixel(outPixels, width, height, col+1, row) +
				getPixel(outPixels, width, height, col, row+1) -
				getPixel(outPixels, width, height, col+1, row+1))/2.0f;
		// 计算振幅与角度
		data[index] = hypot(xg, yg);
		if(xg == 0)
		{
			if(yg > 0)
			{
				magnitudes[index]=90;						
			}
			if(yg < 0)
			{
				magnitudes[index]=-90;
			}
		}
		else if(yg == 0)
		{
			magnitudes[index]=0;
		}
		else
		{
			magnitudes[index] = (float)((Math.atan(yg/xg) * 180)/Math.PI);					
		}
		// make it 0 ~ 180
		magnitudes[index] += 90;
	}
}

在获取了图像每个像素的边缘幅值与角度之后

4.      非最大信号压制

信号压制本来是数字信号处理中经常用的,这里的非最大信号压制主要目的是实现边

缘细化,通过该步处理边缘像素进一步减少。非最大信号压制主要思想是假设3x3的

像素区域,中心像素P(x,y) 根据上一步中计算得到边缘角度值angle,可以将角度分

为四个离散值0、45、90、135分类依据如下:

其中黄色区域取值范围为0~22.5 与157.5~180

绿色区域取值范围为22.5 ~ 67.5

蓝色区域取值范围为67.5~112.5

红色区域取值范围为112.5~157.5

分别表示上述四个离散角度的取值范围。得到角度之后,比较中心像素角度上相邻

两个像素,如果中心像素小于其中任意一个,则舍弃该边缘像素点,否则保留。一

个简单的例子如下:


// 非最大信号压制算法 3x3
Arrays.fill(magnitudes, 0);
for (int row = 0; row < height; row++) {
	for (int col = 0; col < width; col++) {
		index = row * width + col;
		float angle = magnitudes[index];
		float m0 = data[index];
		magnitudes[index] = m0;
		if(angle >=0 && angle < 22.5) // angle 0
		{
			float m1 = getPixel(data, width, height, col-1, row);
			float m2 = getPixel(data, width, height, col+1, row);
			if(m0 < m1 || m0 < m2)
			{
				magnitudes[index] = 0;
			}
		}
		else if(angle >= 22.5 && angle < 67.5) // angle +45
		{
			float m1 = getPixel(data, width, height, col+1, row-1);
			float m2 = getPixel(data, width, height, col-1, row+1);
			if(m0 < m1 || m0 < m2)
			{
				magnitudes[index] = 0;
			}
		}
		else if(angle >= 67.5 && angle < 112.5) // angle 90
		{
			float m1 = getPixel(data, width, height, col, row+1);
			float m2 = getPixel(data, width, height, col, row-1);
			if(m0 < m1 || m0 < m2)
			{
				magnitudes[index] = 0;
			}
		}
		else if(angle >=112.5 && angle < 157.5) // angle 135 / -45
		{
			float m1 = getPixel(data, width, height, col-1, row-1);
			float m2 = getPixel(data, width, height, col+1, row+1);
			if(m0 < m1 || m0 < m2)
			{
				magnitudes[index] = 0;
			}
		}
		else if(angle >=157.5) // 跟零度是一致的,感谢一位网友发现了这个问题
		{
			float m1 = getPixel(data, width, height, col+1, row);
			float m2 = getPixel(data, width, height, col-1, row);
			if(m0 < m1 || m0 < m2)
			{
				magnitudes[index] = 0;
			}
		}
	}
}

1.      双阈值边缘连接

非最大信号压制以后,输出的幅值如果直接显示结果可能会少量的非边缘像素被包

含到结果中,所以要通过选取阈值进行取舍,传统的基于一个阈值的方法如果选择

的阈值较小起不到过滤非边缘的作用,如果选择的阈值过大容易丢失真正的图像边

缘,Canny提出基于双阈值(Fuzzy threshold)方法很好的实现了边缘选取,在实际

应用中双阈值还有边缘连接的作用。双阈值选择与边缘连接方法通过假设两个阈值

其中一个为高阈值TH另外一个为低阈值TL则有

a.      对于任意边缘像素低于TL的则丢弃

b.      对于任意边缘像素高于TH的则保留

c.      对于任意边缘像素值在TL与TH之间的,如果能通过边缘连接到一个像素大于

TH而且边缘所有像素大于最小阈值TL的则保留,否则丢弃。代码实现如下:

Arrays.fill(data, 0);
int offset = 0;
for (int row = 0; row < height; row++) {
	for (int col = 0; col < width; col++) {
		if(magnitudes[offset] >= highThreshold && data[offset] == 0)
		{
			edgeLink(col, row, offset, lowThreshold);
		}
		offset++;
	}
}
基于递归的边缘寻找方法edgeLink的代码如下:

private void edgeLink(int x1, int y1, int index, float threshold) {
	int x0 = (x1 == 0) ? x1 : x1 - 1;
	int x2 = (x1 == width - 1) ? x1 : x1 + 1;
	int y0 = y1 == 0 ? y1 : y1 - 1;
	int y2 = y1 == height -1 ? y1 : y1 + 1;
	
	data[index] = magnitudes[index];
	for (int x = x0; x <= x2; x++) {
		for (int y = y0; y <= y2; y++) {
			int i2 = x + y * width;
			if ((y != y1 || x != x1)
				&& data[i2] == 0 
				&& magnitudes[i2] >= threshold) {
				edgeLink(x, y, i2, threshold);
				return;
			}
		}
	}
}

6.      结果二值化显示 - 不说啦,直接点,自己看吧,太简单啦

// 二值化显示
for(int i=0; i<inPixels.length; i++)
{
	int gray = clamp((int)data[i]);
	outPixels[i] = gray > 0 ? -1 : 0xff000000;     
}
最终运行结果:


四:完整的Canny算法源代码

package com.gloomyfish.filter.study;

import java.awt.image.BufferedImage;
import java.util.Arrays;

public class CannyEdgeFilter extends AbstractBufferedImageOp {
	private float gaussianKernelRadius = 2f;
	private int gaussianKernelWidth = 16;
	private float lowThreshold;
	private float highThreshold;
	// image width, height
	private int width;
	private int height;
	private float[] data;
	private float[] magnitudes;

	public CannyEdgeFilter() {
		lowThreshold = 2.5f;
		highThreshold = 7.5f;
		gaussianKernelRadius = 2f;
		gaussianKernelWidth = 16;
	}

	public float getGaussianKernelRadius() {
		return gaussianKernelRadius;
	}

	public void setGaussianKernelRadius(float gaussianKernelRadius) {
		this.gaussianKernelRadius = gaussianKernelRadius;
	}

	public int getGaussianKernelWidth() {
		return gaussianKernelWidth;
	}

	public void setGaussianKernelWidth(int gaussianKernelWidth) {
		this.gaussianKernelWidth = gaussianKernelWidth;
	}

	public float getLowThreshold() {
		return lowThreshold;
	}

	public void setLowThreshold(float lowThreshold) {
		this.lowThreshold = lowThreshold;
	}

	public float getHighThreshold() {
		return highThreshold;
	}

	public void setHighThreshold(float highThreshold) {
		this.highThreshold = highThreshold;
	}

	@Override
	public BufferedImage filter(BufferedImage src, BufferedImage dest) {
		width = src.getWidth();
		height = src.getHeight();
		if (dest == null)
			dest = createCompatibleDestImage(src, null);
		// 图像灰度化
		int[] inPixels = new int[width * height];
		int[] outPixels = new int[width * height];
		getRGB(src, 0, 0, width, height, inPixels);
		int index = 0;
		for (int row = 0; row < height; row++) {
			int ta = 0, tr = 0, tg = 0, tb = 0;
			for (int col = 0; col < width; col++) {
				index = row * width + col;
				ta = (inPixels[index] >> 24) & 0xff;
				tr = (inPixels[index] >> 16) & 0xff;
				tg = (inPixels[index] >> 8) & 0xff;
				tb = inPixels[index] & 0xff;
				int gray = (int) (0.299 * tr + 0.587 * tg + 0.114 * tb);
				inPixels[index] = (ta << 24) | (gray << 16) | (gray << 8)
						| gray;
			}
		}
		
		// 计算高斯卷积核
		float kernel[][] = new float[gaussianKernelWidth][gaussianKernelWidth];
		for(int x=0; x<gaussianKernelWidth; x++)
		{
			for(int y=0; y<gaussianKernelWidth; y++)
			{
				kernel[x][y] = gaussian(x, y, gaussianKernelRadius);
			}
		}
		// 高斯模糊 -灰度图像
		int krr = (int)gaussianKernelRadius;
		for (int row = 0; row < height; row++) {
			for (int col = 0; col < width; col++) {
				index = row * width + col;
				double weightSum = 0.0;
				double redSum = 0;
				for(int subRow=-krr; subRow<=krr; subRow++)
				{
					int nrow = row + subRow;
					if(nrow >= height || nrow < 0)
					{
						nrow = 0;
					}
					for(int subCol=-krr; subCol<=krr; subCol++)
					{
						int ncol = col + subCol;
						if(ncol >= width || ncol <=0)
						{
							ncol = 0;
						}
						int index2 = nrow * width + ncol;
						int tr1 = (inPixels[index2] >> 16) & 0xff;
						redSum += tr1*kernel[subRow+krr][subCol+krr];
						weightSum += kernel[subRow+krr][subCol+krr];
					}
				}
				int gray = (int)(redSum / weightSum);
				outPixels[index] = gray;
			}
		}
		
		// 计算梯度-gradient, X放与Y方向
		data = new float[width * height];
		magnitudes = new float[width * height];
		for (int row = 0; row < height; row++) {
			for (int col = 0; col < width; col++) {
				index = row * width + col;
				// 计算X方向梯度
				float xg = (getPixel(outPixels, width, height, col, row+1) - 
						getPixel(outPixels, width, height, col, row) + 
						getPixel(outPixels, width, height, col+1, row+1) -
						getPixel(outPixels, width, height, col+1, row))/2.0f;
				float yg = (getPixel(outPixels, width, height, col, row)-
						getPixel(outPixels, width, height, col+1, row) +
						getPixel(outPixels, width, height, col, row+1) -
						getPixel(outPixels, width, height, col+1, row+1))/2.0f;
				// 计算振幅与角度
				data[index] = hypot(xg, yg);
				if(xg == 0)
				{
					if(yg > 0)
					{
						magnitudes[index]=90;						
					}
					if(yg < 0)
					{
						magnitudes[index]=-90;
					}
				}
				else if(yg == 0)
				{
					magnitudes[index]=0;
				}
				else
				{
					magnitudes[index] = (float)((Math.atan(yg/xg) * 180)/Math.PI);					
				}
				// make it 0 ~ 180
				magnitudes[index] += 90;
			}
		}
		
		// 非最大信号压制算法 3x3
		Arrays.fill(magnitudes, 0);
		for (int row = 0; row < height; row++) {
			for (int col = 0; col < width; col++) {
				index = row * width + col;
				float angle = magnitudes[index];
				float m0 = data[index];
				magnitudes[index] = m0;
				if(angle >=0 && angle < 22.5) // angle 0
				{
					float m1 = getPixel(data, width, height, col-1, row);
					float m2 = getPixel(data, width, height, col+1, row);
					if(m0 < m1 || m0 < m2)
					{
						magnitudes[index] = 0;
					}
				}
				else if(angle >= 22.5 && angle < 67.5) // angle +45
				{
					float m1 = getPixel(data, width, height, col+1, row-1);
					float m2 = getPixel(data, width, height, col-1, row+1);
					if(m0 < m1 || m0 < m2)
					{
						magnitudes[index] = 0;
					}
				}
				else if(angle >= 67.5 && angle < 112.5) // angle 90
				{
					float m1 = getPixel(data, width, height, col, row+1);
					float m2 = getPixel(data, width, height, col, row-1);
					if(m0 < m1 || m0 < m2)
					{
						magnitudes[index] = 0;
					}
				}
				else if(angle >=112.5 && angle < 157.5) // angle 135 / -45
				{
					float m1 = getPixel(data, width, height, col-1, row-1);
					float m2 = getPixel(data, width, height, col+1, row+1);
					if(m0 < m1 || m0 < m2)
					{
						magnitudes[index] = 0;
					}
				}
				else if(angle >=157.5) // angle 0
				{
					float m1 = getPixel(data, width, height, col, row+1);
					float m2 = getPixel(data, width, height, col, row-1);
					if(m0 < m1 || m0 < m2)
					{
						magnitudes[index] = 0;
					}
				}
			}
		}
		// 寻找最大与最小值
		float min = 255;
		float max = 0;
		for(int i=0; i<magnitudes.length; i++)
		{
			if(magnitudes[i] == 0) continue;
			min = Math.min(min, magnitudes[i]);
			max = Math.max(max, magnitudes[i]);
		}
		System.out.println("Image Max Gradient = " + max + " Mix Gradient = " + min);

		// 通常比值为 TL : TH = 1 : 3, 根据两个阈值完成二值化边缘连接
		// 边缘连接-link edges
		Arrays.fill(data, 0);
		int offset = 0;
		for (int row = 0; row < height; row++) {
			for (int col = 0; col < width; col++) {
				if(magnitudes[offset] >= highThreshold && data[offset] == 0)
				{
					edgeLink(col, row, offset, lowThreshold);
				}
				offset++;
			}
		}
		
		// 二值化显示
		for(int i=0; i<inPixels.length; i++)
		{
			int gray = clamp((int)data[i]);
			outPixels[i] = gray > 0 ? -1 : 0xff000000;     
		}
		setRGB(dest, 0, 0, width, height, outPixels );
		return dest;
	}
	
	public int clamp(int value) {
		return value > 255 ? 255 :
			(value < 0 ? 0 : value);
	}
	
	private void edgeLink(int x1, int y1, int index, float threshold) {
		int x0 = (x1 == 0) ? x1 : x1 - 1;
		int x2 = (x1 == width - 1) ? x1 : x1 + 1;
		int y0 = y1 == 0 ? y1 : y1 - 1;
		int y2 = y1 == height -1 ? y1 : y1 + 1;
		
		data[index] = magnitudes[index];
		for (int x = x0; x <= x2; x++) {
			for (int y = y0; y <= y2; y++) {
				int i2 = x + y * width;
				if ((y != y1 || x != x1)
					&& data[i2] == 0 
					&& magnitudes[i2] >= threshold) {
					edgeLink(x, y, i2, threshold);
					return;
				}
			}
		}
	}
	
	private float getPixel(float[] input, int width, int height, int col,
			int row) {
		if(col < 0 || col >= width)
			col = 0;
		if(row < 0 || row >= height)
			row = 0;
		int index = row * width + col;
		return input[index];
	}
	
	private float hypot(float x, float y) {
		return (float) Math.hypot(x, y);
	}
	
	private int getPixel(int[] inPixels, int width, int height, int col,
			int row) {
		if(col < 0 || col >= width)
			col = 0;
		if(row < 0 || row >= height)
			row = 0;
		int index = row * width + col;
		return inPixels[index];
	}
	
	private float gaussian(float x, float y, float sigma) {
		float xDistance = x*x;
		float yDistance = y*y;
		float sigma22 = 2*sigma*sigma;
		float sigma22PI = (float)Math.PI * sigma22;
		return (float)Math.exp(-(xDistance + yDistance)/sigma22)/sigma22PI;
	}

}
转载请务必注明出自本博客-gloomyfish

2018-10-02 23:06:53 qq_36554582 阅读数 10225
  • OpenCV3.2 Java图像处理视频学习教程

    OpenCV3.2 Java图像处理视频培训课程:基于OpenCV新版本3.2.0详细讲述Java OpenCV图像处理部分内容,包括Mat对象使用、图像读写、 基于常用核心API讲述基本原理、使用方法、参数、代码演示、图像处理思路与流程讲授。主要内容包括opencv像素操作、滤波、边缘提取、直线与圆检测、形态学操作与分水岭、图像金子塔融合重建、多尺度模板匹配、opencv人脸检测、OpenCV跟Tomcat使用实现服务器端图像处理服务。

    4253 人正在学习 去看看 贾志刚

MATLAB中有几种算法可以对图像进行边缘提取,其中一种就是edge算法,这个edge算法中有好几个算子,每一个算子分别对应着一种边缘提取的原理,接下来就来看一下几种方法的异同

%读取一张图片,并显示
original_picture=imread('D:\SoftWare\matlab2016a\Project\Picture\cat.jpg');
Pic2=im2bw(original_picture,thresh);
figure(1)
subplot(2,2,1);
imshow(original_picture);
title('原始RGB图像')
subplot(222)
imshow(Pic2)
title('二值化图像')

%用edge算法对二值化图像进行边缘提取
PicEdge1=edge(Pic2,'log');
subplot(223);
imshow(PicEdge1);
title('log算子')

PicEdge2 = edge(Pic2,'canny');
subplot(224);
imshow(PicEdge2);
title('canny算子');

PicEdge3=edge(Pic2,'sobel');
figure(2)
subplot(221)
imshow(PicEdge3);
title('sobel算子')

PicEdge4=edge(Pic2,'prewitt');
subplot(222)
imshow(PicEdge4);
title('sprewitt算子')

PicEdge5=edge(Pic2,'zerocross');
subplot(223)
imshow(PicEdge5);
title('zerocross算子')

PicEdge6=edge(Pic2,'roberts');
subplot(224)
imshow(PicEdge6);
title('roberts算子')



虽然我们从提取的结果来看,可能他们的差别不是很明显,但是这几个算子的基本原理还是有区别的,另外由于我采用的原始图片可能图中的猫和背景颜色有的部分很相似,所以会导致有些猫的边缘没有被提取出来,以后还需改善。

2017-04-19 19:54:00 The_star_is_at 阅读数 8510
  • OpenCV3.2 Java图像处理视频学习教程

    OpenCV3.2 Java图像处理视频培训课程:基于OpenCV新版本3.2.0详细讲述Java OpenCV图像处理部分内容,包括Mat对象使用、图像读写、 基于常用核心API讲述基本原理、使用方法、参数、代码演示、图像处理思路与流程讲授。主要内容包括opencv像素操作、滤波、边缘提取、直线与圆检测、形态学操作与分水岭、图像金子塔融合重建、多尺度模板匹配、opencv人脸检测、OpenCV跟Tomcat使用实现服务器端图像处理服务。

    4253 人正在学习 去看看 贾志刚

实验三   图像轮廓提取与边缘检测


一、实验目的:

理解并掌握对二值图像进行轮廓提取的相关算法(比如,掏空内部点法),以及用于图像边缘检测和提取的典型微分算子(梯度算子和拉普拉斯算子)。

二、实验环境:

计算机、Windows XP操作系统,Matlab7.0

二、实验内容:

1、根据掏空内部点算法,运用Matlab编程实现二值图像的轮廓提取。

%以下适用于黑色背景白色前景的二值图像轮廓提取(以二值图像circles为例)

BW=imread('circles.png');     %二值图像circlesuint80黑,255

subplot(1,2,1);  imshow(BW);  title('二值图像');

[M, N]=size(BW);     %M行,N

BW_Buffer=BW;

for i=2: M-1

for j=2: N-1

if (BW(i, j)==255 & BW(i-1, j)==255 & BW(i+1, j)==255 & BW(i, j-1)==255 & BW(i, j+1)==255 & BW(i-1, j-1)==255 & BW(i-1, j+1)==255 & BW(i+1, j-1)==255 & BW(i+1, j+1)==255)    %说明BW(i, j)是前景中的一个内部点

       BW_Buffer(i, j)=0;    %掏空该内部点,即将该内部点置成与背景相同灰度

end

end

end

subplot(1,2,2);  imshow(BW_Buffer);  title('提取轮廓');


%以下适用于白色背景黑色前景的二值图像轮廓提取(以二值图像source为例)

BW=imread('source.bmp');      %二值图像sourceuint80黑,255

subplot(1,2,1);  imshow(BW);  title('二值图像');

[M, N]=size(BW);   %M行,N

BW_Buffer=BW;

for i=2: M-1

for j=2: N-1

if (BW(i, j)==0 & BW(i-1, j)==0 & BW(i+1, j)==0 & BW(i, j-1)==0 & BW(i, j+1)==0 & BW(i-1, j-1)==0 & BW(i-1, j+1)==0 & BW(i+1, j-1)==0 & BW(i+1, j+1)==0)      %说明BW(i, j)是前景中的一个内部点

       BW_Buffer(i, j)=255;   %掏空该内部点,即将该内部点置成与背景相同灰度

end

end

end

subplot(1,2,2);  imshow(BW_Buffer);  title('提取轮廓');

注意:使用掏空内部点的方法来提取二值图像的轮廓时,不能直接在原始二值图像矩阵上判断一个点掏空一个点,否则对前面像素的掏空操作会影响到对后面像素的判断结果。

解决方法:创建原始二值图像矩阵的副本(即图像矩阵BW_Buffer),在原始二值图像矩阵上执行判断操作,即依次判断每个像素点是否为前景中的内部点,如果是,则在图像矩阵BW_Buffer上执行掏空内部点的操作。

2、以灰度图像ricecameraman为例,利用Matlab图像处理工具箱中的edge函数,分别使用Roberts 算子、Sobel算子、Prewitt 算子对其进行边缘检测。

(1)函数格式: BW = edge(I, 'method', thresh)

(2)格式说明:edge函数输入灰度图像矩阵I,输出二值图像矩阵BW;参数'method'用于指定所使用的边缘检测算子,可以是'roberts''sobel''prewitt''log''canny';参数thresh用于指定梯度门限值(也称梯度阈值),图像中梯度值大于等于门限值thresh的像素用白色(1)表示,说明这些地方对应边缘,梯度值小于门限值thresh的像素用黑色(0)表示(edge function will ignore all edges that are not stronger than thresh)。若不指定参数thresh,则edge函数会自动选择阈值。所以edge函数最终将原始灰度图像中的边缘和背景用二值图像的形式展现出来,以突出边缘的位置,达到边缘检测的目的。

(3)程序如下:

I=imread('rice.png');

subplot(2,2,1);  imshow(I);  title('原始图像');

[BW1,thresh1]=edge(I,'roberts');  %进行Roberts算子边缘检测并返回门限值

[BW2,thresh2]=edge(I,'sobel');    %进行Sobel算子边缘检测并返回门限值

[BW3,thresh3]=edge(I,'prewitt');  %进行Prewitt算子边缘检测并返回门限值

subplot(2,2,2);  imshow(BW1);  title('Roberts算子边缘检测结果');

subplot(2,2,3);  imshow(BW2);  title('Sobel算子边缘检测结果');

subplot(2,2,4);  imshow(BW3);  title('Prewitt算子边缘检测结果');

若向原始图像中加入随机噪声(比如高斯噪声),之后再对噪声图像分别运用Roberts 算子、Sobel算子、Prewitt 算子、Log算子(高斯-拉普拉斯算子)进行边缘检测,观察检测结果,试比较4种边缘检测算子的抗噪声干扰能力。

I=imread('rice.png');

subplot(2,3,1);  imshow(I);  title('原始图像');

G=imnoise(I, 'gaussian');  %向原始图像中加入高斯噪声

subplot(2,3,2);  imshow(G);  title('噪声图像');

BW1=edge(G, 'roberts');  %进行Roberts算子边缘检测

BW2=edge(G, 'sobel');    %进行Sobel算子边缘检测

BW3=edge(G, 'prewitt');  %进行Prewitt算子边缘检测

BW4=edge(G, 'log');      %进行Log算子边缘检测

subplot(2,3,3);  imshow(BW1);  title('Roberts算子边缘检测结果');

subplot(2,3,4);  imshow(BW2);  title('Sobel算子边缘检测结果');

subplot(2,3,5);  imshow(BW3);  title('Prewitt算子边缘检测结果');

subplot(2,3,6);  imshow(BW4);  title('Log算子边缘检测结果');

2019-05-03 21:02:44 tyfwin 阅读数 2879
  • OpenCV3.2 Java图像处理视频学习教程

    OpenCV3.2 Java图像处理视频培训课程:基于OpenCV新版本3.2.0详细讲述Java OpenCV图像处理部分内容,包括Mat对象使用、图像读写、 基于常用核心API讲述基本原理、使用方法、参数、代码演示、图像处理思路与流程讲授。主要内容包括opencv像素操作、滤波、边缘提取、直线与圆检测、形态学操作与分水岭、图像金子塔融合重建、多尺度模板匹配、opencv人脸检测、OpenCV跟Tomcat使用实现服务器端图像处理服务。

    4253 人正在学习 去看看 贾志刚

边缘提取算子

一阶:  Roberts算子、Sobel算子、Prewitt算子、Kirsch算子、Robinson算子

二阶: Laplacian算子、Canny算子、Marr-Hildreth(LoG算子)

 

Roberts 算子

在(i+1/2,j+1/2)处差分

转化为模板即为

Roberts算子,又称罗伯茨算子,是一种最简单的算子,是一种利用局部差分算子寻找边缘的算子。他采用对角线方向相邻两象素之差近似梯度幅值检测边缘。检测垂直边缘的效果好于斜向边缘,定位精度高,对噪声敏感,无法抑制噪声的影响。

Roberts算子检测方法对具有陡峭的低噪声的图像处理效果较好,但是提取边缘的结果是边缘比较粗,因此边缘的定位不是很准确。

 

Sobel算子

Sobel算法是一个离散的一阶差分算子,用来计算图像亮度函数的一阶梯度之近似值。在图像的任何一点使用此算子,将会产生该点对应的梯度矢量或是其法矢量。

 

Prewitt算子

Prewitt算子是一种一阶微分算子边缘检测,利用像素点上下、左右邻点的灰度差,在边缘处达到极值检测边缘,去掉部分伪边缘,对噪声具有平滑作用。

其原理是在图像空间利用两个方向模板与图像进行邻域卷积来完成的,这两个方向模板一个检测水平边缘,一个检测垂直边缘。

Prewitt算子检测方法对灰度渐变和噪声较多的图像处理效果较好。但边缘较宽,而且间断点多。

对数字图像f(x,y),Prewitt算子的定义如下:

G(i)={[f(i-1,j-1)+f(i-1,j)+f(i-1,j+1)]-[f(i+1,j-1)+f(i+1,j)+f(i+1,j+1)]}

G(j)={[f(i-1,j+1)+f(i,j+1)+f(i+1,j+1)]-[f(i-1,j-1)+f(i,j-1)+f(i+1,j-1)]}

则P(i,j)=max[G(i),G(j)]或 P(i,j)=G(i)+G(j)

 

Kirsch算子

Kirsch算子是R.Kirsch提出来一种边缘检测算法,它采用8个3*3的模板对图像进行卷积,这8个模板代表8个方向,并取最大值作为图像的边缘输出,8个模板如下,它在保持细节和抗噪声方面都有较好的效果。

 

 

Robinson算子

规则同上,也是8个模板。

 

Laplacian算子

Laplacian算子法对噪声比较敏感,所以很少用该算子检测边缘,而是用来判断边缘像素视为与图像的明区还是暗区。

Laplacian 算子是n维欧几里德空间中的一个二阶微分算子,定义为梯度grad的散度div。可使用运算模板来运算这定理定律。

函数的拉普拉斯算子也是该数的黑塞矩阵的迹,可以证明,它具有各向同性,即与坐标轴方向无关,坐标轴旋转后梯度结果不变。

Laplacian 算子对噪声比较敏感,所以图像一般先经过平滑处理,因为平滑处理也是用模板进行的,所以,通常的分割算法都是把Laplacian 算子和平滑算子结合起来生成一个新的模板。

如果邻域系统是4 邻域,Laplacian 算子的模板为:

如果邻域系统是8 邻域,Laplacian 算子的模板为:

 

Canny算子

Canny方法不容易受噪声干扰,能够检测到真正的弱边缘。

优点在于,使用两种不同的阈值分别检测强边缘和弱边缘,并且当弱边缘和强边缘相连时,才将弱边缘包含在输出图像中。

该算子效果较好,但是它实现起来较为麻烦,Canny算子是一个具有滤波,增强,检测的多阶段的优化算子,在进行处理前,Canny算子先利用高斯平滑滤波器来平滑图像以除去噪声,Canny分割算法采用一阶偏导的有限差分来计算梯度幅值和方向,在处理过程中,Canny算子还将经过一个非极大值抑制的过程,最后Canny算子还采用两个阈值来连接边缘。

Canny边缘检测算法

step1: 用高斯滤波器平滑图象;

step2: 用一阶偏导的有限差分来计算梯度的幅值和方向;

step3: 对梯度幅值进行非极大值抑制

step4: 用双阈值算法检测和连接边缘

参考:

https://blog.csdn.net/fengye2two/article/details/79190759

https://blog.csdn.net/jiachen0212/article/details/80078685

https://www.cnblogs.com/zhuifeng-mayi/p/9563947.html

https://www.cnblogs.com/techyan1990/p/7291771.html

 

Marr-HildrethLoG算子

LoG可以看成是一个高斯模板的拉普拉斯变换   Laplace of Gaussian

(图像的高斯拉普拉斯(LoG),可利用差分高斯(DoG)近似)

参考:

好:LOG高斯-拉普拉斯算子 https://blog.csdn.net/touch_dream/article/details/62237018 

DoG算子和LoG算子 https://blog.csdn.net/dreamguard/article/details/83988814

好:LoG算子 https://blog.csdn.net/songzitea/article/details/12851079

LOG边缘检测--Marr-Hildreth边缘检测算法 https://blog.csdn.net/bettyshasha/article/details/51757185  

好:SIFT特征点检测,特征点描述,特征点匹配理解 https://blog.csdn.net/ds0529/article/details/79617619

matlab代码

 

matlab自带edge函数

matlab中edge函数提供了很多算子,函数原型和算子包括:

matlab自带edge函数 使用代码:

clear;clc;

%读取图像
I=imread('Lena512.png');
%I=rgb2gray(I);
figure;imshow(I,[]);title('Original Image');

%sobel----------------------
sobelBW = edge(I,'sobel',0.06);  %可以省去0.06,系统会默认。
figure;imshow(sobelBW,[]);title('Sobel 0.06')

%roberts----------------------
robertsBW=edge(I,'roberts');
figure;imshow(robertsBW);title('Roberts Edge');

%prewitt----------------------
prewittBW=edge(I,'prewitt');
figure;imshow(prewittBW);title('Prewitt Edge');

%Laplacian----------------------
logBW=edge(I,'log');
figure;imshow(logBW);title('Laplasian of Gaussian Edge');

%canny----------------------
cannyBW=edge(I,'canny');
figure;imshow(cannyBW);title('Canny Edge');

%高斯滤波 再 canny ----------------------
h=fspecial('gaussian',5);%高斯滤波
I2=imfilter(I,h,'replicate');
GcannyBW=edge(I2,'canny');%高斯滤波后使用Canny算子进行边缘检测
figure;imshow(GcannyBW);title('Gaussian & Canny Edge');

 

非matlasb自带 :梯度法、roberts算子、prewitt算子、sobel算子、krisch算子、LoG算子

%% 非matlasb自带 :梯度法、roberts算子、prewitt算子、sobel算子、krisch算子、LoG算子
clear;clc;close all
f=imread('Lena512.png');
T=20;%阈值
[m,n]=size(f);
%------梯度法-------
f_g=zeros(m,n);
for i=2:m-1
    for j=2:n-1
        f_g(i,j)=abs(f(i+1,j)-f(i,j))+abs(f(i,j+1)-f(i,j));
        if f_g(i,j)<T
            f_g(i,j)=0;
        else
            f_g(i,j)=255;
        end
    end
end
figure(1);
subplot(2,3,1);imshow(uint8(f_g));title('梯度法');
%------roberts算子-------
f_r=zeros(m,n);
for i=2:m-1
    for j=2:n-1
        f_r(i,j)=abs(f(i+1,j+1)-f(i,j))+abs(f(i,j+1)-f(i+1,j));
        if f_r(i,j)<T
            f_r(i,j)=0;
        else
            f_r(i,j)=255;
        end
    end
end
%f_r=imbinarize(imfilter(f,r),T);
subplot(2,3,2);imshow(uint8(f_r));title('Roberts算子');
 
%------prewitt算子-------
f_p=zeros(m,n);
for i=2:m-1
    for j=2:n-1
        f_p(i,j)=abs(f(i-1,j-1)+f(i,j-1)+f(i+1,j-1)-f(i-1,j+1)-f(i,j+1)-f(i+1,j+1))+abs(f(i+1,j-1)+f(i+1,j)+f(i+1,j+1)-f(i-1,j-1)-f(i-1,j)-f(i-1,j+1));
        if f_p(i,j)<15
            f_p(i,j)=0;
        else
            f_p(i,j)=255;
        end
    end
end
subplot(2,3,3);imshow(uint8(f_p));title('Prewitt算子');
 
%------sobel算子-------
f_s=zeros(m,n);
for i=2:m-1
    for j=2:n-1
        f_s(i,j)=abs(f(i-1,j-1)+2*f(i,j-1)+f(i+1,j-1)-f(i-1,j+1)-2*f(i,j+1)-f(i+1,j+1))+abs(f(i+1,j-1)+2*f(i+1,j)+f(i+1,j+1)-f(i-1,j-1)-2*f(i-1,j)-f(i-1,j+1));
        if f_s(i,j)<T
            f_s(i,j)=0;
        else
            f_s(i,j)=255;
        end
    end
end
subplot(2,3,4);imshow(uint8(f_s));title('Sobel算子');
 
%------krisch算子-------
k(:,:,1)=[-3,-3,-3;
    -3,0,5;
    -3,5,5];
k(:,:,2)=[-3,-3,5;
    -3,0,5;
    -3,-3,5];
k(:,:,3)=[-3,5,5;
    -3,0,5;
    -3,-3,-3];
k(:,:,4)=[-3,-3,-3;
    -3,0,-3;
    5,5,5];
k(:,:,5)=[5,5,5;
    -3,0,-3;
    -3,-3,-3];
k(:,:,6)=[-3,-3,-3;
    5,0,-3;
    5,5,-3];
k(:,:,7)=[5,-3,-3;
    5,0,-3;
    5,-3,-3];
k(:,:,8)=[5,5,-3;
    5,0,-3;
    -3,-3,-3];
kk=zeros(size(f));
I=double(f);
for i=1:8
    f_k(:,:,i)=conv2(I,k(:,:,i),'same');
    kk=max(kk,f_k(:,:,i));
end
f_kk=imbinarize(kk,600);
subplot(2,3,5);imshow(f_kk);title('Krisch算子');
 
%------LoG算子-------
log1=[0 0 -1 0 0;
    0 -1 -2 -1 0;
    -1 -2 16 -2 -1;
    0 -1 -2 -1 0;
    0 0 -1 0 0];
 
f_l=conv2(f,log1,'same');
f_ll=imbinarize(abs(f_l),300);
subplot(2,3,6);imshow(f_ll);title('LoG算子');

 

非matlab自带   Kirsch算子

(参考https://blog.csdn.net/u014485485/article/details/78339420

clear;clc;close all

%读取图像
bw1=imread('Lena512.png');
%bw1=rgb2gray(bw1);
figure;imshow(bw1,[]);title('Original Image');

%均值滤波
bw2=filter2(fspecial('average',3),bw1);
%高斯滤波
bw3=filter2(fspecial('gaussian'),bw2);
%利用小波变换对图象进行降噪处理
[thr,sorh,keepapp]=ddencmp('den','wv',bw3);     %获得除噪的缺省参数
bw4=wdencmp('gbl',bw3,'sym4',2,thr,sorh,keepapp);%图象进行降噪处理

%提取图象边缘
t=[0.8 1.0 1.5 2.0 2.5].*10^5 ;     %设定阈值
bw5=double(bw4);            
[m,n]=size(bw5);             
g=zeros(m,n); 
d=zeros(1,8);
%利用Kirsch算子进行边缘提取
for i=2:m-1
   for j=2:n-1
       d(1) =(5*bw5(i-1,j-1)+5*bw5(i-1,j)+5*bw5(i-1,j+1)-3*bw5(i,j-1)-3*bw5(i,j+1)-3*bw5(i+1,j-1)-3*bw5(i+1,j)-3*bw5(i+1,j+1))^2; 
       d(2) =((-3)*bw5(i-1,j-1)+5*bw5(i-1,j)+5*bw5(i-1,j+1)-3*bw5(i,j-1)+5*bw5(i,j+1)-3*bw5(i+1,j-1)-3*bw5(i+1,j)-3*bw5(i+1,j+1))^2; 
       d(3) =((-3)*bw5(i-1,j-1)-3*bw5(i-1,j)+5*bw5(i-1,j+1)-3*bw5(i,j-1)+5*bw5(i,j+1)-3*bw5(i+1,j-1)-3*bw5(i+1,j)+5*bw5(i+1,j+1))^2; 
       d(4) =((-3)*bw5(i-1,j-1)-3*bw5(i-1,j)-3*bw5(i-1,j+1)-3*bw5(i,j-1)+5*bw5(i,j+1)-3*bw5(i+1,j-1)+5*bw5(i+1,j)+5*bw5(i+1,j+1))^2; 
       d(5) =((-3)*bw5(i-1,j-1)-3*bw5(i-1,j)-3*bw5(i-1,j+1)-3*bw5(i,j-1)-3*bw5(i,j+1)+5*bw5(i+1,j-1)+5*bw5(i+1,j)+5*bw5(i+1,j+1))^2; 
       d(6) =((-3)*bw5(i-1,j-1)-3*bw5(i-1,j)-3*bw5(i-1,j+1)+5*bw5(i,j-1)-3*bw5(i,j+1)+5*bw5(i+1,j-1)+5*bw5(i+1,j)-3*bw5(i+1,j+1))^2; 
       d(7) =(5*bw5(i-1,j-1)-3*bw5(i-1,j)-3*bw5(i-1,j+1)+5*bw5(i,j-1)-3*bw5(i,j+1)+5*bw5(i+1,j-1)-3*bw5(i+1,j)-3*bw5(i+1,j+1))^2; 
       d(8) =(5*bw5(i-1,j-1)+5*bw5(i-1,j)-3*bw5(i-1,j+1)+5*bw5(i,j-1)-3*bw5(i,j+1)-3*bw5(i+1,j-1)-3*bw5(i+1,j)-3*bw5(i+1,j+1))^2;      
       g(i,j) = max(d);
    end
end 

%显示边缘提取后的图象
figure(2)
for k=1:5
    for i=1:m
        for j=1:n
            if g(i,j)>t(k)
                bw5(i,j)=255;           
            else
                bw5(i,j)=0;
            end
        end
    end
    subplot(1,5,k)
    imshow(bw5,[])
    title(['Kirsch' '   ' num2str(t(k))])
end

最后我之前自己写的梯度法:

I=imread('1.jpg');
I=rgb2gray(I);
T=10;%阈值
[x,y]=gradient(double(I));   %获取梯度
t=sqrt(x.^2+y.^2);  
I(t>=T)=255;           %梯度提取边缘 画黑
I(t<T)=0;
I=uint8(I);
imshow(I);

imwrite(I,'my开方梯度.jpg');  %保存图像到当前目录

本人水平有限,难免有错误疏漏,还望大佬多多批评指正。

图像的边缘提取

阅读数 6116

没有更多推荐了,返回首页