简介这份资源面向计算机视觉与地质工程方向的本科生、研究生及课程设计开发者提供一套基于Python的CT岩芯与岩石裂缝语义分割完整方案可用于期末大作业、课程设计或相关课题的复现与二次开发。压缩包共15个文件约1.15MB包含3个Python脚本数据增强、均值计算等、6张jpg示例图像及对应标注图、若干zbak备份文件、README说明文档与gitignore配置覆盖从数据准备到模型训练的基础环节。资源以岩石、混凝土、CT岩芯等样本图像为对象演示了Pillow、OpenCV与TensorFlow/PyTorch等工具在像素级分类任务中的配合方式帮助读者理解裂隙分布识别与量化分析流程。目前已有72人学习适合希望快速获取可运行代码与配套数据、并在此基础上完成实验报告或项目原型的读者参考。1. 从一张 CT 切片说起岩心裂缝语义分割到底在解决什么一块直径 5 厘米的碳酸盐岩岩心送进工业 CT 扫描一圈出来的是 1000 到 2000 张灰度切片每张 1024×1024。人眼在屏幕上翻能看出裂缝、溶孔、基质的三相边界但翻到第 200 张就眼花。更现实的问题是一个区块几十米岩心按 0.1 毫米层厚扫下来切片数量轻松过万靠地质人员手工勾画裂缝轮廓一个人一周也标不完一根岩心。这就是基于 Python 的 CT 岩心与岩石裂缝语义分割系统要解决的核心痛点——把像素级裂缝识别从人工描边变成模型推理。语义分割在这里的定义很明确给每个像素打一个类别标签裂缝、基质、孔洞各归各的类输出一张和原图同尺寸的掩膜。它和实例分割的区别在于语义分割不区分“第几条裂缝”只区分“是不是裂缝”。对岩心裂缝统计开度、密度、连通性来说语义分割已经够用而且标注成本低得多。这套系统适合三类人做岩石物理实验、需要批量提取裂缝参数的研究生手里有 CT 数据、想快速验证分割效果的地质工程师以及想找一个真实工业数据集练语义分割的算法同学。热词里“语义分割数据集制作”和“语义分割模型”是绕不开的两件事。CT 岩心图像的难点不在模型结构而在数据本身灰度动态范围窄、裂缝和基质对比度低、扫描伪影多、裂缝宽度从 1 个像素到几十个像素不等。所以这套系统的落地路径不是先调模型而是先把数据管线搭稳再选一个对细小目标友好的分割网络最后把推理结果转成可量化的地质参数。下面按这个顺序拆开讲。2. 数据管线从 CT 切片到可训练掩膜的完整链路2.1 为什么 CT 岩心数据不能直接丢进 DataLoader工业 CT 输出的原始格式通常是 DICOM 或 TIFF 序列灰度值范围取决于扫描参数常见 16 位无符号整数实际有效区间可能只占 2000 到 8000 这一小段。直接归一化到 0-255 再送网络裂缝和基质的对比度会被压掉一大半。我一般会先做一次百分位裁剪把 1% 和 99% 分位数之外的值截断再线性拉伸到 0-255。这一步不做后面模型学到的就是一片灰蒙蒙。另一个坑是切片间的相关性。CT 岩心是连续扫描的相邻切片几乎一样如果按随机划分训练集和验证集验证集里会出现和训练集高度相似的切片指标虚高。正确做法是按深度分段划分比如前 70% 深度做训练中间 15% 做验证最后 15% 做测试。这样验证集才能反映模型在未见过的岩心段上的真实表现。标注方面裂缝掩膜是二值的但边缘像素往往处于过渡带。我建议标注时把裂缝边缘向外扩 1 个像素作为“忽略区”训练时不计损失。这样模型不会因为边缘像素的模糊性产生震荡。2.2 用 Python 做切片预处理与掩膜对齐下面这段代码做三件事读取 TIFF 序列、百分位拉伸、把标注掩膜和图像按文件名对齐。假设图像在images/下掩膜在masks/下文件名一一对应。import numpy as np import tifffile as tiff from pathlib import Path from sklearn.model_selection import train_test_split def percentile_stretch(img, low1, high99): 对 16 位 CT 图像做百分位裁剪并拉伸到 0-255 lo, hi np.percentile(img, (low, high)) img np.clip(img, lo, hi) img (img - lo) / (hi - lo 1e-6) * 255.0 return img.astype(np.uint8) def load_pairs(img_dir, mask_dir): 按文件名对齐图像和掩膜返回路径列表 img_paths sorted(Path(img_dir).glob(*.tif)) pairs [] for ip in img_paths: mp Path(mask_dir) / ip.name if mp.exists(): pairs.append((str(ip), str(mp))) return pairs # 按深度分段划分避免相邻切片泄漏 pairs load_pairs(images, masks) n len(pairs) train_pairs pairs[:int(n*0.7)] val_pairs pairs[int(n*0.7):int(n*0.85)] test_pairs pairs[int(n*0.85):] # 示例读取并拉伸一张切片 img tiff.imread(train_pairs[0][0]) img_stretched percentile_stretch(img) print(img_stretched.shape, img_stretched.dtype, img_stretched.min(), img_stretched.max())逻辑说明percentile_stretch里的low和high控制裁剪强度裂缝对比度特别低时可以试 0.5 和 99.5但不要低于 0.1否则噪声会被拉起来。load_pairs用文件名对齐是最稳的方式不要依赖排序索引因为标注软件导出时可能漏掉几张。划分比例 70/15/15 是经验值数据量少于 500 张时验证集可以再小一点但测试集不能省。参数方面tifffile读 16 位 TIFF 不会丢精度如果用 OpenCV 的imread默认会转成 8 位这是常见翻车点。另外拉伸后的图像建议保存成 PNG 而不是 JPGJPG 的块效应会在裂缝边缘产生伪影影响分割精度。2.3 数据增强针对裂缝形态的几何变换CT 岩心裂缝有方向性但方向是地质成因决定的不能随便旋转。我一般只做水平翻转、垂直翻转和 90 度旋转不做任意角度旋转因为任意角度会引入插值模糊细裂缝直接消失。亮度对比度扰动可以做但幅度要小±10% 以内否则会破坏灰度与岩性的对应关系。import albumentations as A train_transform A.Compose([ A.HorizontalFlip(p0.5), A.VerticalFlip(p0.5), A.RandomRotate90(p0.5), A.RandomBrightnessContrast(brightness_limit0.1, contrast_limit0.1, p0.3), ]) val_transform A.Compose([]) # 验证集不做增强RandomRotate90只做 90 度整数倍旋转不插值。RandomBrightnessContrast的brightness_limit和contrast_limit设 0.1 是保守值如果数据本身灰度差异大可以到 0.2但再大就会让模型把亮度变化误判成裂缝。验证集和测试集一律不做增强保证指标可比。3. 模型选型为什么 U-Net 变体仍然是 CT 岩心分割的稳妥起点3.1 裂缝像素占比极低时损失函数比网络结构更关键一根岩心的 CT 切片里裂缝像素占比通常不到 3%有的层位甚至低于 0.5%。这种极端类别不平衡下交叉熵损失会被基质像素主导模型学会全预测背景就能拿到 97% 以上的准确率但裂缝一个都分不出来。我试过直接上交叉熵训练 50 轮后 IoU 停在 0.02血泪经验。解决办法是组合损失Dice Loss 加 Focal Loss。Dice 直接优化预测和真值的重叠度对类别不平衡不敏感Focal Loss 降低易分类样本的权重让模型聚焦在难分的裂缝边缘。权重上我一般 Dice 占 0.6Focal 占 0.4Focal 的 gamma 设 2.0alpha 设 0.75偏向裂缝类。import torch import torch.nn as nn import torch.nn.functional as F class DiceLoss(nn.Module): def __init__(self, smooth1.0): super().__init__() self.smooth smooth def forward(self, logits, targets): probs torch.sigmoid(logits) probs probs.view(-1) targets targets.view(-1) intersection (probs * targets).sum() dice (2. * intersection self.smooth) / (probs.sum() targets.sum() self.smooth) return 1 - dice class FocalLoss(nn.Module): def __init__(self, alpha0.75, gamma2.0): super().__init__() self.alpha alpha self.gamma gamma def forward(self, logits, targets): bce F.binary_cross_entropy_with_logits(logits, targets, reductionnone) pt torch.exp(-bce) focal self.alpha * (1 - pt) ** self.gamma * bce return focal.mean() class CombinedLoss(nn.Module): def __init__(self, dice_w0.6, focal_w0.4): super().__init__() self.dice DiceLoss() self.focal FocalLoss() self.dice_w dice_w self.focal_w focal_w def forward(self, logits, targets): return self.dice_w * self.dice(logits, targets) self.focal_w * self.focal(logits, targets)DiceLoss里的smooth防止分母为零设 1.0 是常规做法。FocalLoss的alpha偏向裂缝类因为裂缝是正类且样本少。CombinedLoss的两个权重可以根据验证集 IoU 微调如果裂缝召回率低就加大 Dice 权重如果误检多就加大 Focal 权重。3.2 U-Net 的编码器换与不换看数据量标准 U-Net 用 VGG 或 ResNet 做编码器参数量在 2000 万到 4000 万。CT 岩心数据集通常只有几百到几千张标注切片这个量级下从零训练大编码器必然过拟合。我的做法是编码器用 ImageNet 预训练的 ResNet34解码器从零初始化训练时编码器学习率设小一点1e-4解码器设大一点1e-3。如果数据少于 300 张直接把编码器冻结只训解码器效果反而更稳。另一个选择是轻量级的 U-Net比如把通道数减半参数量降到 200 万左右。这种适合做快速验证但裂缝细节会丢。我一般先用轻量版跑通全流程确认数据管线没问题再换标准版冲指标。import segmentation_models_pytorch as smp model smp.Unet( encoder_nameresnet34, encoder_weightsimagenet, in_channels1, # CT 灰度图单通道 classes1, # 二分类裂缝/背景 activationNone, # 输出 logits损失函数里做 sigmoid ) # 编码器和解码器分组学习率 encoder_params list(model.encoder.parameters()) decoder_params list(model.decoder.parameters()) list(model.segmentation_head.parameters()) optimizer torch.optim.Adam([ {params: encoder_params, lr: 1e-4}, {params: decoder_params, lr: 1e-3}, ])in_channels1是因为 CT 灰度图是单通道不要为了套用 RGB 预训练权重强行复制成三通道那样第一层卷积的权重会失效。activationNone让模型输出原始 logits配合binary_cross_entropy_with_logits数值更稳定。分组学习率是微调预训练编码器的常规操作如果编码器冻结就把encoder_params的lr设 0 或者直接从优化器里去掉。3.3 训练循环里必须监控的三个指标准确率在裂缝分割里没有意义我只看三个裂缝类的 IoU、裂缝类的召回率、以及裂缝边缘的 F1。IoU 反映整体重叠度召回率反映漏检边缘 F1 反映边界质量。验证时每 5 个 epoch 算一次如果 IoU 连续 10 个 epoch 不涨就降学习率或者停。def compute_metrics(logits, targets, threshold0.5): probs torch.sigmoid(logits) preds (probs threshold).float() targets targets.float() intersection (preds * targets).sum() union preds.sum() targets.sum() - intersection iou (intersection 1e-6) / (union 1e-6) recall (intersection 1e-6) / (targets.sum() 1e-6) return iou.item(), recall.item()threshold默认 0.5但裂缝分割里可以调到 0.4 提高召回代价是误检增加。具体取值看下游任务如果做裂缝统计宁可误检不可漏检阈值调低如果做三维重建误检会生成假裂缝阈值调高。这个参数没有标准答案要在验证集上扫一遍。4. 推理与后处理把概率图变成可量化的裂缝参数4.1 滑窗推理与重叠拼接CT 切片尺寸 1024×1024直接送网络显存吃不消而且 U-Net 的下采样会丢失细裂缝。常规做法是切 256×256 的滑窗步长 128重叠 50%。重叠区域取平均避免拼接缝。步长太小推理慢太大拼接缝明显128 是折中。def sliding_window_inference(model, image, window256, stride128, devicecuda): model.eval() h, w image.shape prob_map np.zeros((h, w), dtypenp.float32) count_map np.zeros((h, w), dtypenp.float32) with torch.no_grad(): for y in range(0, h - window 1, stride): for x in range(0, w - window 1, stride): patch image[y:ywindow, x:xwindow] tensor torch.from_numpy(patch).float().unsqueeze(0).unsqueeze(0).to(device) logits model(tensor) prob torch.sigmoid(logits).squeeze().cpu().numpy() prob_map[y:ywindow, x:xwindow] prob count_map[y:ywindow, x:xwindow] 1 prob_map / np.maximum(count_map, 1) return prob_mapwindow和stride要根据裂缝宽度调。裂缝最宽 20 像素时256 窗口足够如果裂缝贯穿整张切片窗口要加大到 512否则一条裂缝被切成几段后处理时连不起来。count_map保证重叠区域被平均而不是累加。4.2 后处理去小连通域与骨架化模型输出的概率图二值化后会有一些孤立的假阳性小斑点。用连通域面积过滤面积小于 50 像素的直接去掉。然后对裂缝做骨架化提取中心线用来算开度和长度。from skimage.measure import label, regionprops from skimage.morphology import skeletonize def postprocess(prob_map, threshold0.5, min_area50): binary (prob_map threshold).astype(np.uint8) labeled label(binary, connectivity2) cleaned np.zeros_like(binary) for region in regionprops(labeled): if region.area min_area: cleaned[labeled region.label] 1 skeleton skeletonize(cleaned.astype(bool)) return cleaned, skeletonmin_area设 50 是经验值切片分辨率高时可以调到 100分辨率低时降到 20。skeletonize输出的骨架是单像素宽用来统计裂缝条数和分支点。开度可以用距离变换算对二值图做distance_transform_edt骨架位置的距离值乘 2 就是局部开度。4.3 从掩膜到地质参数开度、密度、连通性拿到二值掩膜和骨架后可以算三个常用参数。开度骨架点上的距离变换值乘 2取平均。密度裂缝总长度除以切片面积。连通性骨架的连通域数量越少说明裂缝越连通。这些参数可以直接输出成 CSV配合深度列就能画裂缝密度曲线。from scipy.ndimage import distance_transform_edt import pandas as pd def extract_params(binary, skeleton, depth): dist distance_transform_edt(binary) aperture dist[skeleton].mean() * 2 if skeleton.sum() 0 else 0 length skeleton.sum() area binary.shape[0] * binary.shape[1] density length / area n_components label(skeleton, connectivity2).max() return {depth: depth, aperture: aperture, density: density, components: n_components} # 批量处理 records [] for i, (img_path, _) in enumerate(test_pairs): img percentile_stretch(tiff.imread(img_path)) prob sliding_window_inference(model, img) binary, skeleton postprocess(prob) records.append(extract_params(binary, skeleton, depthi*0.1)) df pd.DataFrame(records) df.to_csv(fracture_params.csv, indexFalse)depth按切片顺序乘层厚得到层厚从 CT 扫描参数里读。aperture单位是像素要乘像素物理尺寸才是毫米。components是骨架连通域数量如果一条裂缝被断开成多段这个值会偏大说明后处理的min_area或推理阈值需要调。5. 避坑与排查CT 岩心分割里最容易翻车的五件事5.1 验证集 IoU 很高测试集一塌糊涂现象训练时验证集 IoU 到 0.85换测试集掉到 0.4。原因按随机划分导致相邻切片泄漏验证集和训练集几乎一样。解决改成按深度分段划分训练/验证/测试之间留至少 50 张切片的间隔。如果数据本身来自多根岩心按岩心编号划分更稳。5.2 裂缝边缘出现锯齿状毛刺现象分割结果里裂缝边界像狗牙不光滑。原因下采样倍数太高U-Net 的 4 次下采样把细裂缝丢了上采样时用最近邻插值。解决把上采样改成双线性插值或者在解码器里加 skip connection 的注意力模块。另一个办法是推理时用 50% 重叠滑窗重叠平均能平滑边界。5.3 模型把扫描伪影当成裂缝现象环形伪影或束硬化伪影被分割成裂缝。原因训练集里伪影区域没有标注为忽略区模型学到了伪影特征。解决在标注阶段把伪影区域标成 255忽略区损失函数里忽略这些像素。如果伪影已经影响训练先做一次 CT 伪影校正用scikit-image的denoise_tv_chambolle做全变分去噪。5.4 显存溢出batch size 只能设 1现象256×256 窗口batch size 设 4 就 OOM。原因U-Net 编码器通道数太大或者用了 1024×1024 全图训练。解决用混合精度训练torch.cuda.amp能省 30% 到 40% 显存。如果还不够把编码器换成 ResNet18 或 MobileNetV2参数量减半。梯度累积也能模拟大 batch累积 4 步等效 batch size 4。5.5 推理速度太慢一根岩心跑一晚上现象1000 张切片每张滑窗推理要 2 秒总共 30 多分钟加上后处理超过 1 小时。原因滑窗步长太小重叠太多。解决步长从 128 调到 192重叠从 50% 降到 25%速度提升近一倍IoU 只掉 0.01 到 0.02。另外把模型转成 TorchScript 或 ONNX推理速度能再提 20% 到 30%。如果只是做参数统计可以隔张推理裂缝在深度方向连续隔张不会漏掉整条裂缝。6. 把系统跑成可复现的流水线一个配置文件的写法整套系统跑通后最怕的是换一台机器就复现不出来。我的习惯是把所有路径、超参、阈值写进一个 YAML 配置文件训练、推理、后处理都读同一个配置。这样别人拿到源码和数据集改一下data_root就能跑。# config.yaml data: image_dir: data/images mask_dir: data/masks split_ratio: [0.7, 0.15, 0.15] percentile: [1, 99] train: window_size: 256 batch_size: 4 epochs: 100 encoder_lr: 1e-4 decoder_lr: 1e-3 loss: dice_weight: 0.6 focal_weight: 0.4 focal_alpha: 0.75 focal_gamma: 2.0 inference: window_size: 256 stride: 192 threshold: 0.45 min_area: 50 output: param_csv: output/fracture_params.csv mask_dir: output/masks配置里threshold设 0.45 而不是 0.5是因为验证集上扫出来 0.45 的召回更高IoU 只低 0.005。stride设 192 对应 25% 重叠速度和质量平衡。min_area设 50 是 1024×1024 切片下的经验值如果切片尺寸变了要按面积比例调。加载配置的代码很简单但要注意路径拼接用pathlib不要用字符串加号Windows 和 Linux 的斜杠不一样。import yaml from pathlib import Path def load_config(pathconfig.yaml): with open(path, r, encodingutf-8) as f: cfg yaml.safe_load(f) root Path(cfg[data][image_dir]).parent.parent cfg[data][image_dir] str(root / cfg[data][image_dir]) cfg[data][mask_dir] str(root / cfg[data][mask_dir]) return cfgroot取image_dir的上两级是为了让配置里的路径相对于项目根目录而不是相对于当前工作目录。这样在train.py和infer.py里都能用同一份配置不用改路径。最后说一个我自己的习惯每次跑完实验把配置文件、训练日志、验证集指标、测试集指标打包存到一个以日期命名的文件夹里。CT 岩心数据标注成本高一次实验可能隔几周才回头对比没有存档就只能重跑。这个习惯帮我省过至少两次返工。希望帮到你。本文还有配套的精品资源点击获取