简介面向数字图像处理课程期末考核的一份完整PDF报告聚焦图像平滑、边缘检测、二值化与Hough变换检测圆形四大核心环节可作为高年级本科生与研究生完成同类大作业的设计参考。资源包为1个PDF文件大小1.09MB已有850人学习浏览。报告以硬币图像为输入完整覆盖实验流程使用方差为1的高斯滤波器对图像进行平滑与拉普拉斯预处理采用Prewitt梯度算子提取水平、垂直、45度和-45度四个方向的边缘通过循环确定全局阈值完成二值化再借助Hough变换在参数空间中检测圆的圆心与半径并将检测到的圆形叠加回原始灰度图像实现边缘增强。正文还给出了算法原理公式、实验流程图、中间结果和最终结果截图以及对应的实现代码方便读者对照复现。整体结构清晰既适合期末大作业的完整参考也适合图像处理课程中Hough变换、边缘检测等知识点的专项复习。1. 从一副硬币图到 Hough 圆检测的完整链路很多人在数字图像处理课程里都能把 Fourier 变换、Prewitt 算子、Hough 变换的公式背得很熟但真拿到一副硬币图像、要求检测其中的圆并叠加回原图时往往会卡在“先做什么、后做什么、参数怎么设”上。这篇博文不是把课程报告复述一遍而是以一个实际的数字图像处理期末大作业为样本把「高斯平滑 → 拉普拉斯变换 → 四方向 Prewitt 边缘检测 → 循环阈值二值化 → Hough 变换检测圆 → 与原图叠加增强」这条完整链路拆开讲清楚。适合正在做图像处理课程设计的学生也适合需要快速上手 MATLAB 实现圆形检测的工程师。文中所有代码都基于 MATLAB参数会逐一解释最后还会给出我在调试时踩过的坑。2. LoG 预处理为什么平滑和拉普拉斯要合并成一次卷积2.1 先平滑再求导的数学逻辑拿到灰度图像后直接做边缘检测噪声会把边缘响应搅得很难看。标准做法是先平滑、再求二阶导也就是 LoGLaplacian of Gaussian算子。这里的关键是把两步合成一步先用高斯函数对图像做卷积再对结果做拉普拉斯变换即[ \nabla^2[h(r)*f(x,y)] [\nabla^2 h(r)] * f(x,y) ]其中 (h(r)\frac{1}{2\pi\sigma^2}e^{-r^2/(2\sigma^2)})而 (\nabla^2 h(r)\frac{r^2-2\sigma^2}{\sigma^4}e^{-r^2/(2\sigma^2)})。也就是说先平滑再拉普拉斯等价于直接用一个 LoG 核与图像卷积省一次卷积运算。作业里取 (\sigma^21)核的尺寸和整个图像一样大然后用傅里叶变换做频域乘法避免时域大核卷积的耗时。2.2 频域实现时最容易错的坐标原点问题直接看原文的 MATLAB 代码h 矩阵的生成方式是t i^2 j^2这是以图像左上角为原点的距离平方。这样做出来的 h 不是真正中心对称的 LoG 核因为 LoG 核应该以核中心为原点。常见做法是先把坐标中心移到图像中点再计算距离S imread(Sweden_coins.bmp); [m, n] size(S); sigma2 1; % 方差 sigma^2 1 h zeros(m, n); for i 1:m for j 1:n x i - (m 1) / 2; % 将坐标原点移到核中心 y j - (n 1) / 2; r2 x^2 y^2; h(i, j) (r2 - 2 * sigma2) * exp(-r2 / (2 * sigma2)) / sigma2^2; end end % 频域滤波逐点相乘再反变换 H fftshift(fft2(h)); G fftshift(fft2(S)); W H .* G; w ifft2(ifftshift(W)); figure; imshow(w, []);这里fftshift把零频移到中心目的是让频域相乘时低频对齐。做完ifft2后要再用ifftshift调回来否则图像会被平移错位。imshow 默认会把负值截断到 0所以显示时最好加[]让图像自动拉伸灰度范围不然洛普拉斯结果的负响应会显示成全黑。我一般还会直接用imagesc(w)配合colormap gray看细节比imshow更直观。2.3 核尺寸与边界处理的影响作业里核和原图等大循环计算量是 (m \times n) 次指数运算对一副几百像素的图像还能接受但再大就会很慢。实际工程中 LoG 核通常取 5×5 或 7×7 就够了因为高斯尾巴衰减快。核太小会让平滑不足核太大则边缘定位偏差明显。边界像素在卷积时会有响应异常MATLAB 的conv2默认补零而频域乘法相当于周期延拓两种方式在图像边缘会有几像素的差异。如果后续 Hough 检测对边缘位置敏感建议先裁掉边界 5 像素再做后续处理。另一个细节是这里的 h 没有归一化LoG 核的累加和接近 0对直流分量几乎没有响应所以滤波结果整体均值也会接近 0这是正常的。3. 四方向 Prewitt 边缘检测与合成3.1 为什么选 Prewitt 而不是 SobelPrewitt 和 Sobel 的区别在于平滑权重Sobel 对中心行/列给了 2 倍权重更强调中心像素Prewitt 所有权重都是 1结构更简单各向异性响应也更均匀。在 Hough 变换检测圆这个场景里我们对边缘的连续性和方向覆盖比单点精度更在意Prewitt 的四方向掩膜可以直接覆盖 0°、45°、90°、135° 四个方向。Sobel 虽然对噪声更鲁棒但它的平滑核只对水平和垂直方向做了对称加权斜方向掩膜在离散网格上的定义不如 Prewitt 直观。3.2 四个方向的掩膜与卷积实现Prewitt 四个方向的 3×3 掩膜如下方向掩膜模板水平检测垂直边缘[-1 -1 -1; 0 0 0; 1 1 1]垂直检测水平边缘[-1 0 1; -1 0 1; -1 0 1]45°[-1 -1 0; -1 0 1; 0 1 1]-45°[0 -1 -1; 1 0 -1; 1 1 0]注意方向命名在不少教材里是反的有的把 H1 叫水平边缘指的其实是水平方向上的灰度变化也就是检测垂直边缘。作业中四个掩膜与图像逐点卷积后直接取绝对值相加。代码写法如下IMG double(w); % 上一步的 LoG 结果转 double 避免溢出 H1 zeros(m, n); % 水平方向梯度垂直边缘 H2 zeros(m, n); % 垂直方向梯度水平边缘 G1 zeros(m, n); % 45 度 G2 zeros(m, n); % -45 度 for i 2:m-1 for j 2:n-1 H1(i,j) -IMG(i-1,j-1) - IMG(i-1,j) - IMG(i-1,j1) ... IMG(i1,j-1) IMG(i1,j) IMG(i1,j1); H2(i,j) -IMG(i-1,j-1) - IMG(i,j-1) - IMG(i1,j-1) ... IMG(i-1,j1) IMG(i,j1) IMG(i1,j1); G1(i,j) -IMG(i-1,j-1) - IMG(i-1,j) - IMG(i,j-1) ... IMG(i,j1) IMG(i1,j) IMG(i1,j1); G2(i,j) -IMG(i-1,j-1) - IMG(i,j-1) - IMG(i1,j) ... IMG(i-1,j1) IMG(i,j1) IMG(i1,j1); end end % 合成总边缘直接累加四个方向响应 Z H1 H2 G1 G2; figure; imshow(Z, []);这里没有取绝对值或开平方直接用带符号的梯度值相加。因为 LoG 结果已经在边缘两侧产生正负交替的响应直接相加会让边缘两侧的极性被保留后续二值化时能通过阈值把边缘像素和背景分开。循环从 2 到 m-1、2 到 n-1是为了跳过图像最外层否则数组下标会越界。四个方向的响应幅度范围不同直接累加可能导致某个方向主导必要的话可以除以 4 做平均。作业里保留原始累加值因为后续二值化用的是全局迭代阈值幅度绝对值大小不会影响分割结果。3.3 为什么四个方向要叠加而不是取最大值有些边缘检测实现会取四个方向响应的最大值这样做可以抑制非边缘方向的干扰但也可能把斜边断裂成多个小段。圆的边缘是连续变化的任意一个切点附近只有一个方向响应最强其余方向响应偏弱。如果取最大值每个边缘像素只保留一个方向信息合成后的圆轮廓仍然是完整的但对于半径较大的圆切点附近的方向变化平缓最大值策略反而会放大噪声。累加策略保留了所有方向的响应边缘更粗后续 Hough 变换对边缘位置的容忍度更高。作业里得到的总边缘图边缘较粗但不破碎就是这个原因。代价是边缘粗会让 Hough 累加器出现多个相邻峰值后面需要通过阈值 p 来过滤非极大值。4. 自适应阈值二值化与 Hough 变换检测圆4.1 循环逼近全局阈值的过程合成边缘图 Z 的灰度分布不均匀用固定阈值切不出来完整圆。作业里用了一个迭代阈值法先以整图均值作为初始阈值 T把像素分成大于 T 和小于等于 T 两类分别求两类的均值 t1、t2新阈值取两个均值的平均数。不断重复直到前后两次阈值的差值绝对值小于 0.5。这个方法本质上是 Otsu 的简化版不需要计算类间方差收敛快但对初始值敏感。均值初始化的优点是鲁棒缺点是如果前景占比极小均值会偏向背景迭代次数增多。MATLAB 代码关键部分如下T sum(Z(:)) / (m * n); % 初始阈值整图均值 T2 20; % 初始化差值随便给个大数 while T2 0.5 Y1 0; Y2 0; N1 0; N2 0; for i 1:m for j 1:n if Z(i,j) T Y1 Y1 double(Z(i,j)); % 高值类累加 N1 N1 1; else Y2 Y2 double(Z(i,j)); % 低值类累加 N2 N2 1; end end end t1 Y1 / N1; t2 Y2 / N2; T1 (t1 t2) / 2; T2 abs(T - T1); T T1; end % 按最终阈值二值化边缘为黑背景为白 GS zeros(m, n); for i 1:m for j 1:n if Z(i,j) T GS(i,j) 0; else GS(i,j) 255; end end end代码里T2初始为 20只要大于 0.5 就继续迭代。终止条件应写成两次阈值的绝对差小于 0.5而不是某个固定迭代次数。循环里Y1/N1如果遇到 N1 为 0 会除零实际中几乎不会发生因为初始阈值是均值任何非均匀图像都有两侧像素但写成防御性代码更好。二值化时把大于阈值的设 0、小于等于阈值的设 255这样边缘像素是白色背景是黑色符合 Hough 变换中常见的输入约定。不过 Hough 变换本身不关心黑白哪边是边缘只要find(BW)能取到边缘点坐标就行。4.2 Hough 圆检测的参数空间映射Hough 变换检测圆的思路是把图像空间中的每个边缘点 ((x,y)) 映射到参数空间 ((a,b,r)) 中的一个三维锥面[ (a-x)^2 (b-y)^2 r^2 ]同一圆周上的所有边缘点在参数空间中对应的锥面会汇聚到同一个 ((a,b,r)) 点。实际实现时采用极坐标形式加快计算[ a x - r\cos\theta,\quad b y - r\sin\theta ]其中 (\theta) 是边缘点梯度方向角。对每个边缘像素、每个候选半径 r 和每个离散角度步长计算出可能的圆心坐标并在累加器hough_space(a,b,r)上加 1。遍历完所有边缘点后累加器中的局部最大值就对应真实圆心的位置和半径。作业中半径范围 10 到 100、半径步长 1、角度步长 (\pi/18)即 10°p0.7 表示累加峰值只有达到最大值的 70% 才被认为是候选圆。4.3 Hough 变换主循环与三种关键参数下面这段是核心实现注意我在原代码基础上加了坐标对齐和边界判断的注释BW double(GS); % 二值图边缘点为非零值 r_min 10; r_max 100; % 硬币半径范围 step_r 1; % 半径步长 step_angle pi / 18; % 角度步长 10 度 p 0.7; % 峰值阈值比例 size_r round((r_max - r_min) / step_r) 1; size_angle round(2 * pi / step_angle); hough_space zeros(m, n, size_r); % 三维累加器 [rows, cols] find(BW); % 所有边缘点坐标 ecount size(rows, 1); % 遍历每个边缘点、每个半径、每个角度方向 for i 1:ecount for r 1:size_r for k 1:size_angle a round(rows(i) - (r_min (r-1)*step_r) * cos(k*step_angle)); b round(cols(i) - (r_min (r-1)*step_r) * sin(k*step_angle)); if a 0 a m b 0 b n hough_space(a, b, r) hough_space(a, b, r) 1; end end end end % 找累加器最大值按 p 比例筛选候选圆心 max_para max(hough_space(:)); index find(hough_space max_para * p); % index 是线性索引需要拆回 a,b,r三个参数对结果的影响分别是参数取值影响step_r1越小检测越精细但累加器维度增大耗时成倍增加step_anglepi/18角度步长越密圆心定位越准但内层循环次数越多p0.7越大越严格只保留最明显的圆越小会输出很多假圆实际调试时先把 p 设到 0.8 以上看能否检到圆如果圆不完整再降 p。用for k1:length(index)把候选参数提取出来后还要做一次边缘点归属判定对每个边缘点计算它与候选圆心距离是否在半径的 ±5 像素范围内是则保留该边缘像素。这个距离容差 (\pm 5) 是经验值和半径步长 1 配合能容忍边缘粗造成的偏差。注意代码里length是内置函数名用它当变量名会在后续循环里覆盖函数工程上应改成num_cand之类。5. 结果叠加与梯度加速验证最后一环是把 Hough 检测到的圆叠加到原灰度图上。检测结果hough_circle是逻辑图像值为 1 的位置就是识别出的圆周像素。叠加策略是逐像素判断如果是圆周点则置 255白否则保留原灰度值这样圆的轮廓会被高亮出来同时底图细节不丢失。LST zeros(m, n); for i 1:m for j 1:n if hough_circle(i, j) 1 LST(i, j) 255; else LST(i, j) S(i, j); end end end figure; imshow(LST, []); title(叠加后的圆检测结果);叠加前最好确认S是 uint8 而LST是 double如果不做类型转换MATLAB 会隐式把 uint8 和 double 略算可能出现超出 255 的截断。上面代码先声明LST zeros(m,n)是 double再把S(i,j)直接赋进去MATLAB 会自动转换但稳妥起见用double(S(i,j))更清晰。叠加后如果圆太多或太少先检查二值化方向是否反了再检查累加器的max_para*p是否把真实圆的峰值卡掉了。当你拿到一副硬币图可以从一个更快的梯度方向法切入来验证结果。本文开头说过圆的梯度方向要么指向圆心、要么背离圆心所以不需要让每个边缘点对全部半径做三维锥面搜索而是先利用梯度角预估计圆心所在的射线方向。把极坐标公式改写成[ a x \pm r\cos\theta(x,y),\quad b y \pm r\sin\theta(x,y) ]这样每个边缘点只需要沿梯度方向正反两个方向对每个半径做一次累加。相比三重循环计算量从 (O(边缘点数 \times 半径数 \times 角度数)) 降到 (O(边缘点数 \times 半径数 \times 2))。作业里半径 10–100、步长 1 有 91 个候选半径角度取 36 个方向三重循环复杂度是 91×36 倍。硬币图像边缘点通常在几千个量级纯 MATLAB 三重循环会跑几分钟而梯度加速版本能在 10 秒内完成。修改的核心是把for k1:size_angle这一层去掉改为读取边缘点处预计算的梯度方向。验证结果时不要只看叠加图打印候选圆的中心和半径检查同一区域是否出现多个相近的圆。比如中心坐标差小于 3、半径差小于 2就认为是同一个圆的重复响应保留累加值最大的那个即可。另外如果检测到的圆明显偏移硬币中心通常是 LoG 核中心偏移导致边缘位置整体平移了几像素回到第 2 章检查h的坐标原点。最后把检测结果保存为二值掩膜用regionprops统计每个连通域的质心、面积和圆形度可以量化评估检测精度。圆形度指标 (4\pi A/P^2) 越接近 1说明检测到的圆越接近理想圆这个指标比肉眼判断可靠得多。本文还有配套的精品资源点击获取 
长河君 · 2026/9/19 14:39:39 · 阅读 1.2k MATLAB实现Hough圆检测：LoG边缘提取与Prewitt边缘融合的完整流程
简介面向数字图像处理课程期末考核的一份完整PDF报告聚焦图像平滑、边缘检测、二值化与Hough变换检测圆形四大核心环节可作为高年级本科生与研究生完成同类大作业的设计参考。资源包为1个PDF文件大小1.09MB已有850人学习浏览。报告以硬币图像为输入完整覆盖实验流程使用方差为1的高斯滤波器对图像进行平滑与拉普拉斯预处理采用Prewitt梯度算子提取水平、垂直、45度和-45度四个方向的边缘通过循环确定全局阈值完成二值化再借助Hough变换在参数空间中检测圆的圆心与半径并将检测到的圆形叠加回原始灰度图像实现边缘增强。正文还给出了算法原理公式、实验流程图、中间结果和最终结果截图以及对应的实现代码方便读者对照复现。整体结构清晰既适合期末大作业的完整参考也适合图像处理课程中Hough变换、边缘检测等知识点的专项复习。1. 从一副硬币图到 Hough 圆检测的完整链路很多人在数字图像处理课程里都能把 Fourier 变换、Prewitt 算子、Hough 变换的公式背得很熟但真拿到一副硬币图像、要求检测其中的圆并叠加回原图时往往会卡在“先做什么、后做什么、参数怎么设”上。这篇博文不是把课程报告复述一遍而是以一个实际的数字图像处理期末大作业为样本把「高斯平滑 → 拉普拉斯变换 → 四方向 Prewitt 边缘检测 → 循环阈值二值化 → Hough 变换检测圆 → 与原图叠加增强」这条完整链路拆开讲清楚。适合正在做图像处理课程设计的学生也适合需要快速上手 MATLAB 实现圆形检测的工程师。文中所有代码都基于 MATLAB参数会逐一解释最后还会给出我在调试时踩过的坑。2. LoG 预处理为什么平滑和拉普拉斯要合并成一次卷积2.1 先平滑再求导的数学逻辑拿到灰度图像后直接做边缘检测噪声会把边缘响应搅得很难看。标准做法是先平滑、再求二阶导也就是 LoGLaplacian of Gaussian算子。这里的关键是把两步合成一步先用高斯函数对图像做卷积再对结果做拉普拉斯变换即[ \nabla^2[h(r)*f(x,y)] [\nabla^2 h(r)] * f(x,y) ]其中 (h(r)\frac{1}{2\pi\sigma^2}e^{-r^2/(2\sigma^2)})而 (\nabla^2 h(r)\frac{r^2-2\sigma^2}{\sigma^4}e^{-r^2/(2\sigma^2)})。也就是说先平滑再拉普拉斯等价于直接用一个 LoG 核与图像卷积省一次卷积运算。作业里取 (\sigma^21)核的尺寸和整个图像一样大然后用傅里叶变换做频域乘法避免时域大核卷积的耗时。2.2 频域实现时最容易错的坐标原点问题直接看原文的 MATLAB 代码h 矩阵的生成方式是t i^2 j^2这是以图像左上角为原点的距离平方。这样做出来的 h 不是真正中心对称的 LoG 核因为 LoG 核应该以核中心为原点。常见做法是先把坐标中心移到图像中点再计算距离S imread(Sweden_coins.bmp); [m, n] size(S); sigma2 1; % 方差 sigma^2 1 h zeros(m, n); for i 1:m for j 1:n x i - (m 1) / 2; % 将坐标原点移到核中心 y j - (n 1) / 2; r2 x^2 y^2; h(i, j) (r2 - 2 * sigma2) * exp(-r2 / (2 * sigma2)) / sigma2^2; end end % 频域滤波逐点相乘再反变换 H fftshift(fft2(h)); G fftshift(fft2(S)); W H .* G; w ifft2(ifftshift(W)); figure; imshow(w, []);这里fftshift把零频移到中心目的是让频域相乘时低频对齐。做完ifft2后要再用ifftshift调回来否则图像会被平移错位。imshow 默认会把负值截断到 0所以显示时最好加[]让图像自动拉伸灰度范围不然洛普拉斯结果的负响应会显示成全黑。我一般还会直接用imagesc(w)配合colormap gray看细节比imshow更直观。2.3 核尺寸与边界处理的影响作业里核和原图等大循环计算量是 (m \times n) 次指数运算对一副几百像素的图像还能接受但再大就会很慢。实际工程中 LoG 核通常取 5×5 或 7×7 就够了因为高斯尾巴衰减快。核太小会让平滑不足核太大则边缘定位偏差明显。边界像素在卷积时会有响应异常MATLAB 的conv2默认补零而频域乘法相当于周期延拓两种方式在图像边缘会有几像素的差异。如果后续 Hough 检测对边缘位置敏感建议先裁掉边界 5 像素再做后续处理。另一个细节是这里的 h 没有归一化LoG 核的累加和接近 0对直流分量几乎没有响应所以滤波结果整体均值也会接近 0这是正常的。3. 四方向 Prewitt 边缘检测与合成3.1 为什么选 Prewitt 而不是 SobelPrewitt 和 Sobel 的区别在于平滑权重Sobel 对中心行/列给了 2 倍权重更强调中心像素Prewitt 所有权重都是 1结构更简单各向异性响应也更均匀。在 Hough 变换检测圆这个场景里我们对边缘的连续性和方向覆盖比单点精度更在意Prewitt 的四方向掩膜可以直接覆盖 0°、45°、90°、135° 四个方向。Sobel 虽然对噪声更鲁棒但它的平滑核只对水平和垂直方向做了对称加权斜方向掩膜在离散网格上的定义不如 Prewitt 直观。3.2 四个方向的掩膜与卷积实现Prewitt 四个方向的 3×3 掩膜如下方向掩膜模板水平检测垂直边缘[-1 -1 -1; 0 0 0; 1 1 1]垂直检测水平边缘[-1 0 1; -1 0 1; -1 0 1]45°[-1 -1 0; -1 0 1; 0 1 1]-45°[0 -1 -1; 1 0 -1; 1 1 0]注意方向命名在不少教材里是反的有的把 H1 叫水平边缘指的其实是水平方向上的灰度变化也就是检测垂直边缘。作业中四个掩膜与图像逐点卷积后直接取绝对值相加。代码写法如下IMG double(w); % 上一步的 LoG 结果转 double 避免溢出 H1 zeros(m, n); % 水平方向梯度垂直边缘 H2 zeros(m, n); % 垂直方向梯度水平边缘 G1 zeros(m, n); % 45 度 G2 zeros(m, n); % -45 度 for i 2:m-1 for j 2:n-1 H1(i,j) -IMG(i-1,j-1) - IMG(i-1,j) - IMG(i-1,j1) ... IMG(i1,j-1) IMG(i1,j) IMG(i1,j1); H2(i,j) -IMG(i-1,j-1) - IMG(i,j-1) - IMG(i1,j-1) ... IMG(i-1,j1) IMG(i,j1) IMG(i1,j1); G1(i,j) -IMG(i-1,j-1) - IMG(i-1,j) - IMG(i,j-1) ... IMG(i,j1) IMG(i1,j) IMG(i1,j1); G2(i,j) -IMG(i-1,j-1) - IMG(i,j-1) - IMG(i1,j) ... IMG(i-1,j1) IMG(i,j1) IMG(i1,j1); end end % 合成总边缘直接累加四个方向响应 Z H1 H2 G1 G2; figure; imshow(Z, []);这里没有取绝对值或开平方直接用带符号的梯度值相加。因为 LoG 结果已经在边缘两侧产生正负交替的响应直接相加会让边缘两侧的极性被保留后续二值化时能通过阈值把边缘像素和背景分开。循环从 2 到 m-1、2 到 n-1是为了跳过图像最外层否则数组下标会越界。四个方向的响应幅度范围不同直接累加可能导致某个方向主导必要的话可以除以 4 做平均。作业里保留原始累加值因为后续二值化用的是全局迭代阈值幅度绝对值大小不会影响分割结果。3.3 为什么四个方向要叠加而不是取最大值有些边缘检测实现会取四个方向响应的最大值这样做可以抑制非边缘方向的干扰但也可能把斜边断裂成多个小段。圆的边缘是连续变化的任意一个切点附近只有一个方向响应最强其余方向响应偏弱。如果取最大值每个边缘像素只保留一个方向信息合成后的圆轮廓仍然是完整的但对于半径较大的圆切点附近的方向变化平缓最大值策略反而会放大噪声。累加策略保留了所有方向的响应边缘更粗后续 Hough 变换对边缘位置的容忍度更高。作业里得到的总边缘图边缘较粗但不破碎就是这个原因。代价是边缘粗会让 Hough 累加器出现多个相邻峰值后面需要通过阈值 p 来过滤非极大值。4. 自适应阈值二值化与 Hough 变换检测圆4.1 循环逼近全局阈值的过程合成边缘图 Z 的灰度分布不均匀用固定阈值切不出来完整圆。作业里用了一个迭代阈值法先以整图均值作为初始阈值 T把像素分成大于 T 和小于等于 T 两类分别求两类的均值 t1、t2新阈值取两个均值的平均数。不断重复直到前后两次阈值的差值绝对值小于 0.5。这个方法本质上是 Otsu 的简化版不需要计算类间方差收敛快但对初始值敏感。均值初始化的优点是鲁棒缺点是如果前景占比极小均值会偏向背景迭代次数增多。MATLAB 代码关键部分如下T sum(Z(:)) / (m * n); % 初始阈值整图均值 T2 20; % 初始化差值随便给个大数 while T2 0.5 Y1 0; Y2 0; N1 0; N2 0; for i 1:m for j 1:n if Z(i,j) T Y1 Y1 double(Z(i,j)); % 高值类累加 N1 N1 1; else Y2 Y2 double(Z(i,j)); % 低值类累加 N2 N2 1; end end end t1 Y1 / N1; t2 Y2 / N2; T1 (t1 t2) / 2; T2 abs(T - T1); T T1; end % 按最终阈值二值化边缘为黑背景为白 GS zeros(m, n); for i 1:m for j 1:n if Z(i,j) T GS(i,j) 0; else GS(i,j) 255; end end end代码里T2初始为 20只要大于 0.5 就继续迭代。终止条件应写成两次阈值的绝对差小于 0.5而不是某个固定迭代次数。循环里Y1/N1如果遇到 N1 为 0 会除零实际中几乎不会发生因为初始阈值是均值任何非均匀图像都有两侧像素但写成防御性代码更好。二值化时把大于阈值的设 0、小于等于阈值的设 255这样边缘像素是白色背景是黑色符合 Hough 变换中常见的输入约定。不过 Hough 变换本身不关心黑白哪边是边缘只要find(BW)能取到边缘点坐标就行。4.2 Hough 圆检测的参数空间映射Hough 变换检测圆的思路是把图像空间中的每个边缘点 ((x,y)) 映射到参数空间 ((a,b,r)) 中的一个三维锥面[ (a-x)^2 (b-y)^2 r^2 ]同一圆周上的所有边缘点在参数空间中对应的锥面会汇聚到同一个 ((a,b,r)) 点。实际实现时采用极坐标形式加快计算[ a x - r\cos\theta,\quad b y - r\sin\theta ]其中 (\theta) 是边缘点梯度方向角。对每个边缘像素、每个候选半径 r 和每个离散角度步长计算出可能的圆心坐标并在累加器hough_space(a,b,r)上加 1。遍历完所有边缘点后累加器中的局部最大值就对应真实圆心的位置和半径。作业中半径范围 10 到 100、半径步长 1、角度步长 (\pi/18)即 10°p0.7 表示累加峰值只有达到最大值的 70% 才被认为是候选圆。4.3 Hough 变换主循环与三种关键参数下面这段是核心实现注意我在原代码基础上加了坐标对齐和边界判断的注释BW double(GS); % 二值图边缘点为非零值 r_min 10; r_max 100; % 硬币半径范围 step_r 1; % 半径步长 step_angle pi / 18; % 角度步长 10 度 p 0.7; % 峰值阈值比例 size_r round((r_max - r_min) / step_r) 1; size_angle round(2 * pi / step_angle); hough_space zeros(m, n, size_r); % 三维累加器 [rows, cols] find(BW); % 所有边缘点坐标 ecount size(rows, 1); % 遍历每个边缘点、每个半径、每个角度方向 for i 1:ecount for r 1:size_r for k 1:size_angle a round(rows(i) - (r_min (r-1)*step_r) * cos(k*step_angle)); b round(cols(i) - (r_min (r-1)*step_r) * sin(k*step_angle)); if a 0 a m b 0 b n hough_space(a, b, r) hough_space(a, b, r) 1; end end end end % 找累加器最大值按 p 比例筛选候选圆心 max_para max(hough_space(:)); index find(hough_space max_para * p); % index 是线性索引需要拆回 a,b,r三个参数对结果的影响分别是参数取值影响step_r1越小检测越精细但累加器维度增大耗时成倍增加step_anglepi/18角度步长越密圆心定位越准但内层循环次数越多p0.7越大越严格只保留最明显的圆越小会输出很多假圆实际调试时先把 p 设到 0.8 以上看能否检到圆如果圆不完整再降 p。用for k1:length(index)把候选参数提取出来后还要做一次边缘点归属判定对每个边缘点计算它与候选圆心距离是否在半径的 ±5 像素范围内是则保留该边缘像素。这个距离容差 (\pm 5) 是经验值和半径步长 1 配合能容忍边缘粗造成的偏差。注意代码里length是内置函数名用它当变量名会在后续循环里覆盖函数工程上应改成num_cand之类。5. 结果叠加与梯度加速验证最后一环是把 Hough 检测到的圆叠加到原灰度图上。检测结果hough_circle是逻辑图像值为 1 的位置就是识别出的圆周像素。叠加策略是逐像素判断如果是圆周点则置 255白否则保留原灰度值这样圆的轮廓会被高亮出来同时底图细节不丢失。LST zeros(m, n); for i 1:m for j 1:n if hough_circle(i, j) 1 LST(i, j) 255; else LST(i, j) S(i, j); end end end figure; imshow(LST, []); title(叠加后的圆检测结果);叠加前最好确认S是 uint8 而LST是 double如果不做类型转换MATLAB 会隐式把 uint8 和 double 略算可能出现超出 255 的截断。上面代码先声明LST zeros(m,n)是 double再把S(i,j)直接赋进去MATLAB 会自动转换但稳妥起见用double(S(i,j))更清晰。叠加后如果圆太多或太少先检查二值化方向是否反了再检查累加器的max_para*p是否把真实圆的峰值卡掉了。当你拿到一副硬币图可以从一个更快的梯度方向法切入来验证结果。本文开头说过圆的梯度方向要么指向圆心、要么背离圆心所以不需要让每个边缘点对全部半径做三维锥面搜索而是先利用梯度角预估计圆心所在的射线方向。把极坐标公式改写成[ a x \pm r\cos\theta(x,y),\quad b y \pm r\sin\theta(x,y) ]这样每个边缘点只需要沿梯度方向正反两个方向对每个半径做一次累加。相比三重循环计算量从 (O(边缘点数 \times 半径数 \times 角度数)) 降到 (O(边缘点数 \times 半径数 \times 2))。作业里半径 10–100、步长 1 有 91 个候选半径角度取 36 个方向三重循环复杂度是 91×36 倍。硬币图像边缘点通常在几千个量级纯 MATLAB 三重循环会跑几分钟而梯度加速版本能在 10 秒内完成。修改的核心是把for k1:size_angle这一层去掉改为读取边缘点处预计算的梯度方向。验证结果时不要只看叠加图打印候选圆的中心和半径检查同一区域是否出现多个相近的圆。比如中心坐标差小于 3、半径差小于 2就认为是同一个圆的重复响应保留累加值最大的那个即可。另外如果检测到的圆明显偏移硬币中心通常是 LoG 核中心偏移导致边缘位置整体平移了几像素回到第 2 章检查h的坐标原点。最后把检测结果保存为二值掩膜用regionprops统计每个连通域的质心、面积和圆形度可以量化评估检测精度。圆形度指标 (4\pi A/P^2) 越接近 1说明检测到的圆越接近理想圆这个指标比肉眼判断可靠得多。本文还有配套的精品资源点击获取