深度剖析图像处理—边缘检测

什么是边缘检测

边缘检测(Edge Detection)就是提取图像中的边缘点(Edge Point)。边缘点是与周围像素相比灰度值有阶跃变化或屋顶状变化的像素。边缘常存在于目标与背景之间、目标与目标之间、目标与其影子之间。

​ 在图像处理和图像分析中，经常要用到边缘(Edge)、边界(Boundary)、轮廓(Contour)等术语。一般来说，边缘指的是边缘点，它不能被称为边缘线。边界指的是图像中不同区域之间的分界线，比如不同灰度、不同颜色的区域之间的分界线，它是线而不是点，可以被称为边界线。轮廓一般是指目标的轮廓，目标就是语义明确的区域，轮廓一般是在二值图像中围绕白色区域的闭合曲线。

​ 边缘检测算法和边界线检测算法一般作用于灰度图像，对于二值图像进行边缘检测是没有意义的。轮廓一定是闭合的，但边界线不一定闭合，比如道路区域与道边植被的边界线；边缘点最多是断断续续的线段，不保证连续，更不保证闭合。掌握边缘、边界、轮廓的准确术语是非常必要的。

边缘类型

边缘检测是一种邻域运算，即一个像素是否是边缘点是由其所处的邻域决定的。在一定大小的邻域内，边缘分为阶跃边缘(step edge)和屋顶状边缘(roof edge)两种类型。下面以一维信号为例，分析这两种不同类型的边缘的导数特征

image-20240421200943004

求导与差分

在边缘检测中，导数的计算通常采用两种方法：

  • 将邻域从离散空间变换到连续空间，得到解析描述，然后进行求导操作**。**

    具体做法是，先将**邻域按照一定的数学模型(曲线拟合、曲面拟合)得到其在连续空间中的解析描述，然后对此解析描述进行求导，得到边缘点。解析描述求导得到的导数位置是有小数位的，比如在位置4.17处取得导数最大值，即边缘点的位置是在4.17而不是像素的整数坐标4。这样得到的边缘点位置精度能够小于1个像素，因此又将此方法称之为亚像素(sub pixel)边缘检测。在已知数学模型的指导下，目前工业界做到的边缘检测最高精度为1/50个像素。**亚像素边缘检测能够在大大节省硬件成本的同时，得到高的边缘检测精度，是图像测量中的常用方法。

  • 直接用差分(difference)代替求导。导数的公式见下式，如果令dx=1，即得到是差分描述，式（4-2）是差分描述的x方向偏导数，式(4-3)是差分描述的y方向偏导数。

image-20240421202909028

image-20240421203632870

边缘强度与边缘方向

导数是有大小也有方向的，因此边缘也有强弱与方向，分别叫做边缘强度(edge intensity)和边缘方向(edge direction)，边缘强度即边缘的幅值(magnitude)。用M(x,y)代表边缘的强度，θ(x,y)代表边缘的方向，有：

image-20240421204859267

实例

我们该如何提取这张图片的边缘呢?

test

首先当然是需要我们写一个函数来把24位彩色图像转化为8位灰度值图像

//24位彩色图像转8位灰度值
//rgbImage原始图像
//grayImage输出灰度图像
//width,height图片的宽和高
void convertToGray(uint8_t* rgbImage, uint8_t* grayImage, int width, int height)
{for (int y = 0; y < height; y++) {for (int x = 0; x < width; x++) {// 获取当前像素的 RGB 分量uint8_t r = rgbImage[3 * (y * width + x) + 0];uint8_t g = rgbImage[3 * (y * width + x) + 1];uint8_t b = rgbImage[3 * (y * width + x) + 2];// 计算灰度值（常用的加权平均法）// 这里使用的加权系数是常见的：R: 0.299, G: 0.587, B: 0.114uint8_t gray = (uint8_t)(0.299 * r + 0.587 * g + 0.114 * b);// 将灰度值写入灰度图像数组grayImage[y * width + x] = gray;}}
}

我们来看看效果

image-20240422174937601

那接下来就需要用边缘检测来提取边缘了

一阶微分算子

根据边缘类型及其导数特征，可以设计不同的检测算法。下面讲述几种常用的边缘检测算法，习惯上称为边缘检测算子(Operator)。当使用差分时，一般写成模板的表示形式。

对于阶跃边缘而言，边缘点处的导数特征是“一阶导数取极值”。若边缘点处的一阶导数为正值，则其为最大值；反之，则为最小值，即在边缘点处的导数绝对值最大。

基于一阶导数的边缘检测算子称为一阶微分算子，常用的一阶微分算子有梯度算子、罗伯特算子、索贝尔算子、Prewitt算子，Robinson算子、Kirsch算子等。

梯度算子

梯度的本意是一个向量（矢量），表示某一函数在该点处的方向导数沿着该方向取得最大值，即函数在该点处沿着该方向（此梯度的方向）变化最快，变化率最大（为该梯度的模）。梯度的含义和边缘点是一致的，因此产生了边缘检测的梯度算子(Gradient Operator)

image-20240422172912964

写成模板的情况就是如下图所示：

image-20240422172958535

但是上面我们只是给出了像素(x,y)的边缘强度，称为梯度值；但是它是不是边缘点，还需要一定的约束条件，比如，设定当Gradient(x,y)≥threshold时，像素(x,y)才是边缘点，threshold称为阈值。

image-20240422173629195

这就是用梯度算子计算的结果，我们上代码看看

void RmwGradientGryImg(uint8_t* pGryImg, int width, int height, uint8_t* pGrdImg)
{uint8_t* pGry, * pGrd;int dx, dy;int x, y;for (y = 0, pGry = pGryImg, pGrd = pGrdImg; y < height - 1; y++){for (x = 0; x < width - 1; x++, pGry++){dx = *pGry - *(pGry + 1);dy = *pGry - *(pGry + width);int gradient = (int)(sqrt(dx * dx * 1.0 + dy * dy));*pGrd++ = (gradient > 255) ? 255 : gradient;}*pGrd++ = 0; //尾列不做,边缘强度赋0pGry++;}memset(pGrd, 0, width); //尾行不做,边缘强度赋0
}

这里简单说一下这个函数

  1. 函数参数：
    • pGryImg：输入的灰度图像数据指针。
    • width：图像的宽度。
    • height：图像的高度。
    • pGrdImg：输出的梯度图像数据指针。
  2. 双层循环：
    • 外层循环 for (y = 0, pGry = pGryImg, pGrd = pGrdImg; y<height-1; y++) 遍历图像的每一行，pGry 和 pGrd 分别指向当前行的灰度图像数据和梯度图像数据。
    • 内层循环 for (x = 0; x<width-1; x++, pGry++) 遍历当前行的每个像素，pGry 指向当前像素的灰度值。
  3. 计算梯度：
    • 梯度的计算采用的是简单的基于像素差值的方法，分别计算水平方向和垂直方向的梯度。
    • dx = *pGry-*(pGry+1)：水平方向的梯度，计算当前像素和右侧像素的灰度差值。
    • dy = *pGry-*(pGry+width)：垂直方向的梯度，计算当前像素和下方像素的灰度差值。
    • sqrt(dx*dx*1.0+dy*dy)：使用欧式距离公式计算梯度幅值。
    • min(255, ...)：确保梯度幅值不超过255，限定在0到255之间。
  4. 内存清零：
    • memset(pGrd, 0, width)：清零输出图像的最后一行，因为最后一行的梯度值未计算。
  5. 返回：
    • 函数返回，处理完成。

看看实际运行的效果是什么样子的

image-20240422180246550

效果并不是很好，我们把原始图像做一下图像增强（灰度值均衡化）试试

image-20240422181123926

换一种图像增强的方法，加一下反相试试

先简单写一下反相的函数

void invertImage(uint8_t *image, int width, int height) {for (int i = 0; i < width * height; i++) {image[i] = 255 - image[i];}
}

这次使用了线性拉伸并简单处理了一下参数

image-20240422181718341

效果还算不错，那如果加一下阈值呢？就像上面说的，我们重新写一下这个函数

//梯度算子加阈值
//pGryImg：输入的灰度图像数据指针。
//width：图像的宽度。
//height：图像的高度。
//pGrdImg：输出的梯度图像数据指针
void RmwGradientGryImgPlus(uint8_t* pGryImg, int width, int height, uint8_t* pGrdImg, int threshold)
{uint8_t* pGry, * pGrd;int dx, dy;int x, y;for (y = 0, pGry = pGryImg, pGrd = pGrdImg; y < height - 1; y++){for (x = 0; x < width - 1; x++, pGry++){dx = *pGry - *(pGry + 1);dy = *pGry - *(pGry + width);int gradient = (int)(sqrt(dx * dx * 1.0 + dy * dy));*(pGrd++) = (gradient > threshold) ? min(255, gradient) : 0;}*(pGrd++) = 0; //尾列不做,边缘强度赋0pGry++;}memset(pGrd, 0, width); //尾行不做,边缘强度赋0
}

看一下效果

image-20240422182112492

画面确实变得更干净了

罗伯特算子

梯度算子的计算只涉及到了3个像素，只在水平和垂直方向做差分。罗伯特算子（Roberts Operator)给出了一个4个像素之间进行运算的算子，分别在两个对角线方向做差分。

image-20240422182433397

​ 由于对角线上2个像素之间的距离为√2，所以罗伯特算子的∆x和∆y采用对角线差分后，不再采用√(∆_x2+∆_y2 )，其描述见下式。罗伯特算子取∆x绝对值与∆y绝对值中的最大值。

image-20240422182518316

我们来看演示

image-20240422182547389

那么罗伯特算子好在哪呢？

罗伯特算子去掉了梯度算子的开方运算，计算复杂度也降低了不少。

void RmwRobertsGryImg(uint8_t *pGryImg, int width, int height, uint8_t *pRbtImg)
{uint8_t *pGry, *pRbt;int dx, dy;int x, y;for (y = 0, pGry = pGryImg, pRbt = pRbtImg; y < height - 1; y++){for (x = 0; x < width - 1; x++, pGry++){dx = *pGry - *(pGry + width + 1);dy = *(pGry + 1) - *(pGry + width);*pRbt++ = (uint8_t)(dx > dy ? dx : dy); // 使用三目运算符选择较大的值}*pRbt++ = 0; // 尾列不做, 边缘强度赋0pGry++;}memset(pRbt, 0, width); // 尾行不做, 边缘强度赋0
}

我们来看结果

image-20240422192132132

为什么看起来的效果好像还没有梯度算子的结果好呢

是因为图像太复杂没有滤波，那么在边缘计算前能不能先滤波？

我们来看下一个算子

索贝尔算子

在一幅噪声较大的图像中，如果不进行图像平滑就进行边缘检测，必然会在边缘图像中产生噪声干扰。

因此，索贝尔算子(Sobel Operator)中，在求∆x和∆y前，先进行滤波。在求∆x前，先执行如下图的所示的高斯均值滤波；在求∆y前，先执行如下图的所示的高斯均值滤波。

image-20240422184902956

另外，索贝尔算子进一步拉大进行差分的2个像素之间的距离，∆x和∆y采用如下的模板形式

image-20240422184942115

我们来看示意图

image-20240422185037120

在图像1中，红色方块代表当前像素(x,y)，先执行图4-8(a)所示的高斯滤波，用D ̅、E ̅代表滤波后的值，则得到：

image-20240422185102632

执行∆x模板，则有∆x=D ̅-E ̅=(A+2D+F)-(C+2E+H)

image-20240422185448932

另外，索贝尔算子在对∆x和∆y的使用上，采用了它们的绝对值相加的形式

image-20240422185539273

看效果：

image-20240422185603424

对原始灰度图像执行索贝尔算子得到的结果，图中虚线框所示的边缘变成了双线宽。

我们来写代码

//索贝尔算子
void RmwSobelGryImg(uint8_t* pGryImg, int width, int height, uint8_t* pSbImg)
{uint8_t* pGry, * pSb;int dx, dy;int x, y;memset(pSbImg, 0, width); // 首行不做, 边缘强度赋0for (y = 1, pGry = pGryImg + width, pSb = pSbImg + width; y < height - 1; y++){*pSb++ = 0; // 首列不做, 边缘强度赋0pGry++;for (x = 1; x < width - 1; x++, pGry++){// 求dxdx = *(pGry - 1 - width) + (*(pGry - 1) * 2) + *(pGry - 1 + width);dx -= *(pGry + 1 - width) + (*(pGry + 1) * 2) + *(pGry + 1 + width);// 求dydy = *(pGry - width - 1) + (*(pGry - width) * 2) + *(pGry - width + 1);dy -= *(pGry + width - 1) + (*(pGry + width) * 2) + *(pGry + width + 1);// 结果*pSb++ = (uint8_t)min(255, abs(dx) + abs(dy));}*pSb++ = 0; // 尾列不做, 边缘强度赋0pGry++;}memset(pSb, 0, width); // 尾行不做, 边缘强度赋0
}

image-20240422191511580

梯度算子、罗伯特算子、索贝尔算子的比较

以下从4个方面对梯度算子、罗伯特算子、索贝尔算子进行比较。

•1. 偏导数∆x和∆y的求取

梯度算子在3个像素之间进行运算，只在水平和垂直方向做差分，做差分的2个像素之间的距离为1。

罗伯特算子在4个像素之间进行运算，分别在两个对角线方向做差分，做差分的2个像素之间的距离为√2。

索贝尔算子在8个像素之间进行运算，只在水平和垂直方向做差分，做差分的2个像素之间的距离为2。

•2. 是否“先平滑后求导”

索贝尔算子在差分之前，进行了加权均值滤波对图像进行平滑（加权函数采用了高斯函数），因此索贝尔具有滤除噪声的效果。梯度算子和罗伯特算子都没有进行平滑。“先平滑后求导”是边缘检测的通用策略，一般在执行梯度算子和罗伯特算子前是需要使用另外的步骤做图像平滑，索贝尔算子则是把平滑写到了算子中。

•3. 边缘强度的大小

按照式(4-4)边缘强度的定义，梯度算子是严格遵守的，罗伯特算子是取∆x绝对值和∆y绝对值中的最大值，索贝尔算子是取∆x与∆y的绝对值之和。而且，索贝尔算子在高斯滤波后没有除以4，所以又相当于∆x、∆y放大了4倍。对于边缘强度，罗伯特算子、梯度算子、索贝尔算子之间的数值关系大致如下：

image-20240422195609623

•4. 邻域与边缘宽度

梯度算子、罗伯特算子的计算只涉及到了2行2列，所以它们得到的边缘宽度是1个像素；索贝尔算子涉及到了3行3列，所以它得到的边缘宽度是2个像素，边缘变成了双线宽。

方向模板

若是能根据边缘的具体走向求偏导数，则边缘强度值应该会更准确。

因此在实际应用中，是先假定了有限的几个边缘方向，再对这些假定的每个边缘方向设置一个特定的模板，计算每个模板的边缘强度，从中选择最大的边缘强度作为边缘强度的结果，而且该最大边缘强度对应模板的方向就认为是边缘的方向。

常用基于方向模板的边缘检测算子有：Prewitt算子，Robinson算子、Kirsch算子。Prewitt算子使用4个方向模板，Robinson算子和Kirsch算子都使用8个方向模板。这些算子都是先采用了均值滤波，然后进行差分计算。

Prewitt算子

Prewitt算子设定了0°、45°、90°和135°，共计4种边缘方向；根据这4种边缘方向，分别设计了4个模板

image-20240422195755759

对每个像素分别计算这个4个模板的值，取绝对值最大者作为该像素的边缘强度Prewitt(x,y)。同时该最大值对应模板的方向作为该像素的边缘方向(与边缘的走向相差90°，因为显然边缘走向的法线方向上的导数最大)。若把这些模板中为“0”点(空白处)的连成一条直线，可以发现这些模板强调了水平线、135°斜线、竖直线、45°斜线的检测。Prewitt算子强调对直线的检测，对于上述走向的直线，总有一个模板的输出值最大。

void RmwPrewittGryImg(uint8_t *pGryImg, int width, int height, uint8_t *pPRTImg)
{uint8_t *pGry, *pPRT;int dx, dy, d45, d135, v1, v2;int x, y;memset(pPRTImg, 0, width); // 首行不做, 边缘强度赋0for (y = 1, pGry = pGryImg + width, pPRT = pPRTImg + width; y < height - 1; y++){*pPRT++ = 0; // 首列不做, 边缘强度赋0pGry++;for (x = 1; x < width - 1; x++, pGry++){// 求dxdx = *(pGry - 1 - width) + *(pGry - 1) + *(pGry - 1 + width);dx -= *(pGry + 1 - width) + *(pGry + 1) + *(pGry + 1 + width);// 求dydy = *(pGry - width - 1) + *(pGry - width) + *(pGry - width + 1);dy -= *(pGry + width - 1) + *(pGry + width) + *(pGry + width + 1);// 求45度d45 = *(pGry - width - 1) + *(pGry - width) + *(pGry - 1);d45 -= *(pGry + width + 1) + *(pGry + width) + *(pGry + 1);// 求135度d135 = *(pGry - width) + *(pGry - width + 1) + *(pGry + 1);d135 -= *(pGry + width - 1) + *(pGry + width) + *(pGry - 1);// 结果v1 = abs(dx) > abs(dy) ? abs(dx) : abs(dy);v2 = abs(d45) > abs(d135) ? abs(d45) : abs(d135);*pPRT++ = (uint8_t)((v1 > v2) ? ((v1 > 255) ? 255 : v1) : ((v2 > 255) ? 255 : v2));}*pPRT++ = 0; // 尾列不做, 边缘强度赋0pGry++;}memset(pPRT, 0, width); // 尾行不做, 边缘强度赋0
}

我们看一下实现效果

image-20240422200701865

Robinson算子

Robinson算子除(a)外的7个模板，都是由其上个模板顺时针旋转1个像素得到的。若是把模板中的负数值合并成一个区域，可以看出该算子强调了对角点的检测；对于各种形状的角点，总有一个模板的输出值最大。

image-20240422200914550

有兴趣大家可以查一下啊，这里就不再过多解释

二阶微分算子

一阶微分算子能够得到边缘强度，但是需要再加上一定的条件约束，比如设置个阈值，才能判定一个像素是不是边缘点。

通过对边缘类型及其导数的分析可知，阶跃边缘的导数特征除了“一阶导数取极值”外，还有“二阶导数过零点”。因此可以采用二阶导数，利用过零点得到边缘点，这样就不需要其他的条件了。

拉普拉斯算子

拉普拉斯算子(Laplacian Operator)是近似给出二阶导数的流行方法，其使用3×3的邻域，给出了4-邻接和8-邻接的邻域的2种模板

image-20240422201917270

对如图所示的原始灰度图像执行4邻域拉普拉斯算子，得到的结果如图所示，图中虚线框所示的位置上发生了过零点（导数由负数变到了正数），此处即边缘。

image-20240422202053683

​ 在拉普拉斯算子的结果图像中，可以发现过零点位置刚好就是边缘的位置。由于过零点是在像素之间，不在整数坐标上，所以在提取边缘点时，往往采取下面的策略：

当一个像素的二阶导数大于0，其邻域内有像素的二阶导数小于0或等于0，则该像素被标记为边缘点。

沈俊算子

唯一一个以国人命名的算子

沈俊教授同样提出了先滤波后求导的边缘检测方法（J. Shen and S. Castan, An optimal linear operator for step edge detection, CVGIP: Graphical Models and Image Processing, Vol. 54 No.2, Mar. 1992, pp.112 – 133），即沈俊算子(ShenJun Edge Operator)。

沈教授在阶跃边缘和可加白噪声模型下，就信噪比最大准则，证明了图像平滑的最佳滤波器是对称的指数函数，形式如下：

image-20240422202508333

显然，当a_0越大时，c2就越小，T(j,i)就越陡越窄，相当于滤波邻域就越小，压制噪声的能力就弱，图像模糊程度就越小，边缘定位的精度就越高。

在算子实现上，沈教授对图像分别按行、按列各进行两次先正方向再反方向的递推滤波实现(|j|、|i|的优点)，等价于用上述指数函数进行图像滤波；证明了滤波结果减去原始灰度值得到的差值乘以2c1〖ln〗^c2，约等于其二阶导数的值。沈俊算子的实现过程如下：

沈俊算子的实现过程
step.1 对每行从左向右进行：
𝑔1(0,𝑦)=𝑔(0,𝑦)，
𝑔1(𝑥,𝑦)=𝑔1(𝑥−1,𝑦)+𝑎_0×(𝑔(𝑥,𝑦)−𝑔1(𝑥−1,𝑦))，𝑥=1,2,⋯𝑤𝑖𝑑𝑡ℎ−1。
step.2 对每行从右向左进行：
𝑔2(𝑤𝑖𝑑𝑡ℎ−1,𝑦)=𝑔1(𝑤𝑖𝑑𝑡ℎ−1,𝑦)，
𝑔2(𝑥,𝑦)=𝑔2(𝑥+1,𝑦)+𝑎_0×(𝑔1(𝑥,𝑦)−𝑔2(𝑥+1,𝑦))，𝑥=𝑤𝑖𝑑𝑡ℎ−2, ⋯, 1,0。
step.3 对每列从上向下进行：
𝑔3(𝑥,0)=𝑔2(𝑥,0)，
𝑔3(𝑥,𝑦)=𝑔3(𝑥,𝑦−1)+𝑎_0×(𝑔2(𝑥,𝑦)−𝑔3(𝑥,𝑦−1))，𝑦=1,2,⋯ℎ𝑒𝑖𝑔ℎ𝑡−1。
step.4 对每列从下向上进行：
𝑔4(𝑥,ℎ𝑖𝑔ℎ𝑡−1)=𝑔3(𝑥,ℎ𝑖𝑔ℎ𝑡−1)，
𝑔4(𝑥,𝑦)=𝑔4(𝑥,𝑦+1)+𝑎_0×(𝑔3(𝑥,𝑦)−𝑔4(𝑥,𝑦+1))，𝑦=ℎ𝑒𝑖𝑔ℎ𝑡−2, ⋯, 1,0。
step.5 对每个像素(𝑥,𝑦)执行𝑆𝐽(𝑥,𝑦)=𝑔4(𝑥,𝑦)−𝑔(𝑥,𝑦)，得到二阶导数𝑆𝐽(𝑥,𝑦)。
step.6 对每个像素𝑆𝐽(𝑥,𝑦)进行过零点检测得到边缘点。

我们来看程序

//沈俊算子
//pGryImg 和 pTmpImg 是指向 uint8_t 类型的指针，它们分别指向原始灰度图像数据和辅助图像数据。
//width 和 height 是整型参数，表示图像的宽度和高度。
//a0 是双精度浮点型参数，表示滤波系数。
//pSJImg 是指向 uint8_t 类型的指针，它指向了输出的图像数据。
void RmwShenJunGryImg(uint8_t* pGryImg,uint8_t* pTmpImg, int width, int height, double a0, uint8_t* pSJImg)
{uint8_t* pGry, * pCur, * pSJ, * pEnd;int LUT[512], * ALUT; // a0查找表int x, y, pre, dif;// Step 1: 初始化查找表a0 = (a0 < 0.01) ? 0.01 : ((a0 > 0.99) ? 0.99 : a0); // 安全性检查// a0查找表, 进行了四舍五入ALUT = LUT + 256;for (ALUT[0] = 0, dif = 1; dif < 256; dif++){ALUT[dif] = (int)(dif * a0 + 0.5);ALUT[-dif] = (int)(-dif * a0 - 0.5);}// Step 2: 递推实现指数滤波// 按行滤波for (y = 0, pGry = pGryImg, pCur = pTmpImg; y < height; y++){// 从左向右: p1(y,x) = p1(y,x-1) + a * [p(y,x) - p1(y,x-1)]*(pCur++) = pre = *(pGry++);for (x = 1; x < width; x++, pGry++)*(pCur++) = pre = pre + ALUT[*pGry - pre];pCur--; // 回到行尾// 从右向左: p2(y,x) = p2(y,x+1) - a * [p1(y,x) - p2(y,x+1)]for (x = width - 2, pCur = pCur - 1; x >= 0; x--)*(pCur--) = pre = pre + ALUT[*pCur - pre];pCur += (width + 1); // 回到下一行的开始}// 按列滤波for (x = 0, pCur = pTmpImg; x < width; x++, pCur = pTmpImg + x){// 从上向下: p3(y,x) = p3(y-1,x) + a * [p2(y,x) - p3(y-1,x)]pre = *pCur;for (y = 1, pCur += width; y < height; y++, pCur += width)*pCur = pre = pre + ALUT[*pCur - pre];pCur -= width; // 回到列尾// 从下向上: p4(i,j) = p4(i+1,j) + a * [p3(i,j) - p4(i+1,j)]for (y = height - 2, pCur -= width; y >= 0; y--, pCur -= width)*pCur = pre = pre + ALUT[*pCur - pre];}// Step 3: 正导数=1，负导数为0，0必须也是0pEnd = pTmpImg + width * height;for (pCur = pTmpImg, pGry = pGryImg; pCur < pEnd; pGry++){*(pCur++) = (*pCur > *pGry);}// Step 4: 过零点检测memset(pSJImg, 0, width * height); // 边缘强度赋0pSJ = pSJImg + width;pCur = pTmpImg + width; // 首行不做 for (y = 1; y < height - 1; y++){pSJ++; pCur++;  // 首列不做for (x = 1; x < width - 1; x++, pGry++, pCur++, pSJ++){if (*pCur) // 正导数{// 下面使用4邻域, 边缘为8连通, 不能保证4连通; 使用8邻域才能保证边缘4连通if ((!*(pCur - 1)) || // 左, 必须<=0, 不能<0(!*(pCur + 1)) || // 右, 必须<=0, 不能<0(!*(pCur - width)) || // 上, 必须<=0, 不能<0(!*(pCur + width)))   // 下, 必须<=0, 不能<0{*pSJ = 255; // 周围有导数小于等于0}}}pSJ++; pCur++;  // 尾列不做}
}

当滤波系数为0.1时，我们来看效果

image-20240422204121096

沈俊算子能够得到闭合的边缘。一种边缘检测方法能够得到一般要用图像分割才能得到的目标轮廓，会具有很高的实用价值。

沈俊算子只需要一个参数a_0，且a_0语义明确，而且沈俊算子的代码量非常小，所以沈俊算子使用起来非常方便。

测试

那么该如何更好的得到我们原始图像的边缘呢

可能一个边缘处理都没有办法得到很好的效果，我们可以结合多个处理办法

image-20240422204713809

我们可以看到效果比较好的是索贝尔算子，那么我们结合一下

//沈俊算子加索贝尔算子
//pGryImg：指向原始灰度图像数据的指针
//pTmpImg：指向辅助图像数据的指针
//width：图像的宽度
//height：图像的高度
//a0：这是沈俊算子的参数，用于控制边缘检测的灵敏度。
//grdThre：这是Sobel算子的梯度阈值
//pEdgeImg：最终边缘图像数据的指针
void RmwExtractRiceEdge(uint8_t* pGryImg,uint8_t* pTmpImg,int width,int height,double a0, int grdThre, uint8_t* pEdgeImg)
{// step.1------------沈俊算子-----------------------//RmwShenJunGryImg(pGryImg, pTmpImg, width, height, a0, pEdgeImg);// step.2------------Sobel算子----------------------//RmwSobelGryImg(pGryImg, width, height, pTmpImg);// step.3------------二者融合-----------------------//for (int i = 0; i < width * height; i++){*(pEdgeImg + i) = (pEdgeImg[i] && (pTmpImg[i] > grdThre)) * 255;}// step.4------------结束---------------------------//return;
}

image-20240422211324490

我们换一张图片试一下

image-20240422212038348

总结

那么到目前为止，我们已经可以成功的进行边缘检测

接下来就是边缘增强，边缘分割等相关内容，会尽快更新

源码

IDP.h

#pragma once#include <stdio.h>
#include <stdint.h>
#include <stdlib.h>
#include <math.h>
#include <nmmintrin.h>uint8_t* readGrayScaleBMP(const char* filename, int* width, int* height);//读取8位灰度图片
void saveGrayScaleBMP(const char* filename, const uint8_t* imageData, int width, int height);// 将8位灰度图像数据保存为BMP文件
uint8_t* readColorBMP(const char* filename, int* width, int* height);//读取24位彩色图像的BMP文件
void saveColorBMP(const char* filename, const uint8_t* imageData, int width, int height);//将24位彩色图像数据保存为BMP文件
void convertToGray(uint8_t* rgbImage, uint8_t* grayImage, int width, int height);//24位彩色图像转8位灰度值
void LinearStretchDemo(uint8_t* pGryImg, int width, int height, double k, double b);//灰度线性拉伸
void GetHistogram(uint8_t* pImg, int width, int height, int* histogram);//统计图像灰度值
void GetBrightContrast(int* histogram, double* bright, double* contrast);//亮度和对比度
void RmwHistogramEqualize(uint8_t* pGryImg, int width, int height);//直方图均衡化
void RmwLogTransform(uint8_t* pGryImg, int width, int height);//对数变换
void RmwAvrFilterBySumCol(uint8_t* pGryImg, int width, int height, int M, int N, uint8_t* pResImg);//基于列积分的快速均值滤波
void RmwDoSumGryImg(uint8_t* pGryImg, int width, int height, int* pSumImg);//基于列积分的积分图实现
void RmwDoSumGryImg_SSE(uint8_t* pGryImg, int width, int height, int* pSumImg);//基于SSE的积分图实现
void RmwAvrFilterBySumImg(int* pSumImg, int width, int height, int M, int N, uint8_t* pResImg);//基于积分图的快速均值滤波  
void GetMedianGry(int* histogram, int N, int* medGry);//求灰度值中值
double RmwMedianFilter(uint8_t* pGryImg, int width, int height, int M, int N, uint8_t* pResImg);//中值滤波
void RmwBinImgFilter(uint8_t* pBinImg, int width, int height, int M, int N, double threshold, uint8_t* pResImg);//二值滤波
void RmwGradientGryImg(uint8_t* pGryImg, int width, int height, uint8_t* pGrdImg);//梯度算子
void RmwGradientGryImgPlus(uint8_t* pGryImg, int width, int height, uint8_t* pGrdImg, int threshold);//梯度算子加阈值
void invertImage(uint8_t* image, int width, int height);//反相
void RmwRobertsGryImg(uint8_t* pGryImg, int width, int height, uint8_t* pRbtImg);//罗伯特算子
void RmwSobelGryImg(uint8_t* pGryImg, int width, int height, uint8_t* pSbImg);//索贝尔算子
void RmwPrewittGryImg(uint8_t* pGryImg, int width, int height, uint8_t* pPRTImg); //Prewitt算子
void RmwShenJunGryImg(uint8_t* pGryImg, uint8_t* pTmpImg, int width, int height, double a0, uint8_t* pSJImg);//沈俊算子
void RmwExtractRiceEdge(uint8_t* pGryImg, uint8_t* pTmpImg, int width, int height, double a0, int grdThre, uint8_t* pEdgeImg);//索贝尔＋沈俊算子

IDP.c

#define _CRT_SECURE_NO_WARNINGS 1#include "IDP.h"//读取8位灰度图片
//filename：字符数组的指针，用于指定要保存的图像文件的名称或路径。
//imageData：无符号 8 位整型数据的指针，代表要保存的图像数据。
//width：图像的宽度。
//height：图像的高度。
uint8_t* readGrayScaleBMP(const char* filename, int* width, int* height) 
{FILE* file = fopen(filename, "rb");if (!file) {fprintf(stderr, "Error opening file %s\n", filename);return NULL;}// 读取BMP文件头部信息uint8_t bmpHeader[54];fread(bmpHeader, 1, 54, file);// 从文件头部提取图像宽度和高度信息*width = *(int*)&bmpHeader[18];*height = *(int*)&bmpHeader[22];// 分配存储图像数据的内存uint8_t* imageData = (uint8_t*)malloc(*width * *height);if (!imageData) {fprintf(stderr, "内存分配失败\n");fclose(file);return NULL;}// 计算调色板的大小int paletteSize = *(int*)&bmpHeader[46];if (paletteSize == 0)paletteSize = 256;// 读取调色板数据uint8_t palette[1024];fread(palette, 1, paletteSize * 4, file);// 读取图像数据fseek(file, *(int*)&bmpHeader[10], SEEK_SET);fread(imageData, 1, *width * *height, file);fclose(file);return imageData;
}// 将8位灰度图像数据保存为BMP文件
//filename：字符数组的指针，用于指定要保存的图像文件的名称或路径。
//imageData：无符号 8 位整型数据的指针，代表要保存的图像数据。
//width：图像的宽度。
//height：图像的高度。
void saveGrayScaleBMP(const char* filename, const uint8_t* imageData, int width, int height) 
{FILE* file = fopen(filename, "wb");if (!file) {fprintf(stderr, "Error creating file %s\n", filename);return;}// BMP文件头部信息uint8_t bmpHeader[54] = {0x42, 0x4D,             // 文件类型标识 "BM"0x36, 0x00, 0x0C, 0x00, // 文件大小（以字节为单位，此处假设图像数据大小不超过4GB）0x00, 0x00,             // 保留字段0x00, 0x00,             // 保留字段0x36, 0x00, 0x00, 0x00, // 位图数据偏移（以字节为单位）0x28, 0x00, 0x00, 0x00, // 位图信息头大小（40字节）0x00, 0x00, 0x00, 0x00, // 图像宽度0x00, 0x00, 0x00, 0x00, // 图像高度0x01, 0x00,             // 目标设备的级别（此处为1，不压缩）0x08, 0x00,             // 每个像素的位数（8位）0x00, 0x00, 0x00, 0x00, // 压缩类型（此处为不压缩）0x00, 0x00, 0x00, 0x00, // 图像数据大小（以字节为单位，此处为0，表示不压缩）0x00, 0x00, 0x00, 0x00, // 水平分辨率（像素/米，此处为0，表示未知）0x00, 0x00, 0x00, 0x00, // 垂直分辨率（像素/米，此处为0，表示未知）0x00, 0x00, 0x00, 0x00, // 使用的颜色索引数（0表示使用所有调色板项）0x00, 0x00, 0x00, 0x00  // 重要的颜色索引数（0表示所有颜色都重要）};// 更新BMP文件头部信息中的宽度和高度*(int*)&bmpHeader[18] = width;*(int*)&bmpHeader[22] = height;// 写入BMP文件头部信息fwrite(bmpHeader, 1, 54, file);// 写入调色板数据for (int i = 0; i < 256; i++) {fputc(i, file);  // 蓝色分量fputc(i, file);  // 绿色分量fputc(i, file);  // 红色分量fputc(0, file);  // 保留字节}// 写入图像数据fwrite(imageData, 1, width * height, file);fclose(file);
}// 读取24位彩色图像的BMP文件
//filename：字符数组的指针，用于指定要读取的 BMP 格式图像文件的名称或路径。
//width：整型变量的指针，用于存储读取的图像的宽度。
//height：整型变量的指针，用于存储读取的图像的高度。
uint8_t* readColorBMP(const char* filename, int* width, int* height) 
{FILE* file = fopen(filename, "rb");if (!file) {fprintf(stderr, "Error opening file %s\n", filename);return NULL;}// 读取BMP文件头部信息uint8_t bmpHeader[54];fread(bmpHeader, 1, 54, file);// 从文件头部提取图像宽度和高度信息*width = *(int*)&bmpHeader[18];*height = *(int*)&bmpHeader[22];// 分配存储图像数据的内存uint8_t* imageData = (uint8_t*)malloc(*width * *height * 3);if (!imageData) {fprintf(stderr, "Memory allocation failed\n");fclose(file);return NULL;}// 读取图像数据fseek(file, *(int*)&bmpHeader[10], SEEK_SET);fread(imageData, 1, *width * *height * 3, file);fclose(file);return imageData;
}//将24位彩色图像数据保存为BMP文件
//filename：字符数组的指针，用于指定要保存的图像文件的名称或路径。
//imageData：无符号 8 位整型数据的指针，代表要保存的图像数据。
//width：图像的宽度。
//height：图像的高度。
void saveColorBMP(const char* filename, const uint8_t* imageData, int width, int height) 
{FILE* file = fopen(filename, "wb");if (!file) {fprintf(stderr, "Error creating file %s\n", filename);return;}// BMP文件头部信息uint8_t bmpHeader[54] = {0x42, 0x4D,             // 文件类型标识 "BM"0x00, 0x00, 0x00, 0x00, // 文件大小（占位，稍后计算）0x00, 0x00,             // 保留字段0x00, 0x00,             // 保留字段0x36, 0x00, 0x00, 0x00, // 位图数据偏移（以字节为单位）0x28, 0x00, 0x00, 0x00, // 位图信息头大小（40字节）0x00, 0x00, 0x00, 0x00, // 图像宽度0x00, 0x00, 0x00, 0x00, // 图像高度0x01, 0x00,             // 目标设备的级别（此处为1，不压缩）0x18, 0x00,             // 每个像素的位数（24位）0x00, 0x00, 0x00, 0x00, // 压缩类型（此处为不压缩）0x00, 0x00, 0x00, 0x00, // 图像数据大小（占位，稍后计算）0x00, 0x00, 0x00, 0x00, // 水平分辨率（像素/米，此处为0，表示未知）0x00, 0x00, 0x00, 0x00, // 垂直分辨率（像素/米，此处为0，表示未知）0x00, 0x00, 0x00, 0x00, // 使用的颜色索引数（0表示使用所有调色板项）0x00, 0x00, 0x00, 0x00  // 重要的颜色索引数（0表示所有颜色都重要）};// 更新BMP文件头部信息中的宽度和高度*(int*)&bmpHeader[18] = width;*(int*)&bmpHeader[22] = height;// 计算图像数据大小uint32_t imageDataSize = width * height * 3 + 54; // 加上文件头部大小bmpHeader[2] = (uint8_t)(imageDataSize & 0xFF);bmpHeader[3] = (uint8_t)((imageDataSize >> 8) & 0xFF);bmpHeader[4] = (uint8_t)((imageDataSize >> 16) & 0xFF);bmpHeader[5] = (uint8_t)((imageDataSize >> 24) & 0xFF);// 写入BMP文件头部信息fwrite(bmpHeader, 1, 54, file);// 写入图像数据fwrite(imageData, width * height * 3, 1, file);fclose(file);
}//24位彩色图像转8位灰度值
//rgbImage原始图像
//grayImage输出灰度图像
//width,height图片的宽和高
void convertToGray(uint8_t* rgbImage, uint8_t* grayImage, int width, int height)
{for (int y = 0; y < height; y++) {for (int x = 0; x < width; x++) {// 获取当前像素的 RGB 分量uint8_t r = rgbImage[3 * (y * width + x) + 0];uint8_t g = rgbImage[3 * (y * width + x) + 1];uint8_t b = rgbImage[3 * (y * width + x) + 2];// 计算灰度值（常用的加权平均法）// 这里使用的加权系数是常见的：R: 0.299, G: 0.587, B: 0.114uint8_t gray = (uint8_t)(0.299 * r + 0.587 * g + 0.114 * b);// 将灰度值写入灰度图像数组grayImage[y * width + x] = gray;}}
}//灰度线性拉伸
//pGryImg：灰度图像数据的指针。
//width：图像的宽度。
//height：图像的高度。
//k：线性拉伸的斜率。它控制着拉伸的速率或程度。当(k) 大于 1 时，图像的对比度增加；当(k) 小于 1 时，对比度降低。
//b：线性拉伸的偏移。它控制着拉伸后灰度值的起始位置。当(b) 大于 0 时，图像的整体亮度增加；当(b) 小于 0 时，整体亮度减小。
void LinearStretchDemo(uint8_t* pGryImg, int width, int height, double k, double b)
{uint8_t* pCur, * pEnd;int LUT[256];    //因为只有[0,255]共256个灰度值//step1. 生成查找表for (int g = 0; g < 256; g++){LUT[g] = max(0, min(255, k * g + b));}//step2. 进行变换for (pCur = pGryImg, pEnd = pGryImg + width * height; pCur < pEnd; pCur++){*pCur = LUT[*pCur];}//step3. 结束return;
}//统计图像灰度值
//pImg：灰度图像数据的指针。
//width：图像的宽度。
//height：图像的高度。
//* histogram：数组首元素地址，需要一个能储存256个变量的整型数组
void GetHistogram(uint8_t* pImg, int width, int height, int* histogram)
{uint8_t* pCur;uint8_t* pEnd = pImg + width * height;// 初始化直方图数组memset(histogram, 0, sizeof(int) * 256);// 直方图统计for (pCur = pImg; pCur < pEnd;){histogram[*pCur]++;pCur++;}// 函数结束return;
}//亮度和对比度
//储存histogram灰度直方图的指针
//接收亮度的变量地址
//接收对比度的变量地址
void GetBrightContrast(int* histogram, double* bright, double* contrast)
{int g;double sum, num; //书上说图像很亮时，int有可能会溢出，所以我这里直接用doubledouble fsum;//step.1 求亮度for (sum = num = 0, g = 0; g < 256; g++){sum += histogram[g] * g;num += histogram[g];}*bright = sum * 1.0 / num;//step.2 求对比度for (fsum = 0.0, g = 0; g < 256; g++){fsum += histogram[g] * (g - *bright) * (g - *bright);}*contrast = sqrt(fsum / (num - 1)); //即Std Dev//step.3 结束return;
}//pGryImg：灰度图像数据的指针。
//width：图像的宽度。
//height：图像的高度。
void RmwHistogramEqualize(uint8_t* pGryImg, int width, int height)
{uint8_t* pCur, * pEnd = pGryImg + width * height; // 指针变量，指向当前像素和图像末尾int histogram[256], LUT[256], A, g; // 直方图数组、查找表数组、累积直方图、灰度级// step.1-------------求直方图--------------------------//memset(histogram, 0, sizeof(int) * 256); // 初始化直方图数组为0for (pCur = pGryImg; pCur < pEnd;)histogram[*(pCur++)]++; // 统计每个灰度级出现的频率// step.2-------------求LUT[g]-------------------------//A = histogram[0]; // 初始化累积直方图的值为第一个灰度级的频率LUT[0] = 255 * A / (width * height); // 计算第一个灰度级对应的均衡化后的灰度值for (g = 1; g < 256; g++) {A += histogram[g]; // 更新累积直方图的值LUT[g] = 255 * A / (width * height); // 计算当前灰度级对应的均衡化后的灰度值}// step.3-------------查表------------------------------//for (pCur = pGryImg; pCur < pEnd;)*(pCur++) = LUT[*pCur]; // 使用查找表对每个像素进行灰度映射// step.4-------------结束------------------------------//return;
}//对数变换
//pGryImg：灰度图像数据的指针。
//width：图像的宽度。
//height：图像的高度。
void RmwLogTransform(uint8_t* pGryImg, int width, int height)
{uint8_t* pCur, * pEnd = pGryImg + width * height; // 指向灰度图像数据的当前指针和结束指针int histogram[256], LUT[256], gmax, g; // 声明直方图数组、查找表数组、最大灰度值、当前灰度值double c; // 声明常数c// step.1-------------求直方图--------------------------//memset(histogram, 0, sizeof(int) * 256); // 初始化直方图数组为0for (pCur = pGryImg; pCur < pEnd;)histogram[*(pCur++)]++; // 遍历图像数据，统计每个灰度级的像素数量// step.2-------------最大值---------------------------//for (gmax = 255; gmax >= 0; gmax++)if (histogram[gmax]) break; // 从最大灰度级开始向低灰度级搜索，找到第一个非零灰度级，即最大灰度值// step.3-------------求LUT[g]-------------------------//c = 255.0 / log(1 + gmax); // 计算常数cfor (g = 0; g < 256; g++){LUT[g] = (int)(c * log(1 + g)); // 根据对数变换公式计算查找表中每个灰度级的映射值}// step.4-------------查表------------------------------//for (pCur = pGryImg; pCur < pEnd;)*(pCur++) = LUT[*pCur]; // 使用查找表将图像数据进行对数变换// step.5-------------结束------------------------------//return; // 函数结束
}//基于列积分的快速均值滤波
//原始灰度图像
//图像的宽度和高度
//滤波邻域：M列N行
//结果图像
void RmwAvrFilterBySumCol(uint8_t* pGryImg,int width, int height,int M, int N,uint8_t* pResImg) 
{uint8_t* pAdd, * pDel, * pRes;int halfx, halfy;int x, y;int sum, c;int sumCol[4096]; // 约定图像宽度不大于4096// step.1------------初始化--------------------------//M = M / 2 * 2 + 1; // 奇数化N = N / 2 * 2 + 1; // 奇数化halfx = M / 2; // 滤波器的半径xhalfy = N / 2; // 滤波器的半径yc = (1 << 23) / (M * N); // 乘法因子memset(sumCol, 0, sizeof(int) * width);for (y = 0, pAdd = pGryImg; y < N; y++) {for (x = 0; x < width; x++) sumCol[x] += *(pAdd++);}// step.2------------滤波----------------------------//for (y = halfy, pRes = pResImg + y * width, pDel = pGryImg; y < height - halfy; y++) {// 初值for (sum = 0, x = 0; x < M; x++) sum += sumCol[x];// 滤波pRes += halfx; // 跳过左侧for (x = halfx; x < width - halfx; x++) {// 求灰度均值// *(pRes++)=sum/(N*M);*(pRes++) = (sum * c) >> 23; // 用整数乘法和移位代替除法// 换列,更新灰度和sum -= sumCol[x - halfx]; // 减左边列sum += sumCol[x + halfx + 1]; // 加右边列}pRes += halfx; // 跳过右侧// 换行,更新sumColfor (x = 0; x < width; x++) {sumCol[x] -= *(pDel++); // 减上一行sumCol[x] += *(pAdd++); // 加下一行}}// step.3------------返回----------------------------//return;
}//基于列积分的积分图实现
//pGryImg, // 原始灰度图像
//width,       // 图像的宽度 
//height,      // 图像的高度
//pSumImg     // 计算得到的积分图
void RmwDoSumGryImg(uint8_t* pGryImg,int width,int height, int* pSumImg)
{uint8_t* pGry;int* pRes;int x, y;int sumCol[4096]; // 约定图像宽度不大于4096memset(sumCol, 0, sizeof(int) * width);for (y = 0, pGry = pGryImg, pRes = pSumImg; y < height; y++){// 最左侧像素的特别处理sumCol[0] += *(pGry++);*(pRes++) = sumCol[0];// 正常处理for (x = 1; x < width; x++){sumCol[x] += *(pGry++);       // 更新列积分int temp = *(pRes - 1);*(pRes++) = temp + sumCol[x];}}return;
}//基于SSE的积分图实现
//pGryImg原始灰度图像
//width图像的宽度，必须是4的倍数
//height图像的高度
//pSumImg计算得到的积分图
void RmwDoSumGryImg_SSE(uint8_t* pGryImg,int width,int height,int* pSumImg)
{int sumCol[4096]; //约定图像宽度不大于4096__m128i* pSumSSE, A;uint8_t* pGry;int* pRes;int x, y;memset(sumCol, 0, sizeof(int) * width);for (y = 0, pGry = pGryImg, pRes = pSumImg; y < height; y++){// 0:需要特别处理sumCol[0] += *(pGry++);*(pRes++) = sumCol[0];// 1sumCol[1] += *(pGry++);*(pRes++) = *(pRes - 1) + sumCol[1];// 2sumCol[2] += *(pGry++);*(pRes++) = *(pRes - 1) + sumCol[2];// 3sumCol[3] += *(pGry++);*(pRes++) = *(pRes - 1) + sumCol[3];// [4...width-1]for (x = 4, pSumSSE = (__m128i*)(sumCol + 4); x < width; x += 4, pGry += 4){// 把变量的低32位(有4个8位整数组成)转换成32位的整数A = _mm_cvtepu8_epi32(_mm_loadl_epi64((__m128i*)pGry));// 4个32位的整数相加*(pSumSSE++) = _mm_add_epi32(*pSumSSE, A);// 递推*(pRes++) = *(pRes - 1) + sumCol[x + 0];*(pRes++) = *(pRes - 1) + sumCol[x + 1];*(pRes++) = *(pRes - 1) + sumCol[x + 2];*(pRes++) = *(pRes - 1) + sumCol[x + 3];}}return;
}//基于积分图的快速均值滤波
//pSumImg计算得到的积分图
//width,height,图像的宽度和高度
//M, N,滤波邻域：M列N行
//pResImg 结果图像
void RmwAvrFilterBySumImg(int* pSumImg,int width, int height,int M, int N,uint8_t* pResImg)
{// 没有对边界上邻域不完整的像素进行处理，可以采用变窗口的策略int* pY1, * pY2;uint8_t* pRes;int halfx, halfy;int y, x1, x2;int sum, c;// step.1------------初始化--------------------------//M = M / 2 * 2 + 1; // 奇数化N = N / 2 * 2 + 1; // 奇数化halfx = M / 2;      // 滤波器的半径xhalfy = N / 2;      // 滤波器的半径yc = (1 << 23) / (M * N); // 乘法因子// step.2------------滤波----------------------------//for (y = halfy + 1, pRes = pResImg + y * width, pY1 = pSumImg, pY2 = pSumImg + N * width;y < height - halfy;y++, pY1 += width, pY2 += width){pRes += halfx + 1; // 跳过左侧for (x1 = 0, x2 = M; x2 < width; x1++, x2++) // 可以简化如此，但不太容易读{sum = *(pY2 + x2) - *(pY2 + x1) - *(pY1 + x2) + *(pY1 + x1);*(pRes++) = (uint8_t)((sum * c) >> 23); // 用整数乘法和移位代替除法}pRes += halfx; // 跳过右侧}// step.3------------返回----------------------------//return;
}void GetMedianGry(int* histogram, int N, int* medGry)
{int g;int num;// step.1-------------求灰度中值------------------------//num = 0;for (g = 0; g < 256; g++){num += histogram[g];if (2 * num > N) break;  //num>N/2}*medGry = g;// step.2-------------结束------------------------------//return;
}//中值滤波
//pGryImg：指向待处理灰度图像数据的指针。
//width、height：表示图像的宽度和高度。
//M、N：分别表示中值滤波器的水平和垂直邻域大小（以像素为单位）。
//pResImg：指向存储结果图像数据的指针。
double RmwMedianFilter(uint8_t* pGryImg, int width, int height, int M, int N, uint8_t* pResImg) 
{uint8_t* pCur, * pRes;int halfx, halfy, x, y, i, j, y1, y2;int histogram[256];int wSize, j1, j2;int num, med, v;int dbgCmpTimes = 0; // 搜索中值所需比较次数的调试M = M / 2 * 2 + 1; // 奇数化N = N / 2 * 2 + 1; // 奇数化halfx = M / 2;      // x半径halfy = N / 2;      // y半径wSize = (halfx * 2 + 1) * (halfy * 2 + 1); // 邻域内像素总个数for (y = halfy, pRes = pResImg + y * width; y < height - halfy; y++) {// step.1----初始化直方图y1 = y - halfy;y2 = y + halfy;memset(histogram, 0, sizeof(int) * 256);for (i = y1, pCur = pGryImg + i * width; i <= y2; i++, pCur += width) {for (j = 0; j < halfx * 2 + 1; j++) {histogram[*(pCur + j)]++;}}// step.2-----初始化中值num = 0; // 记录着灰度值从0到中值的个数for (i = 0; i < 256; i++) {num += histogram[i];if (num * 2 > wSize) {med = i;break;}}// 滤波pRes += halfx; // 没有处理图像左边界侧的像素for (x = halfx; x < width - halfx; x++) {// 赋值*(pRes++) = med;// step.3-----直方图递推: 减去当前邻域最左边的一列,添加邻域右侧的一个新列j1 = x - halfx;     // 最左边列j2 = x + halfx + 1; // 右边的新列for (i = y1, pCur = pGryImg + i * width; i <= y2; i++, pCur += width) {// 减去最左边列v = *(pCur + j1);histogram[v]--;  // 更新直方图if (v <= med) num--; // 更新num// 添加右边的新列v = *(pCur + j2);histogram[v]++; // 更新直方图if (v <= med) num++; // 更新num}// step.4-----更新中值if (num * 2 < wSize) { // 到上次中值med的个数不够了,则med要变大for (med = med + 1; med < 256; med++) {dbgCmpTimes += 2; // 总的比较次数,调试用num += histogram[med];if (num * 2 > wSize) break;}dbgCmpTimes += 1; // 总的比较次数,调试用}else { // 到上次中值med的个数多了,则med要变小while ((num - histogram[med]) * 2 > wSize) { // 若减去后,仍变小dbgCmpTimes++; // 总的比较次数,调试用num -= histogram[med];med--;}dbgCmpTimes += 2; // 总的比较次数,调试用}}pRes += halfx; // 没有处理图像右边界侧的像素}// 返回搜索中值需要的平均比较次数return dbgCmpTimes * 1.0 / ((width - halfx * 2) * (height - halfy * 2));
}//二值滤波
//pBinImg,  原始二值图像
// width, height,图像的宽度和高度
// M, N, 滤波邻域：M列N行
// threshold, 灰度阈值,大于等于该值时结果赋255
// pResImg 结果图像
void RmwBinImgFilter(uint8_t* pBinImg,int width, int height,int M, int N,double threshold,uint8_t* pResImg )
{// 没有对边界上邻域不完整的像素进行处理，可以采用变窗口的策略uint8_t* pAdd, * pDel, * pRes;int halfx, halfy;int x, y, sum, sumThreshold;int sumCol[4096]; //约定图像宽度不大于4096// step.1------------初始化--------------------------//M = M / 2 * 2 + 1; //奇数化N = N / 2 * 2 + 1; //奇数化halfx = M / 2; //滤波器的x半径halfy = N / 2; //滤波器的y半径sumThreshold = max(1, (int)(threshold * M * N)); //转换成邻域内灰度值之和的阈值memset(sumCol, 0, sizeof(int) * width);for (y = 0, pAdd = pBinImg; y < N; y++){for (x = 0; x < width; x++)sumCol[x] += *(pAdd++);}// step.2------------滤波----------------------------//for (y = halfy, pRes = pResImg + y * width, pDel = pBinImg; y < height - halfy; y++){//初值for (sum = 0, x = 0; x < M; x++)sum += sumCol[x];//滤波pRes += halfx; //跳过左侧for (x = halfx; x < width - halfx; x++){//求灰度均值/*if (sum>=sumThreshold){*(pRes++) = 255;}else  *(pRes++) = 0;*/*(pRes++) = (sum >= sumThreshold) * 255; //请理解这个表达式的含义//换列,更新灰度和sum -= sumCol[x - halfx];     //减左边列sum += sumCol[x + halfx + 1]; //加右边列}pRes += halfx; //跳过右侧//换行,更新sumColfor (x = 0; x < width; x++){sumCol[x] -= *(pDel++); //减上一行sumCol[x] += *(pAdd++); //加下一行}}// step.3------------返回----------------------------//return;
}//梯度算子加阈值
//pGryImg：输入的灰度图像数据指针。
//width：图像的宽度。
//height：图像的高度。
//pGrdImg：输出的梯度图像数据指针
void RmwGradientGryImg(uint8_t* pGryImg, int width, int height, uint8_t* pGrdImg)
{uint8_t* pGry, * pGrd;int dx, dy;int x, y;for (y = 0, pGry = pGryImg, pGrd = pGrdImg; y < height - 1; y++){for (x = 0; x < width - 1; x++, pGry++){dx = *pGry - *(pGry + 1);dy = *pGry - *(pGry + width);int gradient = (int)(sqrt(dx * dx * 1.0 + dy * dy));*pGrd++ = (gradient > 255) ? 255 : gradient;}*pGrd++ = 0; //尾列不做,边缘强度赋0pGry++;}memset(pGrd, 0, width); //尾行不做,边缘强度赋0
}//梯度算子加阈值
//pGryImg：输入的灰度图像数据指针。
//width：图像的宽度。
//height：图像的高度。
//pGrdImg：输出的梯度图像数据指针
void RmwGradientGryImgPlus(uint8_t* pGryImg, int width, int height, uint8_t* pGrdImg, int threshold)
{uint8_t* pGry, * pGrd;int dx, dy;int x, y;for (y = 0, pGry = pGryImg, pGrd = pGrdImg; y < height - 1; y++){for (x = 0; x < width - 1; x++, pGry++){dx = *pGry - *(pGry + 1);dy = *pGry - *(pGry + width);int gradient = (int)(sqrt(dx * dx * 1.0 + dy * dy));*(pGrd++) = (gradient > threshold) ? min(255, gradient) : 0;}*(pGrd++) = 0; //尾列不做,边缘强度赋0pGry++;}memset(pGrd, 0, width); //尾行不做,边缘强度赋0
}//反相
void invertImage(uint8_t* image, int width, int height) {for (int i = 0; i < width * height; i++) {image[i] = 255 - image[i];}
}//罗伯特算子
void RmwRobertsGryImg(uint8_t* pGryImg, int width, int height, uint8_t* pRbtImg)
{uint8_t* pGry, * pRbt;int dx, dy;int x, y;for (y = 0, pGry = pGryImg, pRbt = pRbtImg; y < height - 1; y++){for (x = 0; x < width - 1; x++, pGry++){dx = *pGry - *(pGry + width + 1);dy = *(pGry + 1) - *(pGry + width);*pRbt++ = (uint8_t)(dx > dy ? dx : dy); // 使用三目运算符选择较大的值}*pRbt++ = 0; // 尾列不做, 边缘强度赋0pGry++;}memset(pRbt, 0, width); // 尾行不做, 边缘强度赋0
}//索贝尔算子
void RmwSobelGryImg(uint8_t* pGryImg, int width, int height, uint8_t* pSbImg)
{uint8_t* pGry, * pSb;int dx, dy;int x, y;memset(pSbImg, 0, width); // 首行不做, 边缘强度赋0for (y = 1, pGry = pGryImg + width, pSb = pSbImg + width; y < height - 1; y++){*pSb++ = 0; // 首列不做, 边缘强度赋0pGry++;for (x = 1; x < width - 1; x++, pGry++){// 求dxdx = *(pGry - 1 - width) + (*(pGry - 1) * 2) + *(pGry - 1 + width);dx -= *(pGry + 1 - width) + (*(pGry + 1) * 2) + *(pGry + 1 + width);// 求dydy = *(pGry - width - 1) + (*(pGry - width) * 2) + *(pGry - width + 1);dy -= *(pGry + width - 1) + (*(pGry + width) * 2) + *(pGry + width + 1);// 结果*pSb++ = (uint8_t)min(255, abs(dx) + abs(dy));}*pSb++ = 0; // 尾列不做, 边缘强度赋0pGry++;}memset(pSb, 0, width); // 尾行不做, 边缘强度赋0
}//Prewitt算子
void RmwPrewittGryImg(uint8_t* pGryImg, int width, int height, uint8_t* pPRTImg)
{uint8_t* pGry, * pPRT;int dx, dy, d45, d135, v1, v2;int x, y;memset(pPRTImg, 0, width); // 首行不做, 边缘强度赋0for (y = 1, pGry = pGryImg + width, pPRT = pPRTImg + width; y < height - 1; y++){*pPRT++ = 0; // 首列不做, 边缘强度赋0pGry++;for (x = 1; x < width - 1; x++, pGry++){// 求dxdx = *(pGry - 1 - width) + *(pGry - 1) + *(pGry - 1 + width);dx -= *(pGry + 1 - width) + *(pGry + 1) + *(pGry + 1 + width);// 求dydy = *(pGry - width - 1) + *(pGry - width) + *(pGry - width + 1);dy -= *(pGry + width - 1) + *(pGry + width) + *(pGry + width + 1);// 求45度d45 = *(pGry - width - 1) + *(pGry - width) + *(pGry - 1);d45 -= *(pGry + width + 1) + *(pGry + width) + *(pGry + 1);// 求135度d135 = *(pGry - width) + *(pGry - width + 1) + *(pGry + 1);d135 -= *(pGry + width - 1) + *(pGry + width) + *(pGry - 1);// 结果v1 = abs(dx) > abs(dy) ? abs(dx) : abs(dy);v2 = abs(d45) > abs(d135) ? abs(d45) : abs(d135);*pPRT++ = (uint8_t)((v1 > v2) ? ((v1 > 255) ? 255 : v1) : ((v2 > 255) ? 255 : v2));}*pPRT++ = 0; // 尾列不做, 边缘强度赋0pGry++;}memset(pPRT, 0, width); // 尾行不做, 边缘强度赋0
}//沈俊算子
//pGryImg 和 pTmpImg 是指向 uint8_t 类型的指针，它们分别指向原始灰度图像数据和辅助图像数据。
//width 和 height 是整型参数，表示图像的宽度和高度。
//a0 是双精度浮点型参数，表示滤波系数。
//pSJImg 是指向 uint8_t 类型的指针，它指向了输出的图像数据。
void RmwShenJunGryImg(uint8_t* pGryImg,uint8_t* pTmpImg, int width, int height, double a0, uint8_t* pSJImg)
{uint8_t* pGry, * pCur, * pSJ, * pEnd;int LUT[512], * ALUT; // a0查找表int x, y, pre, dif;// Step 1: 初始化查找表a0 = (a0 < 0.01) ? 0.01 : ((a0 > 0.99) ? 0.99 : a0); // 安全性检查// a0查找表, 进行了四舍五入ALUT = LUT + 256;for (ALUT[0] = 0, dif = 1; dif < 256; dif++){ALUT[dif] = (int)(dif * a0 + 0.5);ALUT[-dif] = (int)(-dif * a0 - 0.5);}// Step 2: 递推实现指数滤波// 按行滤波for (y = 0, pGry = pGryImg, pCur = pTmpImg; y < height; y++){// 从左向右: p1(y,x) = p1(y,x-1) + a * [p(y,x) - p1(y,x-1)]*(pCur++) = pre = *(pGry++);for (x = 1; x < width; x++, pGry++)*(pCur++) = pre = pre + ALUT[*pGry - pre];pCur--; // 回到行尾// 从右向左: p2(y,x) = p2(y,x+1) - a * [p1(y,x) - p2(y,x+1)]for (x = width - 2, pCur = pCur - 1; x >= 0; x--)*(pCur--) = pre = pre + ALUT[*pCur - pre];pCur += (width + 1); // 回到下一行的开始}// 按列滤波for (x = 0, pCur = pTmpImg; x < width; x++, pCur = pTmpImg + x){// 从上向下: p3(y,x) = p3(y-1,x) + a * [p2(y,x) - p3(y-1,x)]pre = *pCur;for (y = 1, pCur += width; y < height; y++, pCur += width)*pCur = pre = pre + ALUT[*pCur - pre];pCur -= width; // 回到列尾// 从下向上: p4(i,j) = p4(i+1,j) + a * [p3(i,j) - p4(i+1,j)]for (y = height - 2, pCur -= width; y >= 0; y--, pCur -= width)*pCur = pre = pre + ALUT[*pCur - pre];}// Step 3: 正导数=1，负导数为0，0必须也是0pEnd = pTmpImg + width * height;for (pCur = pTmpImg, pGry = pGryImg; pCur < pEnd; pGry++){*(pCur++) = (*pCur > *pGry);}// Step 4: 过零点检测memset(pSJImg, 0, width * height); // 边缘强度赋0pSJ = pSJImg + width;pCur = pTmpImg + width; // 首行不做 for (y = 1; y < height - 1; y++){pSJ++; pCur++;  // 首列不做for (x = 1; x < width - 1; x++, pGry++, pCur++, pSJ++){if (*pCur) // 正导数{// 下面使用4邻域, 边缘为8连通, 不能保证4连通; 使用8邻域才能保证边缘4连通if ((!*(pCur - 1)) || // 左, 必须<=0, 不能<0(!*(pCur + 1)) || // 右, 必须<=0, 不能<0(!*(pCur - width)) || // 上, 必须<=0, 不能<0(!*(pCur + width)))   // 下, 必须<=0, 不能<0{*pSJ = 255; // 周围有导数小于等于0}}}pSJ++; pCur++;  // 尾列不做}
}//沈俊算子加索贝尔算子
//pGryImg：指向原始灰度图像数据的指针
//pTmpImg：指向辅助图像数据的指针
//width：图像的宽度
//height：图像的高度
//a0：这是沈俊算子的参数，用于控制边缘检测的灵敏度。
//grdThre：这是Sobel算子的梯度阈值
//pEdgeImg：最终边缘图像数据的指针
void RmwExtractRiceEdge(uint8_t* pGryImg,uint8_t* pTmpImg,int width,int height,double a0, int grdThre, uint8_t* pEdgeImg)
{// step.1------------沈俊算子-----------------------//RmwShenJunGryImg(pGryImg, pTmpImg, width, height, a0, pEdgeImg);// step.2------------Sobel算子----------------------//RmwSobelGryImg(pGryImg, width, height, pTmpImg);// step.3------------二者融合-----------------------//for (int i = 0; i < width * height; i++){*(pEdgeImg + i) = (pEdgeImg[i] && (pTmpImg[i] > grdThre)) * 255;}// step.4------------结束---------------------------//return;
}

本文来自互联网用户投稿，该文观点仅代表作者本人，不代表本站立场。本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如若转载，请注明出处：http://www.hqwc.cn/news/641236.html

如若内容造成侵权/违法违规/事实不符，请联系编程知识网进行投诉反馈email:809451989@qq.com，一经查实，立即删除！

相关文章

【YOLOv8改进[注意力]】YOLOv8添加DAT(Vision Transformer with Deformable Attention)助力涨点

【YOLOv8改进[注意力]】YOLOv8添加DAT(Vision Transformer with Deformable Attention)助力涨点

目录 一 DAT 二 YOLOv8添加DAT助力涨点 1 总体修改 2 配置文件 3 训练 其他 一 DAT 官方论文地址&#xff1a;https://openaccess.thecvf.com/content/CVPR2022/papers/Xia_Vision_Transformer_With_Deformable_Attention_CVPR_2022_paper.pdf Transformers最近在各种视…
阅读更多...
学生选课及成绩查询管理系统的设计与开发C#(winform + sqlserver)

学生选课及成绩查询管理系统的设计与开发C#(winform + sqlserver)

源码来自网络 技术栈&#xff1a; C#的窗体程序开发 本系统未采用C#实现MDI——多文档窗口&#xff0c;因为考虑到C#的该技术与java类似&#xff0c;而暑期java实训时&#xff0c;曾用过类似的方法做过停车场管理系统&#xff0c;所以想为这次的系统注入一点新鲜的血液&#x…
阅读更多...
C++修炼之路之多态---多态的原理(虚函数表)

C++修炼之路之多态---多态的原理(虚函数表)

目录 一&#xff1a;多态的原理 1.虚函数表 2.原理分析 3.对于虚表存在哪里的探讨 4.对于是不是所有的虚函数都要存进虚函数表的探讨 二&#xff1a;多继承中的虚函数表 三&#xff1a;常见的问答题 接下来的日子会顺顺利利&#xff0c;万事胜意&#xff0c;生活明朗--…
阅读更多...
五、【易 AI】鼠标事件与目标焦点

五、【易 AI】鼠标事件与目标焦点

所谓辉煌的人生，不过是欲望的囚徒。 ——叔本华 注：自本节开始，所有的示例都以 OpenGLWidget 实现， 一、鼠标事件 重写鼠标事件， #ifndef MYOPENGL_H #define MYOPENGL_H#include <QOpenGLWidget> #include <QTimer> #include <QMouseEvent>class My…
阅读更多...
CUDA编程【2】-（51-78）

CUDA编程【2】-（51-78）

系列文章目录 文章目录 系列文章目录前言51、寄存器溢出51.1 溢出概念51.1 使用控制 52、本地内存和共享内存52.1 本地内存52.2. 共享内存 53. 常量内存53.1 概念53.2 初始化 54. 全局内存54.1 概念54.2 初始化 55. GPU缓存和变量作用域55.1 缓存类型55.2 变量作用域 56. 静态全…
阅读更多...
鲁抗医药专属采购商城上线，携手隆道公司注入数字化采购新动能

鲁抗医药专属采购商城上线，携手隆道公司注入数字化采购新动能

近日&#xff0c;国内领先的医药制造企业——山东鲁抗医药股份有限公司&#xff08;以下简称鲁抗医药&#xff09;与隆道公司联手打造的鲁抗医药专属采购商城上线运行。该商城&#xff0c;通过整合鲁抗医药合作电商和合格供应商资源&#xff0c;创新商城化采供协同模式&#xf…
阅读更多...
【Hello算法】 ＞ 第 3 关 ＞栈与队列

【Hello算法】 ＞ 第 3 关 ＞栈与队列

数据结构 之 数组与链表 1 栈 / 栈的常见操作、实现、应用2 队列 /队列的常见操作、实现、应用3 双向队列4 Tips ———————————————————————————————————————————————————————————- ————————————————…
阅读更多...
动手学大模型应用开发--Chapter 03搭建并使用向量数据库

动手学大模型应用开发--Chapter 03搭建并使用向量数据库

文章目录 前言一、学习目标二、学习知识点概要2.1 什么是词向量2.2 词向量的意义2.3 文本转为词向量的方法 三、总结四、引申阅读 前言 本学习笔记为datawhale动手学大模型应用开发的第三章&#xff0c;学习链接为&#xff1a; https://datawhalechina.github.io/llm-universe…
阅读更多...
系统启动修复和SYSTEM丢失损坏故障处理

系统启动修复和SYSTEM丢失损坏故障处理

系统启动修复和SYSTEM丢失损坏故障处理 一、问题描述 你的电脑/设备需要修复。无法加载应用程序或操作系统&#xff0c;原因是所需文件丢失或包含错误。 文件:\Windows\system32\winload.exe 错误代码: 0xc000000e 二、问题分析 1.查询winload.exe是win7或者win10以上系统…
阅读更多...
算法竞赛相关问题总结记录

算法竞赛相关问题总结记录

前言 日常在校生或者是工作之余的同学或多或少都会参加一些竞赛,参加竞赛一方面可以锻炼自己的理解与实践能力&#xff0c;也能够增加自己的生活费&#xff0c;竞赛中的一些方案也可以后续作为自己论文的base,甚至是横向课题的框架。在算法竞赛中算法的差别个人感觉差距都不大&…
阅读更多...
Transformer - 时间特征的处理

Transformer - 时间特征的处理

Transformer - 时间特征的处理 flyfish ETTm1.csv有如下内容 假如有2016/7/1 0:45:00有这样的时间字符串&#xff0c;如何变成时间特征列表 from typing import Listimport numpy as np import pandas as pd from pandas.tseries import offsets from pandas.tseries.freq…
阅读更多...
携程 Java 暑期实习一面：HashMap 的 key 可以设置为 null 吗？那 ConcurrentHashMap 呢？

携程 Java 暑期实习一面：HashMap 的 key 可以设置为 null 吗？那 ConcurrentHashMap 呢？

更多大厂面试内容可见 -> http://11come.cn 携程 Java 暑期实习一面&#xff1a;HashMap 的 key 可以设置为 null 吗&#xff1f;那 ConcurrentHashMap 呢&#xff1f; Java 基础 1、Java 中有哪些常见的数据结构&#xff1f; 图片来源于&#xff1a;JavaGuide Java 中常…
阅读更多...
推荐文章
最新文章

Copyright @ 2022~2023