基于Matlab的裂缝及长度检测——“算法实现与应用”

基于Matlab的裂缝及其长度检测

工地朋友上周甩我三张照片的时候脸都绿了，说监理要这个季度所有老旧混凝土路的裂缝台账，照片是拍了，但一条一条量尺量腿都跑断，还要拍视频证明量尺的位置，麻烦到想辞职。

基于Matlab的裂缝及其长度检测

我翻了翻照片，果然头大——下午4点的逆光加小石子反光，路面还有风吹沙粒蹭的模糊噪点，细裂缝躲在大光斑边缘若隐若现，差点要瞎。不过转念一想这不是Matlab图像处理的经典入门练手题嘛（当年我大作业差点选这个，后来嫌找不到好看的工地实拍图改做表情包边缘检测了），赶紧打开电脑撸代码。

第一步：先把“脏照片”洗干净

朋友的是1280×720的普通手机直出JPG，先读取进来转灰度图，省得RGB三个通道瞎忙活。

% 读取并转灰度 img = imread('crack_dirty.jpg'); gray = rgb2gray(img); figure, subplot(231), imshow(gray), title('原图灰度版');

转完看直方图：低像素区（暗，可能是裂缝）只有一点点尾巴，中间高像素区（亮，路面）挤成一坨，还有大的尖峰（逆光大光斑）。

去噪组合拳：先中值后高斯

手机拍的野外图一般是椒盐混高斯，沙粒反光是椒盐，手机摄像头CMOS的热噪是高斯。先压椒盐，用medfilt2，别用3×3太碎留不下痕迹，也别7×7把头发丝一样的细裂缝磨没，5×5踩中间刚合适。

% 中值滤波去椒盐 med = medfilt2(gray, [5 5]); subplot(232), imshow(med), title('中值去椒盐');

椒盐没了，再来个低通高斯压热噪，kernel用fspecial('gaussian', [3 3], 0.8)，sigma调0.8刚好够软，别太猛把刚才中值压出来的沙粒轮廓（不对要的是裂缝轮廓）一起糊了。

% 高斯低通去热噪 gauss_kernel = fspecial('gaussian', [3 3], 0.8); gauss = imfilter(med, gauss_kernel, 'symmetric'); % symmetric防边缘变黑 subplot(233), imshow(gauss), title('低通去热噪');

第二步：把躲猫猫的裂缝揪出来

压完噪对比度还是烂，imadjust全局调直方图不行，大光斑过曝路面暗部更暗，adapthisteq自适应的就好用多了（也就是我们常说的CLAHE）。tile选[16 16]，别太碎（比如8×8）把光斑切成马赛克，也别太大（32×32）适配不了不均匀的路面反光；clipLimit调0.02，防止过度放大细裂纹旁边的沙粒反光。

% CLAHE自适应增强 enhance = adapthisteq(gauss, 'NumTiles', [16 16], 'ClipLimit', 0.02); subplot(234), imshow(enhance), title('CLAHE增强');

增强完看，细裂缝终于从大光斑边上冒出来了！

第三步：给裂缝画个黑白漫画

二值化是重头戏，全局Otsu直接阵亡，试一下adaptthresh的局部阈值，选'local'的高斯加权平均（默认是'mean'，加权的对中间暗两边亮的光斑更友好），neighborhood选[25 25]，刚好覆盖大概直径50像素的小石子反光；sensitivity调0.6，别太低漏了大宽缝，别太高把路面沙粒影子全抠出来。

% 局部自适应二值化 level = adaptthresh(enhance, [25 25], 'Statistic', 'gaussian', 'Sensitivity', 0.6); bw = imbinarize(enhance, level); subplot(235), imshow(bw), title('初步二值化');

初步二值化肯定不干净：断缝、小沙粒噪点、刚才加权没完全压的反光斑。

组合除渣大法

1. 先剪小的碎渣：用bwareaopen，把小于20像素的连通域全干掉——朋友照片的细裂纹大概宽度2像素，长度最少也得10像素（2×10=20），先留有用的。

% 去掉小连通域（小于20像素） bw_clean1 = bwareaopen(bw, 20);

1. 再连断的缝：用闭运算（先膨胀后腐蚀），strel选四个主要方向的直线（0°、45°、90°、135°），半径3像素——朋友照片的断缝一般在5像素以内，别太猛不然连到旁边的反光斑。

% 闭运算连断缝（四个方向） se1 = strel('line', 3, 0); se2 = strel('line', 3, 45); se3 = strel('line', 3, 90); se4 = strel('line', 3, 135); bw_clean2 = imclose(bw_clean1, se1); bw_clean2 = imclose(bw_clean2, se2); bw_clean2 = imclose(bw_clean2, se3); bw_clean2 = imclose(bw_clean2, se4);

1. 再剪刚才连错的渣：用bwareaopen，阈值提到50——真裂缝连完应该至少50像素了，两个小沙粒反光连起来一般也超不过40。

% 再次去掉小连通域（小于50像素） bw_final = bwareaopen(bw_clean2, 50); subplot(236), imshow(bw_final), title('最终二值化');

第四步：把粗线条变“骨架”，算长度

算像素长度必须要单像素宽的骨架，用bwmorph的'thin'操作，Inf次直到不变。

% 细化裂缝成单像素骨架 bw_skel = bwmorph(bw_final, 'thin', Inf); figure, imshowpair(bw_final, bw_skel, 'montage'), title('最终二值 vs 单像素骨架');

算像素长度的“笨办法但绝对准”

之前试了regionprops的'Perimeter'除以2，发现对直线型的还行，对折线型的有点误差；试了'MinorAxisMajor'拟合椭圆，完全不对。后来想到regionprops的'PixelList'可以得到所有坐标，但不是连续的轨迹！换个思路：用bwtraceboundary沿着骨架的端点追踪连续的像素点，然后用pdist算相邻两点的欧氏距离，最后求和！

% 找到所有骨架的连通域 stats = regionprops(bw_skel, 'PixelList', 'Centroid'); num_cracks = length(stats); total_pixel_length = 0; figure, imshow(gray), hold on; % 在原图上画骨架，标注长度 % 遍历每个连通域 for i = 1:num_cracks % 找连通域的第一个端点（PixelList的第一个点可能是中间点，bwdist找最远的两个？选第一个就行） % 这里简化一下，直接用bwboundaries的默认追踪？不对bwboundaries是找外边界，单像素线的外边界也是自己？ % 哦更简单的是bwtraceboundary找一个连通域的所有点，选8连通 [x, y] = find(bw_skel, 1); boundary = bwtraceboundary(bw_skel, [x y], 'N', 8); % 算相邻两点的欧氏距离 distances = pdist(boundary, 'euclidean'); pixel_length = sum(distances); total_pixel_length = total_pixel_length + pixel_length; % 在原图上画这个裂缝的骨架，标上序号和像素长度 plot(boundary(:,2), boundary(:,1), 'r-', 'LineWidth', 2); text(stats(i).Centroid(1), stats(i).Centroid(2), ... sprintf('#%d: %.1fpx', i, pixel_length), ... 'Color', 'yellow', 'FontSize', 10, 'FontWeight', 'bold'); % 把这个连通域从bw_skel里去掉，防止下次重复找 bw_skel(sub2ind(size(bw_skel), boundary(:,1), boundary(:,2))) = 0; end hold off; fprintf('总像素长度：%.2fpx\n', total_pixel_length);

这个代码我调了好久，一开始没把追踪过的连通域去掉，每次都找同一个，后来加了bw_skel(sub2ind(...)) = 0就搞定了。

第五步：像素长度转真实长度

最后一步也是朋友最关心的，怎么把px换成cm/m？这个需要一个参考系——朋友拍照片的时候特意在裂缝旁边放了一把50cm的直尺，直尺在照片里的长度是多少px呢？我用ginput手动选了直尺的两个端点，算一下欧氏距离。

% 手动选直尺的两个端点 figure, imshow(gray), title('请点击直尺的两个端点，按Enter结束'); [ruler_pts, ~] = ginput(2); ruler_pixel_length = pdist(ruler_pts, 'euclidean'); ruler_real_length = 50; % 单位cm % 换算比例 scale = ruler_real_length / ruler_pixel_length; % cm/px % 总真实长度 total_real_length = total_pixel_length * scale; fprintf('总真实长度：%.2fcm = %.2fm\n', total_real_length, total_real_length/100);

朋友试了三把照片，手动选直尺和量尺量出来的真实长度误差在5%以内，监理完全接受！他现在逢人就吹我是“工地裂缝检测大师”，还说下次土方开挖的边坡位移检测也要找我——这个我可不敢随便接，位移检测要用连续帧的光流法，难度要大得多。

最后附一张处理后的效果截图，红色的是单像素骨架，黄色的是序号和像素长度，是不是很直观？