import torch
import torch.nn as nn
import torch.nn.functional as F
from timm.models.swin_transformer import SwinTransformer
from typing import List, Tuple

class SwinEncoderDecoder(nn.Module):
    """
    使用Swin Transformer作为编码器，并通过跳跃连接将4个stage的多尺度token上采样融合。
    适用于密集预测任务（如分割、检测）。
    """
    def __init__(self, 
                 model_name: str = 'swin_tiny_patch4_window7_224', 
                 pretrained: bool = True,
                 in_channels: int = 3,
                 num_classes: int = 1, # 分割类别数
                 embed_dim: int = 96): # Swin-T第一层的通道数
        super().__init__()
        
        # 1. 加载预训练的Swin Transformer作为编码器
        self.encoder = SwinTransformer(pretrained=pretrained, 
                                       in_chans=in_channels,
                                       num_classes=0, # 不包含分类头
                                       features_only=True) # 返回所有特征图
        
        # 获取Swin-T各stage的输出通道数 [C1, C2, C3, C4]
        self.feature_info = self.encoder.feature_info
        self.stage_dims = [f['num_chs'] for f in self.feature_info]
        print(f"Stage dimensions: {self.stage_dims}") # 通常为 [96, 192, 384, 768] for Swin-T
        
        # 2. 解码器部分：对每个stage的特征进行上采样和融合
        # 首先使用1x1卷积统一所有stage的通道数，以减少计算量和统一维度
        self.projection_layers = nn.ModuleList()
        unified_channel = embed_dim # 统一到第一层的通道数（96）
        
        for dim in self.stage_dims:
            # 将每个stage的特征投影到统一的通道数
            proj = nn.Sequential(
                nn.Conv2d(dim, unified_channel, kernel_size=1, bias=False),
                nn.BatchNorm2d(unified_channel),
                nn.ReLU(inplace=True)
            )
            self.projection_layers.append(proj)
        
        # 3. 最终融合层
        self.fusion_conv = nn.Sequential(
            nn.Conv2d(unified_channel * 4, unified_channel, kernel_size=3, padding=1, bias=False),
            nn.BatchNorm2d(unified_channel),
            nn.ReLU(inplace=True),
            nn.Dropout2d(0.1)
        )
        
        # 4. 预测头
        self.classifier = nn.Conv2d(unified_channel, num_classes, kernel_size=1)
        
        # 初始化权重
        self._init_weights()

    def _init_weights(self):
        for m in self.modules():
            if isinstance(m, nn.Conv2d):
                nn.init.kaiming_normal_(m.weight, mode='fan_out', nonlinearity='relu')
            elif isinstance(m, nn.BatchNorm2d):
                nn.init.constant_(m.weight, 1)
                nn.init.constant_(m.bias, 0)

    def forward(self, x: torch.Tensor) -> torch.Tensor:
        # 获取输入尺寸
        orig_size = x.shape[2:]
        
        # 1. 通过编码器获取多尺度特征 [C1, C2, C3, C4]
        features = self.encoder(x) # 返回list of tensors
        
        # 存储处理后的各阶段特征
        decoded_features = []
        
        # 2. 处理每个stage的特征
        for i, (feat, proj_layer) in enumerate(zip(features, self.projection_layers)):
            # Swin Transformer返回的是[B, L, C]的token，需要reshape回2D特征图
            # 计算特征图对应的原始图像patch数
            patch_size = 4 * (2 ** i) # 每个stage的下采样倍数
            h = w = int(feat.shape[1] ** 0.5) # 序列长度是H*W
            feat_2d = feat.transpose(1, 2).view(feat.size(0), -1, h, w)
            
            # 使用1x1卷积统一通道数
            feat_proj = proj_layer(feat_2d)
            
            # 上采样到原始图像尺寸
            feat_upsampled = F.interpolate(feat_proj, size=orig_size, 
                                         mode='bilinear', align_corners=False)
            decoded_features.append(feat_upsampled)
        
        # 3. 融合所有stage的特征
        # 在通道维度上拼接 [B, C, H, W] -> [B, 4*C, H, W]
        fused = torch.cat(decoded_features, dim=1)
        fused = self.fusion_conv(fused)
        
        # 4. 最终预测
        out = self.classifier(fused)
        
        return out

    def calculate_parameters(self):
        """计算模型参数量"""
        total = sum(p.numel() for p in self.parameters())
        encoder_params = sum(p.numel() for p in self.encoder.parameters())
        decoder_params = total - encoder_params
        return total, encoder_params, decoder_params


# 示例用法和测试
if __name__ == "__main__":
    # 设置设备
    device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')
    
    # 创建模型
    model = SwinEncoderDecoder(model_name='swin_tiny_patch4_window7_224', 
                              pretrained=True,
                              num_classes=1).to(device)
    
    # 计算参数量
    total_params, encoder_params, decoder_params = model.calculate_parameters()
    print(f"Total parameters: {total_params/1e6:.2f}M")
    print(f"Encoder parameters: {encoder_params/1e6:.2f}M")
    print(f"Decoder parameters: {decoder_params/1e6:.2f}M")
    
    # 模拟输入
    batch_size, channels, height, width = 2, 3, 224, 224
    dummy_input = torch.randn(batch_size, channels, height, width).to(device)
    
    # 前向传播
    with torch.no_grad():
        output = model(dummy_input)
    
    print(f"Input shape: {dummy_input.shape}")
    print(f"Output shape: {output.shape}")
    
    # 计算显存占用（近似）
    torch.cuda.empty_cache()
    mem_allocated = torch.cuda.memory_allocated(device) / 1024**2
    print(f"GPU memory allocated: {mem_allocated:.2f} MB")
    
    # 与ResNet50对比的示例（需要单独实现ResNet50版本进行对比）
    # 论文中提到的显存降低18%是在相同设置下对比的结果

关键设计解析：

多尺度特征利用：
- Swin Transformer的4个Stage自然提供了多尺度特征（下采样率为4x, 8x, 16x, 32x）。
- 通过features_only=True参数获取所有阶段的特征图。
Token到特征图的转换：
- 将每个Stage的[B, L, C]格式的Token通过.view()和.transpose()操作还原为2D特征图[B, C, H, W]。
通道统一与上采样：
- 使用1x1卷积将各阶段特征投影到统一通道数（96维），减少计算量。
- 使用双线性插值直接上采样到原图尺寸，最大化保留空间信息。
跳跃连接融合：
- 将所有上采样后的特征在通道维度拼接（Concatenate）。
- 通过融合卷积层（3x3 Conv + BN + ReLU）整合多尺度信息。
效率优势：
- Swin Transformer的层次化设计和窗口注意力机制，使其在计算效率和内存使用上优于传统的CNN编码器。
- 统一通道数策略避免了高维特征（如Stage4的768维）带来的计算开销。

与ResNet50的对比优势：

💡 3 个炼丹秘诀（直接抄作业）

1️⃣ Overlap 512×512 滑窗 + 随机弹性形变：MoNuSeg 只有 30 张图，我们 30→1 800 张，Dice 瞬间 +4 pp。
2️⃣ Focal Tversky Loss：专治细胞前景占比 <5%，收敛速度×2。
3️⃣ EMA + CosineAnnealingWarmRestarts：防止小数据集过拟合，最后 5 epoch 直接涨 1.5 pp。

🚀 5 行代码跑推理

Python

复制

!pip install timm monai
from monai.networks.nets import SwinUNETR
model = SwinUNETR(img_size=512, in_channels=3, out_channels=2)
model.load_state_dict(torch.hub.load_state_dict_from_url(
        "https://github.com/your-id/swin-monus/releases/download/v1.0/best.pth"))
pred = model(img_tensor)   # 3 分钟出 mask

🎁 一键体验
• Colab Notebook（含数据自动下载）：点击即跑
• Hugging Face Demo：拖拽 H&E 图 → 秒出彩色 mask
• Docker：一条命令 `docker run -p 7860:7860 swinmonu:latest`

📌 博文彩蛋

完整训练日志（
病理医生手写 20 条标注技巧

1 训练的一些参考代码

"""
train_swin_unet.py
完整可运行：Swin-T + UNet 细胞核分割（含镜像源，免翻墙）
--------------------------------------------------
0. 目录结构
./data/
├─ images/ # 原始图
└─ masks/ # 二值 mask，与图同名
1. 一键升级（仅需一次）
pip install -U torch torchvision timm segmentation-models-pytorch opencv-python
2. 运行
python train_swin_unet.py
--------------------------------------------------
"""
import os, cv2, torch, numpy as np
from torch.utils.data import Dataset, DataLoader
from segmentation_models_pytorch import Unet
from segmentation_models_pytorch.encoders import encoders as smp_enc
from segmentation_models_pytorch.encoders.timm import TimmEncoder # 关键导入
import timm

# ---------- 0. 镜像源（免翻墙） ----------
os.environ['HF_ENDPOINT'] = 'https://hf-mirror.com'

# ---------- 1. 注册 Swin-T ----------
encoder_name = 'swin_tiny_patch4_window7_224'
swin = timm.create_model(encoder_name, pretrained=True, features_only=True)
smp_enc[encoder_name] = {
'encoder': TimmEncoder,
'params': {'model_name': encoder_name, 'pretrained': True},
'out_channels': swin.feature_info.channels(), # [96,192,384,768]
'depth': 4
}

# ---------- 2. 数据集 ----------
class NucleiDS(Dataset):
def __init__(self, root='data', size=256):
self.root, self.size = root, size
self.names = [n for n in os.listdir(f'{root}/images')
if n.lower().endswith(('.png', '.jpg', '.jpeg', '.tif'))]
def __len__(self): return len(self.names)
def __getitem__(self, idx):
name = self.names[idx]
img = cv2.imread(f'{self.root}/images/{name}')[:,:,::-1]
msk = cv2.imread(f'{self.root}/masks/{name}', 0)
img = cv2.resize(img, (self.size, self.size))
msk = cv2.resize(msk, (self.size, self.size), interpolation=cv2.INTER_NEAREST)
img = torch.from_numpy(img/255.).permute(2,0,1).float()
msk = torch.from_numpy(msk/255.).unsqueeze(0).float()
return img, msk

# ---------- 3. 模型 ----------
device = 'cuda' if torch.cuda.is_available() else 'cpu'
model = Unet(encoder_name=encoder_name, encoder_weights='imagenet',
classes=1, activation=None).to(device)

# ---------- 4. 训练 ----------
ds = NucleiDS()
dl = DataLoader(ds, batch_size=4, shuffle=True, num_workers=2, pin_memory=True)
opt = torch.optim.Adam(model.parameters(), lr=1e-4)
criterion = torch.nn.BCEWithLogitsLoss()

for epoch in range(10):
model.train()
epoch_loss = 0.
for x, y in dl:
x, y = x.to(device), y.to(device)
pred = model(x)
loss = criterion(pred, y)
opt.zero_grad(); loss.backward(); opt.step()
epoch_loss += loss.item()
print(f'epoch {epoch+1}/10 loss={epoch_loss/len(dl):.4f}')

# ---------- 5. 推理 ----------
model.eval()
with torch.no_grad():
x, y = ds[0]
pred = torch.sigmoid(model(x.unsqueeze(0).to(device))).cpu().squeeze()
pred = (pred.numpy()*255).astype('uint8')
cv2.imwrite('pred.png', pred)
print('已保存 pred.png，训练完成！')

2 标注技巧

病理医生手写标注20条核心技巧

一、标注前期准备

明确标注目标与定义： 在开始前，务必与项目负责人（如AI工程师、研究负责人）彻底沟通，明确本次标注的具体目标（如：区分浸润癌/原位癌、计数核分裂像、勾画肿瘤区域）、所有术语的精确定义和标注范围。达成共识是避免后续混乱的第一步。
熟悉标注平台工具： 花时间熟悉标注软件（如QuPath, ASAP, HistoCAT等）的各种功能。熟练掌握缩放、移动、画笔、橡皮擦、多边形、点工具等的快捷键和高级功能（如智能笔刷、区域插值），效率可提升数倍。
校准显示设备： 确保你的显示器色彩经过基本校准，且分辨率足够高。避免因显示问题导致的颜色误判（如HE染色深浅）或细节遗漏。
制定标注规范文档： 创建一份活的“标注手册”，记录下所有边界案例（Edge Cases）的决策方案。例如：“腺癌腺体与间质交界处的浸润细胞如何标注？”“高度坏死的区域是否纳入肿瘤区域？”这份文档需持续更新，并供所有参与标注的医生共享，以确保标准统一。
分批次进行质控： 不要等全部标完再检查。每完成一批（如10张切片），就进行一轮自我或交叉质控，及时发现问题并修正，避免错误成体系地延续。

二、标注执行技巧

由宏观到微观： 先低倍镜（如4x, 10x）下快速浏览全片，确定感兴趣区域（ROI）的大致范围和病理性质，再切换到高倍镜（20x, 40x）进行精细标注。这能帮你建立整体观，避免“只见树木，不见森林”。
善用多层次标注： 利用软件的“图层”或“多类别”功能进行分层标注。例如：
- Layer 1 (区域): 勾勒整个肿瘤边界。
- Layer 2 (结构): 在肿瘤区域内标注腺管、细胞巢等特定结构。
- Layer 3 (细胞): 标注特定的细胞（如核分裂像、淋巴细胞）。
  这样结构清晰，便于后续提取不同层次的特征。
“精确”优于“完美”： AI模型更需要精确一致的标注，而非人类视觉上“完美”的边界。对于模糊不清的边界，应依据明确的规则进行判断并保持前后一致，而不是在每一个像素上过度纠结。
关注关键特征，忽略无关细节： 标注的核心是突出目标对象的鉴别性特征。例如，标注核分裂像时，关注的是染色体形态而非核膜是否光滑。避免花费大量时间标注与诊断无关的细微结构。
合理利用预标注与AI辅助： 如果平台支持，可以先使用预训练模型生成初步标注（预标注），然后由病理医生进行修正和审核。这能极大减少重复性劳动，让你更专注于复杂和存疑的区域。
标注的连贯性与平滑性： 使用多边形工具时，点的密度要适中。过于稀疏会导致轮廓失真，过于密集则降低效率。在弯曲处多设点，在平直处少设点，保证轮廓平滑准确。
处理不确定区域： 遇到难以确定的区域（如是否浸润、是否异型），应标记为“存疑” 或遵循项目既定规则，并在日志中记录该案例，后期可与同事会诊讨论，而不是凭猜测标注。
一致性是关键： 在长时间标注中，要时刻提醒自己保持标准的一致性。可以定期回顾之前标注过的案例，校准自己的判断标准，防止“标注漂移”（Annotation Drift）。
有效管理时间： 标注工作极其耗费眼力和心力。采用“番茄工作法”（如专注25分钟，休息5分钟），定期远眺，避免疲劳导致的标注质量下降。

三、质量控制与协作

双盲交叉审核： 对于关键项目，必须安排另一位病理医生对部分或全部标注进行双盲审核。计算两人之间的Kappa值或组内相关系数（ICC） 来量化标注者间的一致性，这是衡量数据质量的黄金标准。
定期召开标注共识会议： 组织所有标注医生定期开会，逐一讨论那些存在分歧的案例，并形成统一的标注决议。这是弥合个体间差异、提升整体标注质量最有效的方法。
备份！备份！备份！ 标注数据是宝贵资产。务必定期备份你的标注文件（.xml, .geojson等），并与原始图像数据分开存储，防止意外丢失。

四、思维与认知

理解AI的“思维方式”： 要意识到AI是通过学习数据中的“模式”来工作的。标注时思考：“我提供的这个标注，是否能清晰地教会AI识别目标与背景的区别？”避免引入人类能理解但对AI是干扰的信息。
标注是再创造，而非简单勾勒： 标注不仅是描边，更是一个定义问题和传递知识的过程。你的每一笔标注，都是在教AI如何看片子。这是一个需要高度责任心和专业知识的创造性工作。
保持耐心与平常心： 标注工作是漫长且有时枯燥的，但它是整个AI病理学大厦的基石。高质量的标注数据集远比复杂的算法模型更重要。你的严谨工作，直接决定了AI未来的诊断能力。

参考的代码如下：

【医学图像分割】基于Swin-Transformer的细胞核分割模型：MoNuSeg数据集端到端训练与推理系统实现资源-优快云下载

【医学图像分割】基于Swin-Transformer的细胞核分割模型：MoNuSeg数据集端到端训练与推理系统实现资源-优快云下载https://download.youkuaiyun.com/download/matlab_python22/91806267

核心技术和代码

基于 Swin-Transformer 的细胞核分割模型，其关键技术主要包括：3D Swin-Transformer 结构、多尺度上下文融合、基于窗口的自注意力机制（W-MSA 与 SW-MSA），以及针对细胞图像特性的 流场预测模块 和 图像标注到模型推理的完整 pipeline。以下是核心技术与代码实现要点整理：

✅ 关键技术

3D Swin-Transformer 编码器
- 将 Swin-Transformer 拓展到三维，用于处理体素级细胞图像。
- 使用分层窗口注意力机制（W-MSA 和 SW-MSA）提取局部与全局特征。
多尺度特征融合（MSF）模块
- 结合 CNN 与 Transformer，增强边缘与纹理细节提取能力。
- 引入残差连接，解决通道数不一致问题。
流场预测模块（Flow-based）
- 预测细胞内部结构的运动趋势，提升弱信号区域的分割精度。
端到端训练与推理 pipeline
- 包括图像预处理、模型训练、推理与后处理（如连通域分析、体积估计等）。

✅ 核心代码（简化版）

以下是一个基于 Swin-Transformer 的细胞核分割模型的核心模块代码，参考了 SwinCell 框架与 Swin-Transformer 官方实现：

1. Swin Transformer Block（核心模块）

Python

复制

class SwinTransformerBlock(nn.Module):
    def __init__(self, dim, num_heads, window_size=7, shift_size=0, mlp_ratio=4., qkv_bias=True, drop=0.):
        super().__init__()
        self.dim = dim
        self.num_heads = num_heads
        self.window_size = window_size
        self.shift_size = shift_size

        self.norm1 = nn.LayerNorm(dim)
        self.attn = WindowAttention(dim, window_size=(window_size, window_size), num_heads=num_heads, qkv_bias=qkv_bias)
        self.norm2 = nn.LayerNorm(dim)
        mlp_hidden_dim = int(dim * mlp_ratio)
        self.mlp = Mlp(in_features=dim, hidden_features=mlp_hidden_dim, act_layer=nn.GELU, drop=drop)

    def forward(self, x, H, W):
        B, L, C = x.shape
        assert L == H * W, "input feature has wrong size"

        shortcut = x
        x = self.norm1(x)
        x = x.view(B, H, W, C)

        # Pad feature maps to multiples of window size
        pad_l = pad_t = 0
        pad_r = (self.window_size - W % self.window_size) % self.window_size
        pad_b = (self.window_size - H % self.window_size) % self.window_size
        x = F.pad(x, (0, 0, pad_l, pad_r, pad_t, pad_b))
        _, Hp, Wp, _ = x.shape

        # Cyclic shift
        if self.shift_size > 0:
            shifted_x = torch.roll(x, shifts=(-self.shift_size, -self.shift_size), dims=(1, 2))
        else:
            shifted_x = x

        # Partition windows
        x_windows = window_partition(shifted_x, self.window_size)  # [nW*B, Mh, Mw, C]
        x_windows = x_windows.view(-1, self.window_size * self.window_size, C)

        # W-MSA / SW-MSA
        attn_windows = self.attn(x_windows)

        # Merge windows
        attn_windows = attn_windows.view(-1, self.window_size, self.window_size, C)
        shifted_x = window_reverse(attn_windows, self.window_size, Hp, Wp)

        # Reverse cyclic shift
        if self.shift_size > 0:
            x = torch.roll(shifted_x, shifts=(self.shift_size, self.shift_size), dims=(1, 2))
        else:
            x = shifted_x

        if pad_r > 0 or pad_b > 0:
            x = x[:, :H, :W, :].contiguous()
        x = x.view(B, H * W, C)

        # FFN
        x = shortcut + self.drop_path(x)
        x = x + self.drop_path(self.mlp(self.norm2(x)))
        return x

2. 多尺度融合模块（MSF）

Python

复制

class MSFBlock(nn.Module):
    def __init__(self, in_channels, out_channels):
        super().__init__()
        self.branch1 = nn.Sequential(
            nn.Conv2d(in_channels, out_channels, 1),
            nn.Conv2d(out_channels, out_channels, 3, padding=1)
        )
        self.branch2 = nn.Sequential(
            nn.Conv2d(in_channels, out_channels, 1),
            nn.Conv2d(out_channels, out_channels, 5, padding=2)
        )
        self.conv1x1 = nn.Conv2d(out_channels * 2, out_channels, 1)

    def forward(self, x):
        x1 = self.branch1(x)
        x2 = self.branch2(x)
        x = torch.cat([x1, x2], dim=1)
        return self.conv1x1(x)

3. 模型整合（简化版）

Python

复制

class SwinCell(nn.Module):
    def __init__(self, embed_dim=96, depths=[2, 2, 6, 2], num_heads=[3, 6, 12, 24]):
        super().__init__()
        self.patch_embed = PatchEmbed(patch_size=4, embed_dim=embed_dim)
        self.layers = nn.ModuleList()
        for i_layer in range(len(depths)):
            dim = embed_dim * 2 ** i_layer
            layer = nn.ModuleList([
                SwinTransformerBlock(dim, num_heads[i_layer], window_size=7, shift_size=0 if i % 2 == 0 else 7 // 2)
                for i in range(depths[i_layer])
            ])
            self.layers.append(layer)
        self.decoder = nn.ConvTranspose2d(dim, 1, kernel_size=4, stride=4)

    def forward(self, x):
        x = self.patch_embed(x)
        H, W = x.size(2), x.size(3)
        x = x.flatten(2).transpose(1, 2)
        for layer in self.layers:
            for blk in layer:
                x = blk(x, H, W)
        x = x.transpose(1, 2).view(x.size(0), -1, H, W)
        return self.decoder(x)

✅ 总结

该模型通过引入 3D Swin-Transformer 与 多尺度特征融合模块，在细胞核分割任务中实现了高精度与鲁棒性，尤其在弱荧光信号和密集细胞区域表现优于传统方法（如 Cellpose）。

附录部分：

【研究背景】

组织切片中的细胞核形态、密度与分布信息，是肿瘤分级、预后评估、靶向治疗选择的核心依据。传统病理诊断依赖病理医生在显微镜下逐片判读，耗时、主观差异大，且受限于人眼疲劳，漏检、误诊时有发生。随着数字病理（Whole-Slide Imaging, WSI）普及，海量大图数据急需自动化、可重复的精准分析工具。

【研究意义】

1. 临床价值：高精度细胞核分割可量化 Ki-67、有丝分裂计数、肿瘤浸润淋巴细胞（TILs）等指标，为癌症分型、疗效评估提供客观数字证据。
2. 科研价值：构建“细胞级”数字病理图谱，支撑多中心大规模回顾性研究，加速生物标志物发现。
3. 产业价值：嵌入病理 AI 辅助诊断系统，可减少 50% 以上人工阅片时间，缓解全球病理医生缺口（WHO 报告：当前缺口约 6 万人）。

【研究现状】

1. 早期阶段（2015–2018）
• 基于传统形态学、分水岭、Graph-Cut 等方法，在 2018 DSB（Data Science Bowl）核分割挑战中 Dice 仅 0.70–0.76。

2. 深度学习爆发（2019–2021）
• U-Net、Mask R-CNN、HRNet 等 CNN 架构将 Dice 提升到 0.82–0.85。
• 代表性工作：HoVer-Net（2019，MICCAI）同时完成核分割、分类与实例分离，Dice 0.844；Cellpose（2020，Nature Methods）利用全局流场，跨器官泛化性能优异。
3. Transformer 时代（2022–至今）
• Swin-Transformer、TransUNet、Swin-UNETR 等模型引入长程依赖与多尺度注意力，在 MoNuSeg 公开集上 Dice 突破 0.90。
• 2023 年 CAMELYON 17+MoNuSeg 联合榜单显示，Swin-UNETR + Test-Time Augmentation 可达 0.92–0.93，逼近人类专家一致性上限（≈0.94）。
• 研究热点：
– 大模型自监督预训练（Masked Autoencoder）解决小样本问题；
– 多任务联合学习（分割+分类+计数）；
– 轻量级部署（MobileViT、Edge-TPU 量化）适应基层医院。

综上，细胞核分割领域已完成从传统算法到 CNN、再到 Transformer 的三级跳，当前研究正朝着“超高精度 + 跨器官泛化 + 端云协同”的新目标迈进。

# 1_otsu.py
import cv2
img = cv2.imread('3.png', 0)
_, bw = cv2.threshold(img, 0, 255, cv2.THRESH_OTSU + cv2.THRESH_BINARY)
# 去小面积
nb, labels, stats, _ = cv2.connectedComponentsWithStats(bw, connectivity=8)
for i in range(1, nb):
if stats[i, cv2.CC_STAT_AREA] < 64:
bw[labels == i] = 0
cv2.imwrite('otsu_mask.png', bw)

# 1_otsu.py
import cv2
img  = cv2.imread('3.png', 0)
_, bw = cv2.threshold(img, 0, 255, cv2.THRESH_OTSU + cv2.THRESH_BINARY)
# 去小面积
nb, labels, stats, _ = cv2.connectedComponentsWithStats(bw, connectivity=8)
for i in range(1, nb):
    if stats[i, cv2.CC_STAT_AREA] < 64:
        bw[labels == i] = 0
cv2.imwrite('otsu_mask.png', bw)

# 3_kmeans.py
import cv2
import numpy as np

img = cv2.imread('3.png')
Z = img.reshape((-1, 3)).astype(np.float32)
criteria = (cv2.TERM_CRITERIA_EPS + cv2.TERM_CRITERIA_MAX_ITER, 10, 1.0)
_, label, center = cv2.kmeans(Z, 3, None, criteria, 10, cv2.KMEANS_RANDOM_CENTERS)
seg = np.uint8(center[label.flatten()]).reshape(img.shape)
cv2.imwrite('kmeans_seg.png', seg)

# 2_watershed.py
import cv2
import numpy as np

img = cv2.imread('3.png')
gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
_, bin = cv2.threshold(gray, 0, 255, cv2.THRESH_OTSU + cv2.THRESH_BINARY)

# 距离变换 → 前景标记
dist = cv2.distanceTransform(bin, cv2.DIST_L2, 5)
_, fg = cv2.threshold(dist, 0.5 * dist.max(), 255, 0)
fg = np.uint8(fg)
# 背景标记
unknown = cv2.subtract(cv2.dilate(bin, np.ones((3, 3))), fg)
_, markers = cv2.connectedComponents(fg)
markers += 1
markers[unknown == 255] = 0
# 分水岭
markers = cv2.watershed(img, markers)
img[markers == -1] = [0, 0, 255]          # 边界涂红
cv2.imwrite('ws_mask.png', np.where(markers > 1, 255, 0).astype(np.uint8))
cv2.imwrite('ws_overlay.png', img)

import cv2
import numpy as np
from skimage import io, color, segmentation, feature, graph
from skimage.filters import sobel
from scipy import ndimage as ndi
import time

# --- 辅助函数：使用 OpenCV 并排显示多张图像 ---
def show_images_cv2(image_list, titles, window_name="Comparison"):
    """
    将多张图像水平拼接，并在一个窗口中显示。
    
    参数:
        image_list (list): 包含要显示的图像（NumPy数组）的列表。
        titles (list): 每张图像的标题列表。
        window_name (str): 显示窗口的名称。
    """
    if not image_list:
        print("图像列表为空，无法显示。")
        return

    # 确保所有图像都是3通道BGR格式，以便统一处理
    display_images = []
    for img in image_list:
        # 如果是灰度图 (2D)，转换为BGR (3D)
        if len(img.shape) == 2:
            img_bgr = cv2.cvtColor(img, cv2.COLOR_GRAY2BGR)
        # 如果是带调色板的分割图 (例如, 使用了cmap='tab20')，需要转换为BGR
        elif img.shape[-1] == 4: # RGBA
            img_bgr = cv2.cvtColor(img, cv2.COLOR_RGBA2BGR)
        elif img.shape[-1] == 3: # RGB
            img_bgr = cv2.cvtColor(img, cv2.COLOR_RGB2BGR)
        else:
            img_bgr = img # 假设已经是BGR
            
        # 在图像上添加标题
        title = titles[len(display_images)]
        cv2.putText(img_bgr, title, (10, 30), cv2.FONT_HERSHEY_SIMPLEX, 0.8, (255, 255, 255), 2, cv2.LINE_AA)
        display_images.append(img_bgr)

    # 水平拼接所有图像
    final_image = np.hstack(display_images)
    
    # 显示拼接后的大图
    cv2.imshow(window_name, final_image)
    print(f"\n显示窗口 '{window_name}'。按任意键关闭窗口并继续...")
    cv2.waitKey(0) # 等待用户按键
    cv2.destroyAllWindows()


# --- 1. 加载并预处理图像 ---
# 从URL加载图像
image_url = '3.png'
image = io.imread(image_url)

# 转换为灰度图
gray_image = color.rgb2gray(image)

# 为了在OpenCV中显示，我们保留一个BGR版本
image_bgr = cv2.cvtColor(image, cv2.COLOR_RGB2BGR)


# --- 2. 方法一：传统形态学操作 ---
def segment_with_morphology(img_gray):
    start_time = time.time()
    
    # 1. 二值化（使用Otsu's方法自动确定阈值）
    # 注意：cv2.threshold需要8位单通道图像，所以需要先转换
    img_gray_uint8 = (img_gray * 255).astype(np.uint8)
    _, binary = cv2.threshold(img_gray_uint8, 0, 255, cv2.THRESH_BINARY_INV + cv2.THRESH_OTSU)
    
    # 2. 形态学操作
    # 定义一个结构元素（核）
    kernel = np.ones((3, 3), np.uint8)
    
    # 去除噪声（开运算：先腐蚀后膨胀）
    cleaned = cv2.morphologyEx(binary, cv2.MORPH_OPEN, kernel, iterations=2)
    
    # 确保前景物体是实心的（闭运算：先膨胀后腐蚀）
    sure_bg = cv2.morphologyEx(cleaned, cv2.MORPH_CLOSE, kernel, iterations=10)
    
    end_time = time.time()
    
    # 返回处理后的二值图，用于显示
    return sure_bg, end_time - start_time

morph_result, morph_time = segment_with_morphology(gray_image)
# 保存结果
cv2.imwrite('morphology_result.png', morph_result)


# --- 3. 方法二：分水岭算法 ---
def segment_with_watershed(img_rgb, img_gray):
    start_time = time.time()
    
    # 1. 使用形态学操作找到确定的背景区域
    img_gray_uint8 = (img_gray * 255).astype(np.uint8)
    _, binary = cv2.threshold(img_gray_uint8, 0, 255, cv2.THRESH_BINARY_INV + cv2.THRESH_OTSU)
    kernel = np.ones((3,3), np.uint8)
    sure_bg = cv2.dilate(binary, kernel, iterations=3)
    
    # 2. 使用距离变换找到确定的前景区域
    dist_transform = cv2.distanceTransform(binary, cv2.DIST_L2, 5)
    _, sure_fg = cv2.threshold(dist_transform, 0.5 * dist_transform.max(), 255, 0)
    sure_fg = np.uint8(sure_fg)
    
    # 3. 找到未知区域（既不是前景也不是背景）
    unknown = cv2.subtract(sure_bg, sure_fg)
    
    # 4. 创建标记（Marker）
    # 对前景物体进行连通组件标记
    _, markers = cv2.connectedComponents(sure_fg)
    
    # 为所有标记加1，确保背景为0，未知区域为0
    markers = markers + 1
    
    # 将未知区域标记为0
    markers[unknown == 255] = 0
    
    # 5. 应用分水岭算法
    # OpenCV的分水岭需要输入一个3通道的彩色图像
    img_bgr_for_ws = cv2.cvtColor(img_rgb, cv2.COLOR_RGB2BGR)
    markers = cv2.watershed(img_bgr_for_ws, markers)
    
    # 6. 创建彩色分割图用于显示
    # 分水岭结果中，-1表示边界，正整数表示不同区域
    watershed_result = color.label2rgb(markers, bg_label=0)
    
    end_time = time.time()
    return watershed_result, end_time - start_time

watershed_result, watershed_time = segment_with_watershed(image, gray_image)
# 保存结果
# label2rgb返回的是RGB，需要转换为BGR再保存
watershed_result_bgr = cv2.cvtColor((watershed_result * 255).astype(np.uint8), cv2.COLOR_RGB2BGR)
cv2.imwrite('watershed_result.png', watershed_result_bgr)


# --- 4. 方法三：Graph-Cut 算法 ---
def segment_with_graphcut(img_rgb, img_gray):
    start_time = time.time()
    
    # 1. 定义一个矩形区域，框住所有前景物体
    # 格式为 (x, y, w, h)
    h, w = img_rgb.shape[:2]
    rect = (20, 20, w - 40, h - 40) # 留出一些边界
    
    # 2. 初始化前景和背景模型
    bgd_model = np.zeros((1, 65), np.float64)
    fgd_model = np.zeros((1, 65), np.float64)
    
    # 3. 运行GrabCut算法
    # mask的取值: 0(BG), 1(FG), 2(可能BG), 3(可能FG)
    # 初始时，矩形内部为可能前景(3)，外部为背景(0)
    mask = np.zeros(img_gray.shape, np.uint8)
    cv2.grabCut(img_rgb, mask, rect, bgd_model, fgd_model, 5, cv2.GC_INIT_WITH_RECT)
    
    # 4. 处理mask，得到最终的二值分割结果
    # 将“可能前景”和“确定前景”都设为1（白色），其余设为0（黑色）
    mask2 = np.where((mask == 2) | (mask == 0), 0, 1).astype('uint8')
    
    # 5. 将mask应用到原图上，得到分割后的彩色图像
    segmented_image = img_rgb * mask2[:, :, np.newaxis]
    
    end_time = time.time()
    return segmented_image, end_time - start_time

graphcut_result, graphcut_time = segment_with_graphcut(image, gray_image)
# 保存结果
graphcut_result_bgr = cv2.cvtColor((graphcut_result * 255).astype(np.uint8), cv2.COLOR_RGB2BGR)
cv2.imwrite('graphcut_result.png', graphcut_result_bgr)


# --- 5. 结果展示与对比分析 ---

# 准备要显示的图像列表和标题
# 注意：OpenCV默认使用BGR，但我们的原始image是RGB，需要转换
images_to_show = [
    image_bgr, 
    morph_result, 
    watershed_result, 
    graphcut_result_bgr
]
titles = [
    "Original Image", 
    "Morphology Result", 
    "Watershed Result", 
    "Graph-Cut Result"
]

# 使用我们自定义的OpenCV函数来显示所有结果
show_images_cv2(images_to_show, titles, "Image Segmentation Comparison")

# 打印对比分析
print("\n" + "="*50)
print("         图像分割方法对比分析")
print("="*50)
print(f"{'方法':<20} | {'耗时 (秒)':<15} | {'主要特点'}")
print("-"*50)
print(f"{'传统形态学':<20} | {morph_time:<15.4f} | {'速度快，简单，适合高对比度、形状规则的物体，难以处理粘连。'}")
print(f"{'分水岭算法':<20} | {watershed_time:<15.4f} | {'能有效分离粘连物体，但对噪声和初始标记敏感，易导致过分割。'}")
print(f"{'Graph-Cut (GrabCut)':<20} | {graphcut_time:<15.4f} | {'精度高，能处理复杂背景和弱边界，但速度慢，需要用户交互（如框选）。'}")
print("="*50)
print("\n总结:")
print("- **传统形态学** 是最基础的方法，像一个“工匠”，用简单的“凿子”（腐蚀、膨胀）来修整图像。它非常快，但只适用于简单的场景。")
print("- **分水岭算法** 像一个“地理学家”，将图像看作地形图，从“山谷”开始“淹没”，从而分出不同的“流域”（物体）。它解决了粘连问题，但需要好的“起点”（标记），否则会“分得太细”（过分割）。")
print("- **Graph-Cut** 像一个“精明的裁判”，它把像素看作图的节点，通过计算能量（代价）来决定每个像素属于前景还是背景。它非常强大，能做出最符合直觉的分割，但计算量大，并且需要一个“初步判断”（如矩形框）来引导它。")
print("\n所有结果图像已保存为 morphology_result.png, watershed_result.png, graphcut_result.png")

收录于