简介:这份资源面向计算机视觉入门者、地质图像分析方向的学生以及需要完成期末大作业或课程设计的学习者,提供一套基于Python的CT岩芯与岩石裂缝语义分割完整方案,用于解决像素级裂隙识别与量化分析问题。压缩包共15个文件,约1.15MB,包含3个py脚本用于数据增强与均值计算,6张jpg图像作为岩石、混凝土及CT扫描的原始图与标注掩码,另有zbak备份、md说明文档和gitignore等辅助文件,结构紧凑便于直接运行与二次修改。资源覆盖从图像读取、形态学处理到深度学习模型搭建的端到端流程,可帮助读者理解语义分割在地质非破坏性检测中的落地方式,掌握裂隙分布模式识别与定量分析思路。目前已有72人学习下载,适合作为课程实践与算法入门的参考素材。
1. CT岩心裂缝语义分割:从灰度切片到可训练掩码的完整链路
拿到一批工业 CT 扫描的岩心切片,想把里面的裂缝自动勾出来,这件事在石油地质和岩石力学圈子里需求很实在。人工在几百张切片上逐像素描裂缝,一张图少则十几分钟,多则半小时,标注一致性还差。基于 Python 的 CT 岩心与岩石裂缝语义分割系统,本质就是把「切片预处理 → 裂缝像素级标注 → 语义分割模型训练 → 推理出掩码」这条链路用代码串起来,让裂缝从灰度图里被自动分离成二值掩码。它适合两类人:一类是手里有 CT 数据、想快速搭一套能跑的分割流程的地质或岩土工程师;另一类是刚接触语义分割、想找一个真实工业场景练手的 Python 开发者。这套东西不追求 SOTA 指标,追求的是数据能进、模型能训、掩码能出、结果能复核。下面按我实际搭过的顺序,把每个环节的参数和坑讲清楚。
2. 数据准备:CT 切片怎么变成语义分割能吃的格式
CT 岩心数据通常是 16 位灰度 TIFF 或 DICOM 序列,灰度范围能到 0–65535,而裂缝在灰度上往往只比基质暗一点点,对比度极低。直接丢给模型,网络学不到东西。所以第一步不是写模型,是把数据整理成「图像 + 掩码」成对、灰度归一、尺寸统一的格式。语义分割和实例分割的区别在这里很关键:裂缝分割只关心「这个像素是不是裂缝」,不关心「这是第几条裂缝」,所以用语义分割的标签体系,一张图对应一张单通道掩码,像素值 0 是背景、1 是裂缝,不需要给每条裂缝编号。这也是为什么标题里写的是语义分割而不是实例分割——岩心裂缝的工程诉求是统计裂缝面积占比和走向,不是数裂缝条数。
2.1 从 CT 原始切片到 PNG 的批量转换
工业 CT 导出的切片常见是 16 位无符号整型,Python 里用tifffile或pydicom读。转成 8 位 PNG 是为了后续标注工具和大部分分割框架好处理,但直接线性压缩会丢掉裂缝和基质的微小灰度差,所以要先做对比度拉伸。下面这段是我常用的转换脚本。
import numpy as np import tifffile from PIL import Image import os def ct_slice_to_png(src_dir, dst_dir, low_pct=1, high_pct=99): os.makedirs(dst_dir, exist_ok=True) for name in sorted(os.listdir(src_dir)): if not name.lower().endswith(('.tif', '.tiff')): continue img = tifffile.imread(os.path.join(src_dir, name)).astype(np.float32) # 按百分位裁剪,避免个别极亮/极暗像素拉垮整体对比度 lo, hi = np.percentile(img, [low_pct, high_pct]) img = np.clip(img, lo, hi) img = (img - lo) / (hi - lo + 1e-6) * 255.0 img = img.astype(np.uint8) Image.fromarray(img).save(os.path.join(dst_dir, name.rsplit('.', 1)[0] + '.png')) ct_slice_to_png('./ct_raw', './ct_png')逻辑说明:np.percentile取 1% 和 99% 分位作为拉伸上下限,比直接用 min/max 稳,因为 CT 里常有金属矿物或扫描伪影造成的极端亮斑。参数low_pct和high_pct是可调的,裂缝对比度还是不够时可以收到 5/95,让拉伸更激进。转换后务必抽查几张,确认裂缝没有在拉伸中被压没。
2.2 掩码标注规范与目录结构
标注工具用 LabelMe 或 CVAT 都行,导出成 PNG 掩码。关键是标注规范要统一:裂缝边缘怎么算、宽度小于 2 像素的细缝标不标、裂缝和孔洞连在一起时怎么切分。我一般定三条规则:宽度小于 2 像素的裂缝不标(模型学不稳,标注也不一致);裂缝与孔洞连通时,只标裂缝主体,孔洞归背景;掩码边缘允许 1 像素误差。目录结构按下面组织,训练脚本直接按文件名配对。
dataset/ images/ slice_0001.png slice_0002.png masks/ slice_0001.png slice_0002.png掩码必须是单通道、像素值只有 0 和 1(或 0 和 255),如果标注工具导出的是彩色索引图,要转成灰度二值。这一步不做,训练时 loss 会直接报类别数不匹配。
2.3 数据集划分与增强的边界
按切片划分训练/验证/测试,比例 7:2:1。注意不能随机打乱所有切片再分,因为相邻切片高度相似,随机分会导致验证集里出现和训练集几乎一样的图,指标虚高。正确做法是按深度区间分块,比如前 70% 深度做训练,中间 20% 验证,最后 10% 测试。增强用水平翻转、垂直翻转、±10° 旋转、亮度抖动就够了。CT 岩心是各向异性的,过度旋转会破坏裂缝的真实走向分布,所以旋转角度别开太大。
3. 语义分割模型选型与训练:U-Net 为什么还是裂缝分割的稳妥起点
裂缝分割是典型的「细长目标 + 低对比度」任务,目标像素占比往往不到 5%,属于强类别不平衡。这个场景下,U-Net 系列依然是性价比最高的起点,编码器-解码器加跳跃连接的结构对细结构的保留能力好,训练数据需求也比 Transformer 类模型低。如果数据量上千张、算力充足,可以上 SegFormer 或 DeepLabV3+,但对大多数岩心项目,U-Net 加一个预训练 ResNet 编码器就够用。下面给一套能直接跑的 PyTorch 训练代码。
3.1 U-Net 模型定义与损失函数
import torch import torch.nn as nn import torchvision class UNet(nn.Module): def __init__(self, pretrained=True): super().__init__() resnet = torchvision.models.resnet34(weights='IMAGENET1K_V1' if pretrained else None) self.enc0 = nn.Sequential(resnet.conv1, resnet.bn1, resnet.relu) # 1/2 self.enc1 = nn.Sequential(resnet.maxpool, resnet.layer1) # 1/4 self.enc2 = resnet.layer2 # 1/8 self.enc3 = resnet.layer3 # 1/16 self.enc4 = resnet.layer4 # 1/32 self.up4 = nn.ConvTranspose2d(512, 256, 2, stride=2) self.dec4 = nn.Sequential(nn.Conv2d(512, 256, 3, padding=1), nn.BatchNorm2d(256), nn.ReLU()) self.up3 = nn.ConvTranspose2d(256, 128, 2, stride=2) self.dec3 = nn.Sequential(nn.Conv2d(256, 128, 3, padding=1), nn.BatchNorm2d(128), nn.ReLU()) self.up2 = nn.ConvTranspose2d(128, 64, 2, stride=2) self.dec2 = nn.Sequential(nn.Conv2d(128, 64, 3, padding=1), nn.BatchNorm2d(64), nn.ReLU()) self.up1 = nn.ConvTranspose2d(64, 32, 2, stride=2) self.dec1 = nn.Sequential(nn.Conv2d(64, 32, 3, padding=1), nn.BatchNorm2d(32), nn.ReLU()) self.head = nn.Conv2d(32, 1, 1) def forward(self, x): e0 = self.enc0(x) e1 = self.enc1(e0) e2 = self.enc2(e1) e3 = self.enc3(e2) e4 = self.enc4(e3) d4 = self.dec4(torch.cat([self.up4(e4), e3], dim=1)) d3 = self.dec3(torch.cat([self.up3(d4), e2], dim=1)) d2 = self.dec2(torch.cat([self.up2(d3), e1], dim=1)) d1 = self.dec1(torch.cat([self.up1(d2), e0], dim=1)) return self.head(d1)逻辑说明:编码器用 ResNet34 预训练权重,跳跃连接把编码器各层特征拼到解码器对应层,这是 U-Net 保留细裂缝的关键。head输出单通道 logits,配合下面的损失函数。参数上,pretrained=True在数据少于 500 张时强烈建议开,能明显加快收敛。
损失函数用 Dice + BCE 组合,Dice 负责应对类别不平衡,BCE 负责像素级稳定梯度。
class DiceBCELoss(nn.Module): def __init__(self, dice_weight=0.5): super().__init__() self.dice_weight = dice_weight self.bce = nn.BCEWithLogitsLoss() def forward(self, logits, targets): bce = self.bce(logits, targets) probs = torch.sigmoid(logits) intersection = (probs * targets).sum() dice = 1 - (2 * intersection + 1e-6) / (probs.sum() + targets.sum() + 1e-6) return self.dice_weight * dice + (1 - self.dice_weight) * bcedice_weight默认 0.5,如果裂缝像素占比低于 2%,可以提到 0.7,让模型更关注裂缝。但别设成 1.0,纯 Dice 在训练初期梯度不稳,容易震荡。
3.2 训练循环与关键超参
from torch.utils.data import Dataset, DataLoader from PIL import Image import numpy as np class CrackDataset(Dataset): def __init__(self, img_dir, mask_dir, size=512): self.img_dir, self.mask_dir, self.size = img_dir, mask_dir, size self.names = sorted(os.listdir(img_dir)) def __len__(self): return len(self.names) def __getitem__(self, idx): name = self.names[idx] img = Image.open(os.path.join(self.img_dir, name)).convert('L').resize((self.size, self.size)) mask = Image.open(os.path.join(self.mask_dir, name)).convert('L').resize((self.size, self.size)) img = np.array(img, dtype=np.float32) / 255.0 mask = (np.array(mask) > 127).astype(np.float32) return torch.from_numpy(img)[None], torch.from_numpy(mask)[None] train_ds = CrackDataset('./dataset/images', './dataset/masks') loader = DataLoader(train_ds, batch_size=4, shuffle=True, num_workers=2) model = UNet().cuda() criterion = DiceBCELoss(dice_weight=0.6) optimizer = torch.optim.AdamW(model.parameters(), lr=1e-4, weight_decay=1e-4) scheduler = torch.optim.lr_scheduler.CosineAnnealingLR(optimizer, T_max=50) for epoch in range(50): model.train() for img, mask in loader: img, mask = img.cuda(), mask.cuda() optimizer.zero_grad() loss = criterion(model(img), mask) loss.backward() optimizer.step() scheduler.step() print(f'epoch {epoch}, loss {loss.item():.4f}')逻辑说明:输入统一 resize 到 512×512,convert('L')保证单通道。掩码用>127二值化,兼容 0/255 和 0/1 两种导出。优化器用 AdamW,学习率 1e-4,配合余弦退火。batch_size 受显存限制,4 是 8GB 显存的稳妥值,显存够可以上 8。训练轮数 50 是起点,看验证集 Dice 曲线决定要不要加。
3.3 评估指标:别只看准确率
裂缝像素占比低,准确率(accuracy)会骗人——全预测背景也能到 95% 以上。必须看 Dice 和 IoU。Dice 对裂缝这种小目标更敏感,IoU 更严格。验证时按切片算 Dice 再平均,不要把所有像素混在一起算,否则大图会主导指标。
def dice_score(pred, target, eps=1e-6): pred = (torch.sigmoid(pred) > 0.5).float() inter = (pred * target).sum() return (2 * inter + eps) / (pred.sum() + target.sum() + eps)阈值 0.5 是默认,实际部署时可以调,后面第 5 章会讲怎么调。
4. 推理与后处理:让掩码从「能出」到「能用」
模型训完,推理出的原始掩码往往有毛刺、断线、小噪点。直接拿去统计裂缝面积,误差会很大。这一章讲推理脚本和三种后处理,把掩码修到能进分析流程。
4.1 批量推理脚本
import torch from PIL import Image import numpy as np import os def predict(model, img_path, size=512, threshold=0.5): model.eval() img = Image.open(img_path).convert('L').resize((size, size)) x = torch.from_numpy(np.array(img, dtype=np.float32) / 255.0)[None, None].cuda() with torch.no_grad(): logits = model(x) prob = torch.sigmoid(logits)[0, 0].cpu().numpy() mask = (prob > threshold).astype(np.uint8) * 255 return mask, prob model = UNet().cuda() model.load_state_dict(torch.load('best_unet.pth')) for name in os.listdir('./test_images'): mask, prob = predict(model, os.path.join('./test_images', name)) Image.fromarray(mask).save(f'./pred_masks/{name}')逻辑说明:推理时保持和训练一致的 resize 尺寸,否则尺度不匹配会掉点。threshold是二值化阈值,默认 0.5,后面会讲怎么调。保存概率图prob是为了后处理时能重新选阈值,不用重跑模型。
4.2 三种后处理:去噪、断线连接、阈值调优
去噪用连通域面积过滤,小于 30 像素的连通域直接删掉,这些多半是伪影。
from scipy import ndimage def remove_small(mask, min_area=30): labeled, n = ndimage.label(mask > 0) for i in range(1, n + 1): if (labeled == i).sum() < min_area: mask[labeled == i] = 0 return mask断线连接用形态学闭运算,3×3 或 5×5 核,把裂缝断开的地方接上。核别开太大,否则会把两条平行细缝粘成一条。
from scipy.ndimage import binary_closing def connect_cracks(mask, kernel_size=3): structure = np.ones((kernel_size, kernel_size), dtype=bool) return binary_closing(mask > 0, structure=structure).astype(np.uint8) * 255阈值调优:在验证集上扫 0.3 到 0.7,每 0.05 一档,算 Dice,选最高的。裂缝分割里 0.4 左右往往比 0.5 好,因为模型对裂缝边缘的置信度偏低,降阈值能把边缘捞回来,代价是噪点变多,配合面积过滤正好。
4.3 后处理顺序与参数联动
顺序是:先阈值二值化 → 闭运算连断线 → 面积过滤去噪。顺序反了会出问题:先面积过滤再闭运算,可能把刚连上的细缝又当成小连通域删掉。参数联动上,闭运算核越大,面积过滤的min_area也要相应调大,否则连出来的大块伪影删不掉。我一般固定闭运算核 3×3,min_area在 20–50 之间试。
5. 避坑与排查:裂缝分割里最容易翻车的五件事
这一章是我踩过的坑,按「现象 → 原因 → 解决」写,每条都能对上前面章节的操作。
现象一:训练 loss 一直不降,Dice 卡在 0.1 左右。原因多半是掩码没二值化,或者掩码和图像文件名没对上,模型在学噪声。解决:写个检查脚本,随机抽 10 对图,把掩码叠加到原图上可视化,确认裂缝位置对得上;再确认掩码像素值只有 0 和 1(或 0 和 255)。
现象二:验证集 Dice 很高,测试集一塌糊涂。原因是数据集按随机划分,相邻切片泄漏。解决:改成按深度区间分块划分,训练/验证/测试的切片在深度上不重叠。这个坑最隐蔽,指标虚高会让人误以为模型能用。
现象三:推理掩码全是噪点,裂缝反而没出来。原因通常是推理时的 resize 尺寸和训练不一致,或者归一化方式不同(训练除了 255,推理没除)。解决:把预处理封装成一个函数,训练和推理共用,杜绝两套代码。
现象四:细裂缝断成一段一段,统计面积偏小。原因是模型对细结构召回不足,加上阈值 0.5 偏高。解决:降阈值到 0.4,加闭运算连断线,再面积过滤。如果还断,考虑在损失里提高 Dice 权重,或在训练时加细裂缝的过采样。
现象五:换一批新岩心的 CT 数据,模型直接失效。原因是不同扫描设备的灰度分布、分辨率、伪影特征都不一样,模型过拟合了旧数据。解决:新数据上做少量标注(几十张),用低学习率 1e-5 微调,别从头训。如果新数据灰度分布差异大,重新做百分位拉伸再微调。
6. 把分割结果接进裂缝定量分析:一个可复现的统计脚本
模型出掩码只是中间产物,地质上真正要的是裂缝面积占比、裂缝密度、走向分布这些量。这一章给一个从掩码算裂缝面积占比和等效宽度的脚本,并讲怎么验证统计结果可信。
import numpy as np from PIL import Image from scipy import ndimage def crack_stats(mask_path, pixel_size_mm=0.05): mask = np.array(Image.open(mask_path).convert('L')) > 127 labeled, n = ndimage.label(mask) total_pixels = mask.size crack_pixels = mask.sum() area_ratio = crack_pixels / total_pixels # 每条裂缝的等效宽度 = 面积 / 骨架长度,这里用周长近似 widths = [] for i in range(1, n + 1): comp = labeled == i area = comp.sum() if area < 30: continue perimeter = np.logical_xor(comp, ndimage.binary_erosion(comp)).sum() if perimeter > 0: widths.append(2 * area / perimeter * pixel_size_mm) return { 'area_ratio': round(area_ratio, 4), 'crack_count': n, 'mean_width_mm': round(float(np.mean(widths)), 3) if widths else 0.0 } print(crack_stats('./pred_masks/slice_0001.png', pixel_size_mm=0.05))逻辑说明:pixel_size_mm是 CT 扫描的体素实际尺寸,必须从扫描参数里拿到,否则算出来的宽度没有物理意义。等效宽度用2 * 面积 / 周长近似,对细长裂缝够用,对分叉裂缝会偏大,所以统计时最好按连通域分别看。area_ratio是裂缝面积占比,直接对应地质上的裂缝孔隙度贡献。
验证统计可信度,我一般做两件事。一是抽 20 张测试切片,把模型掩码和人工标注掩码分别跑这个脚本,对比area_ratio的相对误差,控制在 10% 以内算可用。二是把掩码叠加回原图,肉眼扫一遍,重点看有没有把孔洞误判成裂缝、有没有漏掉大裂缝。这两步做完,统计结果才敢往报告里写。
最后说个习惯:这套流程里,模型结构可以换、损失可以调,但数据预处理和掩码规范一旦定下来就别轻易动。我见过太多项目,模型换了三四版,指标上上下下,最后发现是标注规范中途改了,前后数据不可比。把预处理和标注规范固化成脚本和文档,比追 SOTA 模型值钱得多。希望帮到你。
本文还有配套的精品资源,点击获取