欢迎来到《计算机视觉实战》系列教程的第四章。在前三章中,我们学习了图像基础、OpenCV核心操作以及滤波与边缘检测。本章我们将深入探索形态学操作和图像分割技术,这是计算机视觉中连接底层处理与高层理解的关键桥梁。

想象一下,当你用边缘检测提取出物体的轮廓后,如何才能把这个轮廓变成一个完整的区域?如何从复杂的场景中分离出我们感兴趣的目标?这些就是形态学操作和图像分割要解决的问题。


1. 环境声明

  • Python版本Python 3.12+
  • PyTorch版本PyTorch 2.2+
  • OpenCV版本OpenCV 4.10+
  • scikit-image版本0.22+
  • NumPy版本1.26+
  • 操作系统:Windows / macOS / Linux (通用)

2. 形态学操作基础

2.1 什么是形态学

形态学(Morphology)一词源于希腊语,本意是研究生物体形式和结构的学科。在图像处理中,形态学是一种基于形状的信息处理方法,它使用预定义的结构元素(Structuring Element)在图像中滑动,根据结构元素与图像区域的匹配程度来执行各种操作。

核心思想:把图像看作一个集合,把结构元素看作一个"探针"。当探针在图像中移动时,记录探针与图像的契合程度,这就是形态学操作的本质。

类比理解:就像用不同形状的印章去盖在纸上,观察印章与纸面的接触情况。如果印章是圆形的,我们就能知道纸上哪些区域是圆形的;如果是方形的,就能检测方形区域。

2.2 结构元素详解

结构元素是形态学操作的核心,它定义了哪些像素需要被考虑以及它们之间的空间关系。

import cv2
import numpy as np
import matplotlib.pyplot as plt

# 创建不同形状的结构元素
# 1. 矩形结构元素
rect_kernel = cv2.getStructuringElement(cv2.MORPH_RECT, (5, 5))
print(f"矩形结构元素:\n{rect_kernel}")

# 2. 椭圆结构元素
ellipse_kernel = cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (5, 5))
print(f"\n椭圆结构元素:\n{ellipse_kernel}")

# 3. 十字形结构元素
cross_kernel = cv2.getStructuringElement(cv2.MORPH_CROSS, (5, 5))
print(f"\n十字形结构元素:\n{cross_kernel}")

# 自定义结构元素
custom_kernel = np.array([
    [0, 1, 0],
    [1, 1, 1],
    [0, 1, 0]
], dtype=np.uint8)
print(f"\n自定义十字结构元素:\n{custom_kernel}")

# 可视化结构元素
fig, axes = plt.subplots(1, 3, figsize=(12, 4))

axes[0].imshow(rect_kernel, cmap='gray')
axes[0].set_title('矩形结构元素 (5x5)')
axes[0].axis('off')

axes[1].imshow(ellipse_kernel, cmap='gray')
axes[1].set_title('椭圆结构元素 (5x5)')
axes[1].axis('off')

axes[2].imshow(cross_kernel, cmap='gray')
axes[2].set_title('十字形结构元素 (5x5)')
axes[2].axis('off')

plt.tight_layout()
plt.show()

结构元素的选择原则

  • 平滑区域:使用较大结构元素,可以快速去除噪声
  • 细节丰富:使用较小结构元素,保留更多细节
  • 形状匹配:根据目标形状选择相似的结构元素
  • 十字形:适合检测对角线或T形交叉
  • 矩形:适合去除方形噪声

3. 基本形态学运算

3.1 腐蚀(Erosion)

腐蚀是最基本的形态学操作之一。它的原理是:当结构元素完全包含在前景像素中时,保留中心像素;否则去除中心像素。

通俗理解:就像侵蚀作用,边缘被一点点"腐蚀"掉。图像中的白色区域会收缩,黑色区域会扩大。

应用场景

  • 去除小的白色噪声(如盐粒噪声)
  • 分离两个接触的物体
  • 去除微小突起
import cv2
import numpy as np
import matplotlib.pyplot as plt

# 创建测试图像:包含一个矩形和一个圆形
img = np.zeros((300, 400), dtype=np.uint8)

# 画一个矩形
cv2.rectangle(img, (50, 50), (180, 150), 255, -1)

# 画一个圆形
cv2.circle(img, (280, 200), 60, 255, -1)

# 添加噪声
noise = np.random.randint(0, 2, img.shape, dtype=np.uint8) * 255
img_with_noise = cv2.add(img, noise)

# 应用腐蚀
kernel = np.ones((5, 5), dtype=np.uint8)

# 不同迭代次数的腐蚀效果
erosion_1 = cv2.erode(img, kernel, iterations=1)
erosion_2 = cv2.erode(img, kernel, iterations=2)
erosion_3 = cv2.erode(img, kernel, iterations=3)

# 可视化
fig, axes = plt.subplots(2, 3, figsize=(15, 10))

axes[0, 0].imshow(img, cmap='gray')
axes[0, 0].set_title('原始图像')
axes[0, 0].axis('off')

axes[0, 1].imshow(img_with_noise, cmap='gray')
axes[0, 1].set_title('带噪声图像')
axes[0, 1].axis('off')

axes[0, 2].imshow(erosion_1, cmap='gray')
axes[0, 2].set_title('腐蚀 1次')
axes[0, 2].axis('off')

axes[1, 0].imshow(erosixon_2, cmap='gray')
axes[1, 0].set_title('腐蚀 2次')
axes[1, 0].axis('off')

axes[1, 1].imshow(erosion_3, cmap='gray')
axes[1, 1].set_title('腐蚀 3次')
axes[1, 1].axis('off')

# 用噪声图像测试去噪效果
noise_erosion = cv2.erode(img_with_noise, kernel, iterations=1)
axes[1, 2].imshow(noise_erosion, cmap='gray')
axes[1, 2].set_title('噪声图像腐蚀去噪')
axes[1, 2].axis('off')

plt.tight_layout()
plt.show()

print(f"原始像素数: {np.sum(img > 0)}")
print(f"腐蚀1次像素数: {np.sum(erosion_1 > 0)}")
print(f"腐蚀2次像素数: {np.sum(erosion_2 > 0)}")
print(f"腐蚀3次像素数: {np.sum(erosion_3 > 0)}")

3.2 膨胀(Dilation)

膨胀是腐蚀的反操作。当结构元素与前景像素有任何重叠时,保留中心像素。

通俗理解:就像膨胀作用,边缘向外扩张。图像中的白色区域会扩大,黑色区域会缩小。

应用场景

  • 填补物体内部的小孔洞
  • 连接相邻的物体
  • 增强物体边缘
import cv2
import numpy as np
import matplotlib.pyplot as plt

# 创建测试图像:包含一个环形
img_ring = np.zeros((300, 400), dtype=np.uint8)
cv2.circle(img_ring, (200, 150), 80, 255, -1)
cv2.circle(img_ring, (200, 150), 40, 0, -1)  # 挖空中心

# 应用膨胀
kernel = np.ones((5, 5), dtype=np.uint8)

# 不同迭代次数的膨胀效果
dilation_1 = cv2.dilate(img_ring, kernel, iterations=1)
dilation_2 = cv2.dilate(img_ring, kernel, iterations=2)
dilation_3 = cv2.dilate(img_ring, kernel, iterations=3)

# 可视化
fig, axes = plt.subplots(2, 2, figsize=(12, 12))

axes[0, 0].imshow(img_ring, cmap='gray')
axes[0, 0].set_title('原始环形图像')
axes[0, 0].axis('off')

axes[0, 1].imshow(dilation_1, cmap='gray')
axes[0, 1].set_title('膨胀 1次 - 环变粗')
axes[0, 1].axis('off')

axes[1, 0].imshow(dilation_2, cmap='gray')
axes[1, 0].set_title('膨胀 2次 - 环更粗')
axes[1, 0].axis('off')

axes[1, 1].imshow(dilation_3, cmap='gray')
axes[1, 1].set_title('膨胀 3次 - 中心被填充')
axes[1, 1].axis('off')

plt.tight_layout()
plt.show()

# 验证膨胀如何填补孔洞
print(f"原始图像中心是否为空: {img_ring[150, 200] == 0}")
print(f"膨胀3次后中心是否为空: {dilation_3[150, 200] == 0}")

3.3 开运算与闭运算

开运算和闭运算是腐蚀和膨胀的组合,它们是形态学中最常用的操作。

开运算(Opening) = 腐蚀 → 膨胀

  • 效果:去除小的亮点(噪声),保持物体的大致形状
  • 特性:平滑物体轮廓,断开细小的连接

闭运算(Closing) = 膨胀 → 腐蚀

  • 效果:填补小的暗点(孔洞),保持物体的大致形状
  • 特性:平滑物体轮廓,连接细小的断裂
import cv2
import numpy as np
import matplotlib.pyplot as plt

# 创建测试图像
img = np.zeros((300, 400), dtype=np.uint8)

# 画一个有缺口的矩形
cv2.rectangle(img, (50, 50), (180, 150), 255, -1)
# 添加白色噪声点
cv2.circle(img, (100, 80), 5, 0, -1)  # 小暗点

# 画一个有孔洞的圆形
cv2.circle(img, (280, 150), 70, 255, -1)
cv2.circle(img, (280, 150), 20, 0, -1)  # 中心孔洞

kernel = np.ones((7, 7), dtype=np.uint8)

# 开运算:先腐蚀后膨胀 - 去除小白点
opening = cv2.morphologyEx(img, cv2.MORPH_OPEN, kernel)

# 闭运算:先膨胀后腐蚀 - 填补小暗点和孔洞
closing = cv2.morphologyEx(img, cv2.MORPH_CLOSE, kernel)

# 可视化
fig, axes = plt.subplots(2, 2, figsize=(12, 12))

axes[0, 0].imshow(img, cmap='gray')
axes[0, 0].set_title('原始图像\n(有噪声点和孔洞)')
axes[0, 0].axis('off')

axes[0, 1].imshow(opening, cmap='gray')
axes[0, 1].set_title('开运算\n(去除了小白点)')
axes[0, 1].axis('off')

axes[1, 0].imshow(closing, cmap='gray')
axes[1, 0].set_title('闭运算\n(填补了孔洞)')
axes[1, 0].axis('off')

# 梯度运算:膨胀 - 腐蚀,得到轮廓
gradient = cv2.morphologyEx(img, cv2.MORPH_GRADIENT, kernel)
axes[1, 1].imshow(gradient, cmap='gray')
axes[1, 1].set_title('形态学梯度\n(提取轮廓)')
axes[1, 1].axis('off')

plt.tight_layout()
plt.show()

print("开运算原理:腐蚀 -> 膨胀")
print("闭运算原理:膨胀 -> 腐蚀")
print("梯度运算:膨胀 - 腐蚀")

3.4 其他形态学运算

import cv2
import numpy as np

# 创建测试图像
img = np.zeros((200, 200), dtype=np.uint8)
cv2.rectangle(img, (40, 40), (160, 160), 255, -1)

kernel = np.ones((11, 11), dtype=np.uint8)

# 顶帽变换(Top Hat):原图 - 开运算
# 提取小目标
tophat = cv2.morphologyEx(img, cv2.MORPH_TOPHAT, kernel)

# 黑帽变换(Black Hat):闭运算 - 原图
# 提取小孔洞
blackhat = cv2.morphologyEx(img, cv2.MORPH_BLACKHAT, kernel)

print(f"顶帽变换结果(亮点区域): {np.sum(tophat > 0)} 像素")
print(f"黑帽变换结果(暗点区域): {np.sum(blackhat > 0)} 像素")

4. 阈值分割

阈值分割是最简单也是最有效的图像分割方法。它的核心思想是将图像划分为前景和背景两部分,基于像素值的阈值来进行分类。

4.1 全局阈值法

固定阈值法是最简单的阈值方法,但需要手动设定阈值,不适合光照不均匀的图像。

import cv2
import numpy as np
import matplotlib.pyplot as plt

# 读取图像并转为灰度
img = cv2.imread('image.jpg')
gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)

# 固定阈值分割
threshold_value = 127  # 手动设定的阈值
_, binary_fixed = cv2.threshold(gray, threshold_value, 255, cv2.THRESH_BINARY)

# 自动计算阈值
_, binary_auto = cv2.threshold(gray, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)

# 可视化
fig, axes = plt.subplots(2, 3, figsize=(15, 10))

axes[0, 0].imshow(cv2.cvtColor(img, cv2.COLOR_BGR2RGB))
axes[0, 0].set_title('原始图像')
axes[0, 0].axis('off')

axes[0, 1].imshow(gray, cmap='gray')
axes[0, 1].set_title('灰度图')
axes[0, 1].axis('off')

axes[0, 2].hist(gray.ravel(), 256, [0, 256])
axes[0, 2].axvline(x=127, color='r', linestyle='--', label='固定阈值 127')
axes[0, 2].set_title('灰度直方图')
axes[0, 2].legend()

axes[1, 0].imshow(binary_fixed, cmap='gray')
axes[1, 0].set_title(f'固定阈值分割 (T=127)')
axes[1, 0].axis('off')

axes[1, 1].imshow(binary_auto, cmap='gray')
axes[1, 1].set_title('Otsu自动阈值分割')
axes[1, 1].axis('off')

# 对比两种方法
diff = cv2.absdiff(binary_fixed, binary_auto)
axes[1, 2].imshow(diff, cmap='gray')
axes[1, 2].set_title('两种方法差异')
axes[1, 2].axis('off')

plt.tight_layout()
plt.show()

4.2 Otsu阈值法

Otsu算法是一种自动确定阈值的方法,它的核心思想是找到一个阈值,使得前景和背景的类间方差最大。

数学原理
设图像有L个灰度级,前景和背景的分割阈值为T。类间方差定义为:
σ2=P1⋅P2⋅(μ1−μ2)2\sigma^2 = P_1 \cdot P_2 \cdot (\mu_1 - \mu_2)^2σ2=P1P2(μ1μ2)2

其中P1和P2分别是前景和背景的概率,μ1和μ2是它们的均值。Otsu算法遍历所有可能的T,选择使σ²最大的T作为最优阈值。

import cv2
import numpy as np

def otsu_threshold_manual(gray_img):
    """手动实现Otsu阈值算法"""
    # 计算直方图
    hist = cv2.calcHist([gray_img], [0], None, [256], [0, 256]).flatten()
    hist = hist / hist.sum()  # 归一化

    # 总均值
    total_mean = np.sum(np.arange(256) * hist)

    # 遍历所有阈值,找到使类间方差最大的阈值
    max_variance = 0
    optimal_threshold = 0

    sum_bg = 0  # 背景灰度加权和
    weight_bg = 0  # 背景概率

    for t in range(256):
        # 背景部分 [0, t]
        weight_bg += hist[t]
        if weight_bg == 0:
            continue

        # 前景部分 [t+1, 255]
        weight_fg = 1 - weight_bg
        if weight_fg == 0:
            break

        sum_bg += t * hist[t]

        # 计算均值
        mean_bg = sum_bg / weight_bg
        mean_fg = (total_mean - sum_bg) / weight_fg

        # 计算类间方差
        variance = weight_bg * weight_fg * (mean_bg - mean_fg) ** 2

        if variance > max_variance:
            max_variance = variance
            optimal_threshold = t

    return optimal_threshold

# 读取图像
gray = cv2.imread('image.jpg', cv2.IMREAD_GRAYSCALE)

# 使用手动实现的Otsu算法
manual_otsu_t = otsu_threshold_manual(gray)

# 使用OpenCV的Otsu算法
_, cv2_otsu = cv2.threshold(gray, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)

print(f"手动Otsu阈值: {manual_otsu_t}")
print(f"OpenCV Otsu阈值: {cv2_otsu}")
print(f"两种方法差异: {abs(manual_otsu_t - cv2_otsu)}")

4.3 自适应阈值法

当图像光照不均匀时,全局阈值法效果会很差。自适应阈值法通过计算局部区域的阈值来解决这个问题。

import cv2
import numpy as np
import matplotlib.pyplot as plt

# 读取图像
gray = cv2.imread('uneven_lighting.jpg', cv2.IMREAD_GRAYSCALE)

# 全局阈值(效果不好)
_, global_thresh = cv2.threshold(gray, 127, 255, cv2.THRESH_BINARY)

# 自适应阈值 - 均值方法
adaptive_mean = cv2.adaptiveThreshold(
    gray, 255,
    cv2.ADAPTIVE_THRESH_MEAN_C,
    cv2.THRESH_BINARY,
    blockSize=11,  # 必须是奇数
    C=2  # 偏移量
)

# 自适应阈值 - 高斯方法
adaptive_gaussian = cv2.adaptiveThreshold(
    gray, 255,
    cv2.ADAPTIVE_THRESH_GAUSSIAN_C,
    cv2.THRESH_BINARY,
    blockSize=11,
    C=2
)

# 可视化
fig, axes = plt.subplots(2, 2, figsize=(12, 12))

axes[0, 0].imshow(gray, cmap='gray')
axes[0, 0].set_title('光照不均匀的灰度图')
axes[0, 0].axis('off')

axes[0, 1].imshow(global_thresh, cmap='gray')
axes[0, 1].set_title('全局阈值 (效果差)')
axes[0, 1].axis('off')

axes[1, 0].imshow(adaptive_mean, cmap='gray')
axes[1, 0].set_title('自适应阈值 (均值)')
axes[1, 0].axis('off')

axes[1, 1].imshow(adaptive_gaussian, cmap='gray')
axes[1, 1].set_title('自适应阈值 (高斯)')
axes[1, 1].axis('off')

plt.tight_layout()
plt.show()

print("自适应阈值的blockSize选择原则:")
print("- 太小:容易受噪声影响")
print("- 太大:会丢失细节")
print("- 通常选择11、15、21等奇数")
print("\nC值的选择原则:")
print("- C太小:阈值偏低,亮区增多")
print("- C太大:阈值偏高,暗区增多")

4.4 直方图分析法

import cv2
import numpy as np
import matplotlib.pyplot as plt

def analyze_histogram(gray_img):
    """分析直方图特征"""
    hist = cv2.calcHist([gray_img], [0], None, [256], [0, 256]).flatten()

    # 计算基本统计量
    mean_val = np.mean(gray_img)
    std_val = np.std(gray_img)
    median_val = np.median(gray_img)

    # 找到峰值
    peaks, _ = find_peaks(hist, height=np.max(hist) * 0.1, distance=20)

    return {
        'histogram': hist,
        'mean': mean_val,
        'std': std_val,
        'median': median_val,
        'peaks': peaks
    }

# 使用skimage的峰值检测
from scipy.signal import find_peaks

gray = cv2.imread('image.jpg', cv2.IMREAD_GRAYSCALE)
analysis = analyze_histogram(gray)

# 可视化
fig, axes = plt.subplots(1, 2, figsize=(15, 5))

axes[0].imshow(gray, cmap='gray')
axes[0].set_title('灰度图')
axes[0].axis('off')

axes[1].plot(analysis['histogram'])
axes[1].axvline(x=analysis['mean'], color='r', linestyle='--', label=f'均值: {analysis["mean"]:.1f}')
axes[1].axvline(x=analysis['median'], color='g', linestyle='--', label=f'中位数: {analysis["median"]:.1f}')
for peak in analysis['peaks']:
    axes[1].axvline(x=peak, color='b', linestyle=':', alpha=0.7)
axes[1].set_title('灰度直方图与峰值分析')
axes[1].set_xlabel('灰度值')
axes[1].set_ylabel('像素数量')
axes[1].legend()

plt.tight_layout()
plt.show()

print(f"峰值位置: {analysis['peaks']}")
print(f"均值: {analysis['mean']:.2f}")
print(f"标准差: {analysis['std']:.2f}")
print(f"中位数: {analysis['median']:.2f}")

5. 连通域分析

连通域分析是分割后处理的重要步骤,它可以将二值图像中的不同区域标记出来,便于后续的目标识别和计数。

5.1 连通性的定义

4连通(4-connectivity):只考虑上下左右四个邻居
8连通(8-connectivity):考虑上下左右和四个对角线共八个邻居

import cv2
import numpy as np
import matplotlib.pyplot as plt

# 创建测试图像:包含两个接触的物体
img = np.zeros((100, 100), dtype=np.uint8)

# 物体1(左上)
cv2.rectangle(img, (10, 10), (45, 45), 255, -1)

# 物体2(右下)
cv2.rectangle(img, (55, 55), (90, 90), 255, -1)

# 两个物体在角落接触
cv2.circle(img, (50, 50), 8, 255, -1)

# 连通域标记
num_labels_4, labels_4, stats_4, centroids_4 = cv2.connectedComponentsWithStats(img, connectivity=4)
num_labels_8, labels_8, stats_8, centroids_8 = cv2.connectedComponentsWithStats(img, connectivity=8)

# 可视化
fig, axes = plt.subplots(1, 3, figsize=(15, 5))

axes[0].imshow(img, cmap='gray')
axes[0].set_title('原始二值图像')
axes[0].axis('off')

# 使用不同颜色显示标签
labels_4_color = np.zeros((*labels_4.shape, 3), dtype=np.uint8)
labels_4_color[labels_4 == 1] = [255, 0, 0]
labels_4_color[labels_4 == 2] = [0, 255, 0]

axes[1].imshow(labels_4_color)
axes[1].set_title(f'4连通 - {num_labels_4 - 1}个物体')
axes[1].axis('off')

labels_8_color = np.zeros((*labels_8.shape, 3), dtype=np.uint8)
labels_8_color[labels_8 == 1] = [255, 0, 0]

axes[2].imshow(labels_8_color)
axes[2].set_title(f'8连通 - {num_labels_8 - 1}个物体')
axes[2].axis('off')

plt.tight_layout()
plt.show()

print(f"4连通检测到的物体数: {num_labels_4 - 1}")
print(f"8连通检测到的物体数: {num_labels_8 - 1}")
print("\n4连通会把接触的角点算作两个物体")
print("8连通会把接触的角点算作同一个物体")

5.2 连通域属性分析

import cv2
import numpy as np

def analyze_connected_components(binary_img):
    """分析连通域的各种属性"""
    num_labels, labels, stats, centroids = cv2.connectedComponentsWithStats(
        binary_img, connectivity=8
    )

    results = []

    for i in range(1, num_labels):  # 跳过背景(label=0)
        x, y, w, h, area = stats[i]
        cx, cy = centroids[i]

        # 计算边界框
        bbox = (x, y, x + w, y + h)

        # 计算周长(近似)
        perimeter = cv2.arcLength(
            cv2.findContours(
                (labels == i).astype(np.uint8),
                cv2.RETR_EXTERNAL,
                cv2.CHAIN_APPROX_SIMPLE
            )[0], True
        )

        # 计算紧凑度 = 4π * 面积 / 周长²
        compactness = 4 * np.pi * area / (perimeter ** 2) if perimeter > 0 else 0

        results.append({
            'label': i,
            'bbox': bbox,
            'area': area,
            'perimeter': perimeter,
            'centroid': (cx, cy),
            'compactness': compactness
        })

    return results

# 创建测试图像
img = np.zeros((300, 400), dtype=np.uint8)

# 添加不同形状的物体
cv2.rectangle(img, (30, 30), (80, 80), 255, -1)  # 方形
cv2.circle(img, (180, 60), 30, 255, -1)  # 圆形
cv2.circle(img, (320, 60), 30, 255, -1)  # 圆形(接近)

# 细长物体
cv2.rectangle(img, (50, 150), (60, 250), 255, -1)

# 不规则物体
pts = np.array([[200, 130], [250, 150], [270, 200], [230, 220], [180, 180]], np.int32)
cv2.fillPoly(img, [pts], 255)

# 分析连通域
components = analyze_connected_components(img)

# 打印结果
print("检测到的物体属性:")
print("-" * 60)
for comp in components:
    print(f"物体 {comp['label']}:")
    print(f"  边界框: {comp['bbox']}")
    print(f"  面积: {comp['area']} 像素")
    print(f"  周长: {comp['perimeter']:.2f} 像素")
    print(f"  质心: ({comp['centroid'][0]:.1f}, {comp['centroid'][1]:.1f})")
    print(f"  紧凑度: {comp['compactness']:.3f} (1=完美圆形, 0=细长)")
    print("-" * 60)

5.3 去除小目标与区域填充

import cv2
import numpy as np
import matplotlib.pyplot as plt

# 创建测试图像
img = np.zeros((400, 400), dtype=np.uint8)

# 添加大物体
cv2.rectangle(img, (50, 50), (180, 150), 255, -1)
cv2.circle(img, (300, 300), 60, 255, -1)

# 添加小的噪声点
for _ in range(50):
    x, y = np.random.randint(0, 400, 2)
    cv2.circle(img, (x, y), 3, 255, -1)

# 定义最小面积阈值
min_area = 200

# 方法1:去除小目标
num_labels, labels, stats, centroids = cv2.connectedComponentsWithStats(img, connectivity=8)

# 创建掩膜,只保留大于最小面积的目标
mask = np.zeros_like(img)
for i in range(1, num_labels):
    if stats[i, cv2.CC_STAT_AREA] >= min_area:
        mask[labels == i] = 255

# 方法2:填充孔洞
filled_img = img.copy()
h, w = img.shape[:2]
mask_flood = np.zeros((h + 2, w + 2), dtype=np.uint8)
cv2.floodFill(filled_img, mask_flood, (0, 0), 255)  # 从左上角填充

# 反转填充结果,得到原物体的孔洞
holes = cv2.bitwise_not(filled_img)
# 合并原物体和孔洞
filled_img = cv2.bitwise_or(img, holes)

# 可视化
fig, axes = plt.subplots(2, 2, figsize=(12, 12))

axes[0, 0].imshow(img, cmap='gray')
axes[0, 0].set_title('原始图像(包含噪声点)')
axes[0, 0].axis('off')

axes[0, 1].imshow(labels, cmap='nipy_spectral')
axes[0, 1].set_title(f'连通域标记({num_labels - 1}个区域)')
axes[0, 1].axis('off')

axes[1, 0].imshow(mask, cmap='gray')
axes[1, 0].set_title(f'去除小目标后(面积<{min_area}被去除)')
axes[1, 0].axis('off')

axes[1, 1].imshow(filled_img, cmap='gray')
axes[1, 1].set_title('填充孔洞后')
axes[1, 1].axis('off')

plt.tight_layout()
plt.show()

6. 分水岭算法

分水岭算法是一种基于地形学的分割算法,它将图像看作三维地形(x、y是空间坐标,z是灰度值),然后在灰度值的"山谷"处放置水,让水从高处流向低处,最终形成分水岭分割。

6.1 分水岭原理

核心思想

  1. 将图像视为地形学表面,灰度值代表海拔高度
  2. 在灰度最低点(谷底)打孔,水开始涌入
  3. 水逐渐填满山谷,不同水系相遇时建立分水岭
  4. 最终的分水岭线就是分割边界

问题:噪声会导致过度分割,解决方法是使用标记图(Marker)来控制分水岭的起始位置。

import cv2
import numpy as np
import matplotlib.pyplot as plt

def watershed_segmentation(image_path):
    """分水岭算法分割示例"""

    # 读取图像并转为灰度
    img = cv2.imread(image_path)
    gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)

    # 1. 阈值分割得到初始标记
    _, thresh = cv2.threshold(gray, 0, 255, cv2.THRESH_BINARY_INV + cv2.THRESH_OTSU)

    # 2. 形态学操作去除噪声
    kernel = np.ones((3, 3), dtype=np.uint8)
    opening = cv2.morphologyEx(thresh, cv2.MORPH_OPEN, kernel, iterations=2)

    # 3. 膨胀操作确定背景
    sure_bg = cv2.dilate(opening, kernel, iterations=3)

    # 4. 距离变换确定前景
    dist_transform = cv2.distanceTransform(opening, cv2.DIST_L2, 5)
    _, sure_fg = cv2.threshold(dist_transform, 0.7 * dist_transform.max(), 255, 0)

    # 5. 未知区域 = 背景 - 前景
    sure_fg = np.uint8(sure_fg)
    unknown = cv2.subtract(sure_bg, sure_fg)

    # 6. 创建标记图
    _, markers = cv2.connectedComponents(sure_fg)

    # 7. 给背景标记加1,背景为1,未知区域为0
    markers = markers + 1
    markers[unknown == 255] = 0

    # 8. 应用分水岭算法
    markers = cv2.watershed(img, markers)

    # 9. 绘制分割边界
    img[markers == -1] = [0, 0, 255]  # 红色标记边界

    return img, markers

# 使用硬币图像测试
img, markers = watershed_segmentation('coins.jpg')

# 可视化过程
fig, axes = plt.subplots(2, 3, figsize=(15, 10))

axes[0, 0].imshow(cv2.cvtColor(cv2.imread('coins.jpg'), cv2.COLOR_BGR2RGB))
axes[0, 0].set_title('原始图像')
axes[0, 0].axis('off')

_, thresh = cv2.threshold(cv2.cvtColor(cv2.imread('coins.jpg'), cv2.COLOR_BGR2GRAY),
                          0, 255, cv2.THRESH_BINARY_INV + cv2.THRESH_OTSU)
axes[0, 1].imshow(thresh, cmap='gray')
axes[0, 1].set_title('Otsu阈值分割')
axes[0, 1].axis('off')

kernel = np.ones((3, 3), dtype=np.uint8)
opening = cv2.morphologyEx(thresh, cv2.MORPH_OPEN, kernel, iterations=2)
axes[0, 2].imshow(opening, cmap='gray')
axes[0, 2].set_title('形态学开运算去噪')
axes[0, 2].axis('off')

dist_transform = cv2.distanceTransform(opening, cv2.DIST_L2, 5)
axes[1, 0].imshow(dist_transform, cmap='jet')
axes[1, 0].set_title('距离变换')
axes[1, 0].axis('off')

axes[1, 1].imshow(markers, cmap='nipy_spectral')
axes[1, 1].set_title('分水岭标记图')
axes[1, 1].axis('off')

axes[1, 2].imshow(cv2.cvtColor(img, cv2.COLOR_BGR2RGB))
axes[1, 2].set_title('分水岭分割结果')
axes[1, 2].axis('off')

plt.tight_layout()
plt.show()

6.2 交互式分水岭分割

import cv2
import numpy as np

class InteractiveWatershed:
    """交互式分水岭分割工具"""

    def __init__(self, image_path):
        self.img = cv2.imread(image_path)
        self.gray = cv2.cvtColor(self.img, cv2.COLOR_BGR2GRAY)
        self.markers = np.zeros(self.gray.shape, dtype=np.int32)
        self.current_marker = 1  # 当前标记值
        self.drawing = False

    def mouse_callback(self, event, x, y, flags, param):
        """鼠标回调函数"""
        if event == cv2.EVENT_LBUTTONDOWN:
            self.drawing = True
            cv2.circle(self.markers, (x, y), 5, self.current_marker, -1)
            cv2.circle(self.img, (x, y), 5, (0, 0, 255), -1)

        elif event == cv2.EVENT_MOUSEMOVE:
            if self.drawing:
                cv2.circle(self.markers, (x, y), 5, self.current_marker, -1)
                cv2.circle(self.img, (x, y), 5, (0, 0, 255), -1)

        elif event == cv2.EVENT_LBUTTONUP:
            self.drawing = False

    def add_background_marker(self):
        """添加背景标记(按2)"""
        self.current_marker = 1

    def add_foreground_marker(self):
        """添加前景标记(按3)"""
        self.current_marker = 2

    def run(self):
        """运行交互式分割"""
        cv2.namedWindow('Interactive Watershed')
        cv2.setMouseCallback('Interactive Watershed', self.mouse_callback)

        print("交互式分水岭分割:")
        print("  按 '1' - 添加背景标记")
        print("  按 '2' - 添加前景标记")
        print("  按 'w' - 运行分水岭")
        print("  按 'r' - 重置")
        print("  按 'q' - 退出")

        while True:
            cv2.imshow('Interactive Watershed', self.img)

            key = cv2.waitKey(1) & 0xFF

            if key == ord('1'):
                self.current_marker = 1
                print("背景标记模式")

            elif key == ord('2'):
                self.current_marker = 2
                print("前景标记模式")

            elif key == ord('w'):
                # 运行分水岭
                self.markers = self.markers + 1
                self.markers = cv2.watershed(self.img, self.markers)
                print(f"分水岭完成,检测到 {len(np.unique(self.markers)) - 2} 个区域")

            elif key == ord('r'):
                # 重置
                self.img = cv2.imread('image.jpg')
                self.markers = np.zeros(self.gray.shape, dtype=np.int32)
                print("重置完成")

            elif key == ord('q'):
                break

        cv2.destroyAllWindows()
        return self.markers

# 使用示例
# watershed = InteractiveWatershed('image.jpg')
# markers = watershed.run()

7. GrabCut算法

GrabCut是一种基于图割(Graph Cut)的交互式前景提取算法,它只需要少量的用户输入就能获得很好的分割效果。

7.1 GrabCut原理

算法流程

  1. 用户初始化:定义一个包含前景的矩形框
  2. 迭代优化:
    • 估计:基于当前标记,像素被分类为前景或背景
    • 更新:使用高斯混合模型(GMM)更新前景和背景的颜色分布
    • 割:使用最小割算法找到最优分割
  3. 重复迭代直到收敛
import cv2
import numpy as np
import matplotlib.pyplot as plt

def grabcut_segmentation(image_path, rect):
    """GrabCut前景提取"""

    img = cv2.imread(image_path)
    mask = np.zeros(img.shape[:2], dtype=np.uint8)

    # 创建GrabCut所需的数据结构
    bgd_model = np.zeros((1, 65), dtype=np.float64)
    fgd_model = np.zeros((1, 65), dtype=np.float64)

    # rect格式:(x, y, w, h)
    # 运行GrabCut算法
    cv2.grabCut(img, mask, rect, bgd_model, fgd_model,
                5,  # 迭代次数
                cv2.GC_INIT_WITH_RECT)  # 使用矩形初始化

    # 创建最终掩膜
    # 0 = 确定的背景
    # 1 = 确定的前景
    # 2 = 可能的背景
    # 3 = 可能的前景
    final_mask = np.where((mask == 2) | (mask == 0), 0, 1).astype('uint8')

    # 应用掩膜提取前景
    foreground = img * final_mask[:, :, np.newaxis]

    return img, final_mask, foreground

# 读取图像
img = cv2.imread('person.jpg')

# 定义包含前景的矩形 (x, y, w, h)
# 前景物体大致在这个矩形内
h, w = img.shape[:2]
rect = (int(w * 0.1), int(h * 0.1), int(w * 0.8), int(h * 0.8))

# 执行GrabCut
original, mask, foreground = grabcut_segmentation('person.jpg', rect)

# 可视化
fig, axes = plt.subplots(2, 2, figsize=(12, 12))

axes[0, 0].imshow(cv2.cvtColor(original, cv2.COLOR_BGR2RGB))
axes[0, 0].set_title('原始图像')
axes[0, 0].axis('off')
# 绘制矩形框
x, y, rw, rh = rect
cv2.rectangle(original, (x, y), (x + rw, y + rh), (0, 255, 0), 2)
axes[0, 0].imshow(cv2.cvtColor(original, cv2.COLOR_BGR2RGB))
axes[0, 0].set_title('初始化矩形')
axes[0, 0].axis('off')

axes[0, 1].imshow(mask, cmap='gray')
axes[0, 1].set_title('GrabCut掩膜')
axes[0, 1].axis('off')

axes[1, 0].imshow(cv2.cvtColor(foreground, cv2.COLOR_BGR2RGB))
axes[1, 0].set_title('提取的前景')
axes[1, 0].axis('off')

# 创建透明背景
bg_removed = original.copy()
bg_removed[mask == 0] = 255  # 白色背景
axes[1, 1].imshow(cv2.cvtColor(bg_removed, cv2.COLOR_BGR2RGB))
axes[1, 1].set_title('背景移除效果')
axes[1, 1].axis('off')

plt.tight_layout()
plt.show()

7.2 迭代式GrabCut refinement

import cv2
import numpy as np

def grabcut_with_iteration(image_path, init_rect, num_iterations=10):
    """迭代式GrabCut以获得更精细的结果"""

    img = cv2.imread(image_path)
    mask = np.zeros(img.shape[:2], dtype=np.uint8)
    bgd_model = np.zeros((1, 65), dtype=np.float64)
    fgd_model = np.zeros((1, 65), dtype=np.float64)

    # 第一次GrabCut
    cv2.grabCut(img, mask, init_rect, bgd_model, fgd_model,
                5, cv2.GC_INIT_WITH_RECT)

    # 获取当前的前景区域
    for _ in range(num_iterations - 1):
        # 创建更精确的初始掩膜
        new_mask = mask.copy()

        # 手动标记确定的前景和背景(基于当前结果)
        # 0 = 背景, 1 = 前景, 2 = 可能背景, 3 = 可能前景
        certain_bg = (mask == 0) | (mask == 2)
        certain_fg = (mask == 1) | (mask == 3)

        new_mask[certain_bg] = cv2.GC_BGD
        new_mask[certain_fg] = cv2.GC_FGD

        # 再次运行GrabCut
        cv2.grabCut(img, new_mask, init_rect, bgd_model, fgd_model,
                    3, cv2.GC_INIT_WITH_MASK)

        mask = new_mask

    # 创建最终掩膜
    final_mask = np.where((mask == 2) | (mask == 0), 0, 1).astype('uint8')

    return img, final_mask

# 使用示例
h, w = 480, 640  # 假设图像大小
rect = (50, 50, w - 100, h - 100)
img, mask = grabcut_with_iteration('person.jpg', rect, num_iterations=5)

print(f"前景像素数: {np.sum(mask == 1)}")
print(f"背景像素数: {np.sum(mask == 0)}")

8. 实战项目:文档扫描仪

结合所学知识,创建一个实用的文档扫描仪程序:

import cv2
import numpy as np
import matplotlib.pyplot as plt
from typing import Tuple, List, Optional

class DocumentScanner:
    """文档扫描仪:自动检测并校正文档"""

    def __init__(self, image_path: str):
        self.original = cv2.imread(image_path)
        if self.original is None:
            raise ValueError(f"无法读取图像: {image_path}")
        self.image = self.original.copy()
        self.gray = cv2.cvtColor(self.image, cv2.COLOR_BGR2GRAY)

    def preprocess(self) -> np.ndarray:
        """预处理:去噪和增强对比度"""
        # 高斯模糊去噪
        blurred = cv2.GaussianBlur(self.gray, (5, 5), 0)

        # 自适应阈值增强
        adaptive_thresh = cv2.adaptiveThreshold(
            blurred, 255,
            cv2.ADAPTIVE_THRESH_GAUSSIAN_C,
            cv2.THRESH_BINARY,
            11, 2
        )

        # 形态学操作去噪
        kernel = np.ones((3, 3), dtype=np.uint8)
        cleaned = cv2.morphologyEx(adaptive_thresh, cv2.MORPH_CLOSE, kernel)

        return cleaned

    def find_document_contour(self, binary_img: np.ndarray) -> Optional[np.ndarray]:
        """找到文档边缘"""

        # 轮廓检测
        contours, _ = cv2.findContours(
            binary_img, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE
        )

        if not contours:
            return None

        # 找到最大轮廓(应该是文档)
        max_area = 0
        max_contour = None

        for contour in contours:
            area = cv2.contourArea(contour)
            if area > max_area:
                max_area = area
                max_contour = contour

        if max_contour is None:
            return None

        # 多边形逼近
        peri = cv2.arcLength(max_contour, True)
        approx = cv2.approxPolyDP(max_contour, 0.02 * peri, True)

        # 应该是四边形
        if len(approx) == 4:
            return approx.reshape(4, 2)
        else:
            # 如果不是四边形,尝试找最接近的四边形
            hull = cv2.convexHull(max_contour)
            hull_approx = cv2.approxPolyDP(hull, 0.02 * cv2.arcLength(hull, True), True)

            if len(hull_approx) >= 4:
                # 返回最接近四边形的四个角点
                return hull_approx[:4].reshape(4, 2)

        return approx.reshape(4, 2) if len(approx) >= 4 else None

    def order_points(self, pts: np.ndarray) -> np.ndarray:
        """将四个角点按左上、右上、右下、左下顺序排列"""

        rect = np.zeros((4, 2), dtype=np.float32)

        # 按x坐标排序
        s = pts.sum(axis=1)
        rect[0] = pts[np.argmin(s)]      # 左上
        rect[2] = pts[np.argmax(s)]      # 右下

        # 按y坐标排序
        diff = np.diff(pts, axis=1)
        rect[1] = pts[np.argmin(diff)]   # 右上
        rect[3] = pts[np.argmax(diff)]   # 左下

        return rect

    def perspective_transform(self, pts: np.ndarray,
                               target_size: Tuple[int, int] = (800, 600)) -> np.ndarray:
        """透视变换校正文档"""

        rect = self.order_points(pts.astype(np.float32))

        # 计算目标顶点
        (tl, tr, br, bl) = rect

        # 计算宽度
        width_a = np.sqrt(((br[0] - bl[0]) ** 2) + ((br[1] - bl[1]) ** 2))
        width_b = np.sqrt(((tr[0] - tl[0]) ** 2) + ((tr[1] - tl[1]) ** 2))
        max_width = max(int(width_a), int(width_b))

        # 计算高度
        height_a = np.sqrt(((tr[0] - br[0]) ** 2) + ((tr[1] - br[1]) ** 2))
        height_b = np.sqrt(((tl[0] - bl[0]) ** 2) + ((tl[1] - bl[1]) ** 2))
        max_height = max(int(height_a), int(height_b))

        # 目标点
        dst = np.array([
            [0, 0],
            [max_width - 1, 0],
            [max_width - 1, max_height - 1],
            [0, max_height - 1]
        ], dtype=np.float32)

        # 透视变换矩阵
        M = cv2.getPerspectiveTransform(rect, dst)

        # 应用变换
        warped = cv2.warpPerspective(self.original, M, (max_width, max_height))

        return warped

    def enhance_document(self, doc: np.ndarray) -> np.ndarray:
        """增强文档可读性"""
        # 转为灰度
        gray = cv2.cvtColor(doc, cv2.COLOR_BGR2GRAY)

        # 自适应阈值
        enhanced = cv2.adaptiveThreshold(
            gray, 255,
            cv2.ADAPTIVE_THRESH_GAUSSIAN_C,
            cv2.THRESH_BINARY,
            11, 2
        )

        # 转回BGR
        enhanced = cv2.cvtColor(enhanced, cv2.COLOR_GRAY2BGR)

        return enhanced

    def scan(self) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
        """执行完整扫描流程"""

        # 1. 预处理
        binary = self.preprocess()

        # 2. 找到文档边缘
        corners = self.find_document_contour(binary)

        if corners is None:
            print("警告:未检测到文档边缘")
            return self.original, self.original, self.original

        # 3. 透视变换
        scanned = self.perspective_transform(corners)

        # 4. 增强
        enhanced = self.enhance_document(scanned)

        return self.original, scanned, enhanced

    def visualize(self):
        """可视化扫描结果"""
        original, scanned, enhanced = self.scan()

        fig, axes = plt.subplots(1, 3, figsize=(18, 6))

        axes[0].imshow(cv2.cvtColor(original, cv2.COLOR_BGR2RGB))
        axes[0].set_title('原始图像')
        axes[0].axis('off')

        axes[1].imshow(cv2.cvtColor(scanned, cv2.COLOR_BGR2RGB))
        axes[1].set_title('透视校正后')
        axes[1].axis('off')

        axes[2].imshow(cv2.cvtColor(enhanced, cv2.COLOR_BGR2RGB))
        axes[2].set_title('增强处理后')
        axes[2].axis('off')

        plt.tight_layout()
        plt.show()

    def save(self, output_path: str):
        """保存扫描结果"""
        _, scanned, enhanced = self.scan()
        cv2.imwrite(output_path, enhanced)
        print(f"已保存到: {output_path}")


# 使用示例
if __name__ == "__main__":
    scanner = DocumentScanner('document.jpg')
    scanner.visualize()
    scanner.save('scanned_document.jpg')

9. 避坑小贴士

常见错误1:结构元素大小选择不当

现象:腐蚀/膨胀后图像失真严重

原因

  • 结构元素太大:会丢失重要的细节
  • 结构元素太小:无法有效去除噪声

正确做法

# 先检查噪声大小,选择合适的结构元素
noise_size = estimate_noise_size(gray_img)
kernel_size = max(3, noise_size * 2 + 1)  # 通常取噪声大小的2倍+1
kernel = np.ones((kernel_size, kernel_size), dtype=np.uint8)

常见错误2:Otsu阈值在光照不均时失效

现象:Otsu自动计算的阈值在光照不均匀的图像上效果很差

原因:Otsu假设前景和背景的灰度分布是双峰的,光照不均会破坏这个假设

解决方案

# 方法1:使用自适应阈值
adaptive_thresh = cv2.adaptiveThreshold(
    gray, 255,
    cv2.ADAPTIVE_THRESH_GAUSSIAN_C,
    cv2.THRESH_BINARY,
    blockSize=11,
    C=2
)

# 方法2:先做直方图均衡化
equalized = cv2.equalizeHist(gray)
_, otsu_after_eq = cv2.threshold(equalized, 0, 255, cv2.THRESH_OTSU)

# 方法3:分割图像为多个区域,分别阈值

常见错误3:分水岭过度分割

现象:分水岭将一个物体分割成多个区域

原因:噪声和细小的纹理被当作"山谷"

解决方案

# 方法1:使用标记图控制分水岭
# 先用形态学操作创建标记
opening = cv2.morphologyEx(thresh, cv2.MORPH_OPEN, kernel, iterations=2)

# 方法2:提高阈值以减少小的"山谷"
_, sure_bg = cv2.dilate(opening, kernel, iterations=3)

# 方法3:使用距离变换创建更精确的前景标记
dist_transform = cv2.distanceTransform(opening, cv2.DIST_L2, 5)
_, sure_fg = cv2.threshold(dist_transform, 0.3 * dist_transform.max(), 255, 0)

常见错误4:连通域分析忽略边界情况

现象:遗漏小目标或计数错误

原因:没有处理边界情况,如空图像、单像素目标等

正确做法

num_labels, labels, stats, centroids = cv2.connectedComponentsWithStats(binary_img)

# 必须检查
if num_labels <= 1:
    print("没有检测到任何前景")
    return []

for i in range(1, num_labels):
    area = stats[i, cv2.CC_STAT_AREA]
    if area < min_area_threshold:
        continue  # 过滤太小的目标
    # 处理每个有效目标

10. 本章小结

通过本章的学习,你应该已经掌握了:

  1. 形态学基础:理解了结构元素的概念和不同形状的特点
  2. 基本形态学运算:掌握了腐蚀、膨胀、开运算、闭运算的原理和应用
  3. 阈值分割:学会了固定阈值、Otsu自动阈值、自适应阈值的适用场景
  4. 连通域分析:能够使用连通域标记进行目标计数和属性提取
  5. 分水岭算法:理解了基于地形学的分水岭原理和标记图控制方法
  6. GrabCut算法:掌握了交互式前景提取技术
  7. 实战项目:完成了一个文档扫描仪的完整实现

一句话总结:形态学操作是图像分割的瑞士军刀,不同的组合可以解决不同的问题;而阈值分割则是最简单有效的快速分割工具,选择合适的方法取决于具体的应用场景。


11. 练习与思考

  1. 形态学设计:设计一个结构元素和运算组合,用于提取车牌字符
  2. 阈值对比实验:对比固定阈值、Otsu、自适应阈值在同一图像上的效果差异
  3. 连通域应用:使用连通域分析实现血细胞计数
  4. 分水岭改进:思考如何改进分水岭算法以减少过度分割
  5. 综合项目:实现一个名片扫描识别系统

下一章预告:第5章《特征提取与描述子》将带你学习SIFT、ORB、Harris等经典特征算法,这是图像匹配和目标识别的核心技术。


如果本章内容对你有帮助,欢迎点赞、收藏和关注。有任何问题可以在评论区留言。

Logo

腾讯云面向开发者汇聚海量精品云计算使用和开发经验,营造开放的云计算技术生态圈。

更多推荐