图像处理基础 第2课/共35课

第02课:滤波与卷积

📖 课程概述

本课学习滤波与卷积的核心原理与实现。我们将从数学基础出发,通过代码实践真正理解每一个算法的运作机制。

📑 本课目录

1. 卷积的数学原理2. 均值滤波3. 高斯滤波4. 中值滤波5. 双边滤波6. 自定义卷积核7. 滤波器选择指南

1. 卷积的数学原理

卷积(Convolution)是图像处理最核心的操作。它用一个小的核(kernel)在图像上滑动,每个位置计算核与对应像素的加权和。

(f * g)(x, y) = ∑st f(x-s, y-t) · g(s, t)
其中 f 是输入图像,g 是卷积核

🔑 卷积的本质理解

import cv2, numpy as np

def manual_conv2d(image, kernel, padding='zero'):
    kh, kw = kernel.shape
    ph, pw = kh // 2, kw // 2
    if padding == 'zero':
        padded = np.pad(image, ((ph,ph),(pw,pw)), mode='constant', constant_values=0)
    elif padding == 'reflect':
        padded = np.pad(image, ((ph,ph),(pw,pw)), mode='reflect')
    else:
        padded = np.pad(image, ((ph,ph),(pw,pw)), mode='edge')
    h, w = image.shape
    output = np.zeros_like(image, dtype=np.float64)
    for i in range(h):
        for j in range(w):
            output[i,j] = np.sum(padded[i:i+kh, j:j+kw] * kernel)
    return np.clip(output, 0, 255).astype(np.uint8)

img = np.zeros((64,64), dtype=np.uint8)
img[20:44, 20:44] = 200
noise = np.random.normal(0, 30, img.shape)
img_noisy = np.clip(img + noise, 0, 255).astype(np.uint8)

mean_kernel = np.ones((5,5)) / 25.0
mean_result = manual_conv2d(img_noisy, mean_kernel)
cv_mean = cv2.filter2D(img_noisy, -1, mean_kernel)
diff = np.mean(np.abs(mean_result.astype(float) - cv_mean.astype(float)))
print(f"手动 vs OpenCV 均值差: {diff:.4f}")
print(f"均值滤波MSE: {np.mean((mean_result.astype(float)-img.astype(float))**2):.2f}")
print(f"原图MSE: {np.mean((img_noisy.astype(float)-img.astype(float))**2):.2f}")
cv2.imwrite('/var/www/ttl/cv/l02_conv.png', np.hstack([img_noisy, mean_result]))

2. 均值滤波

均值滤波是最简单的线性滤波器,用邻域内所有像素的算术平均值替代中心像素。

h(x,y) = (1/mn) ∑(i,j)∈W f(i,j)
其中 W 是 m×n 的邻域窗口
import cv2, numpy as np
img = np.zeros((200,200), dtype=np.uint8)
cv2.rectangle(img, (50,50), (150,150), 200, -1)
cv2.circle(img, (100,100), 20, 255, -1)
noise = np.random.normal(0, 25, img.shape)
img_noisy = np.clip(img.astype(np.float64) + noise, 0, 255).astype(np.uint8)

for ksize in [3, 5, 9, 15]:
    blurred = cv2.blur(img_noisy, (ksize, ksize))
    mse = np.mean((blurred.astype(float) - img.astype(float))**2)
    print(f"核 {ksize}x{ksize}: MSE={mse:.2f}")
vis = [img_noisy]
for ksize in [3, 5, 9, 15]:
    vis.append(cv2.blur(img_noisy, (ksize, ksize)))
cv2.imwrite('/var/www/ttl/cv/l02_mean.png', np.hstack(vis))

3. 高斯滤波

高斯滤波是最常用的平滑滤波器,核权重按高斯分布分配。

G(x,y) = (1/2πσ²) · exp(-(x²+y²)/(2σ²))
σ 越大,平滑范围越广
import cv2, numpy as np
img = np.zeros((200,200), dtype=np.uint8)
cv2.rectangle(img, (40,40), (160,160), 180, -1)
noise_img = img.copy()
sp = np.random.random(img.shape)
noise_img[sp < 0.02] = 0; noise_img[sp > 0.98] = 255

for sigma in [0.5, 1.0, 2.0, 5.0]:
    blurred = cv2.GaussianBlur(noise_img, (0,0), sigmaX=sigma)
    mse = np.mean((blurred.astype(float) - img.astype(float))**2)
    print(f"sigma={sigma}: MSE={mse:.2f}")

kernel_3x3 = cv2.getGaussianKernel(3, 1.0)
gauss_2d = kernel_3x3 @ kernel_3x3.T
print("3x3高斯核归一化:"); print(np.round(gauss_2d/gauss_2d.sum(), 4))

gauss_r = cv2.GaussianBlur(noise_img, (5,5), 1.5)
mean_r = cv2.blur(noise_img, (5,5))
print(f"高斯MSE: {np.mean((gauss_r.astype(float)-img.astype(float))**2):.2f}")
print(f"均值MSE: {np.mean((mean_r.astype(float)-img.astype(float))**2):.2f}")
cv2.imwrite('/var/www/ttl/cv/l02_gauss.png', np.hstack([noise_img, gauss_r, mean_r]))

4. 中值滤波

中值滤波是非线性滤波器,对椒盐噪声效果极佳。

h(x,y) = median{ f(i,j) | (i,j) ∈ W }
import cv2, numpy as np
np.random.seed(42)
img = np.zeros((200,200), dtype=np.uint8)
cv2.rectangle(img, (30,30), (170,170), 180, -1)
cv2.putText(img, 'CV', (70,120), cv2.FONT_HERSHEY_SIMPLEX, 2, 255, 2)
sp_img = img.copy()
mask = np.random.random(img.shape)
sp_img[mask < 0.05] = 0; sp_img[mask > 0.95] = 255
actual_rate = np.sum(sp_img != img) / img.size
print(f"实际噪声比例: {actual_rate:.2%}")

median_r = cv2.medianBlur(sp_img, 3)
gauss_r = cv2.GaussianBlur(sp_img, (5,5), 1.5)
mean_r = cv2.blur(sp_img, (5,5))
for name, result in [('中值', median_r), ('高斯', gauss_r), ('均值', mean_r)]:
    mse = np.mean((result.astype(float) - img.astype(float))**2)
    print(f"{name}滤波 MSE: {mse:.2f}")
cv2.imwrite('/var/www/ttl/cv/l02_median.png', np.hstack([sp_img, median_r, gauss_r]))

5. 双边滤波

双边滤波同时考虑空间距离像素值差异,保边去噪。

BF(I) = (1/Wp) ∑ Gσs(‖p-q‖) · Gσr(|Ip-Iq|) · Iq
import cv2, numpy as np
img = np.zeros((200,200), dtype=np.uint8)
img[:,:100] = 50; img[:,100:] = 200
noise = np.random.normal(0, 20, img.shape)
noisy = np.clip(img.astype(np.float64) + noise, 0, 255).astype(np.uint8)
bilateral = cv2.bilateralFilter(noisy, d=9, sigmaColor=75, sigmaSpace=75)
gauss = cv2.GaussianBlur(noisy, (9,9), 2.0)
edge_b = np.abs(np.diff(bilateral[:,100].astype(float))).sum()
edge_g = np.abs(np.diff(gauss[:,100].astype(float))).sum()
print(f"边缘梯度 - 双边: {edge_b:.0f}, 高斯: {edge_g:.0f}")
var_b = np.var(bilateral[50:150,20:80].astype(float))
var_g = np.var(gauss[50:150,20:80].astype(float))
print(f"平坦方差 - 双边: {var_b:.1f}, 高斯: {var_g:.1f}")
cv2.imwrite('/var/www/ttl/cv/l02_bilateral.png', np.hstack([noisy, bilateral, gauss]))

6. 自定义卷积核

import cv2, numpy as np
img = np.zeros((200,200), dtype=np.uint8)
cv2.rectangle(img, (50,50), (150,150), 200, -1)
cv2.circle(img, (100,100), 25, 255, -1)
kernels = {
    '浮雕': np.array([[-2,-1,0],[-1,1,1],[0,1,2]]),
    '锐化': np.array([[0,-1,0],[-1,5,-1],[0,-1,0]]),
    '拉普拉斯': np.array([[0,1,0],[1,-4,1],[0,1,0]]),
    'Sobel_X': np.array([[-1,0,1],[-2,0,2],[-1,0,1]]),
    'Sobel_Y': np.array([[-1,-2,-1],[0,0,0],[1,2,1]]),
}
for name, kernel in kernels.items():
    filtered = cv2.filter2D(img, -1, kernel)
    diff = np.mean(np.abs(filtered.astype(float) - img.astype(float)))
    print(f"{name}: 平均像素变化={diff:.2f}")

7. 滤波器选择指南

噪声类型最佳选择次优选择避免
高斯噪声高斯滤波均值滤波中值滤波
椒盐噪声中值滤波双边滤波均值滤波
需保边去噪双边滤波中值滤波高斯
需增强边缘Unsharp Mask拉普拉斯均值滤波
实战建议:先用中值滤波处理椒盐噪声,再用高斯/双边做进一步平滑。组合使用效果往往优于单一滤波器。

🔑 滤波与卷积核心要点

📝 课后练习

  1. 实现本课所有算法并对比效果
  2. 在不同参数下测试算法的鲁棒性
  3. 将本课方法应用到自己的数据集
  4. 分析算法的计算复杂度和优化方向
延伸阅读:推荐阅读《Digital Image Processing》(Gonzalez)相关章节,以及OpenCV官方文档中的详细API说明。实际项目中,建议先在简单数据上验证算法,再迁移到复杂场景。

📊 方法对比总结

方法优点缺点适用场景
方法A简单高效精度有限快速原型
方法B精度高计算量大离线处理
方法C平衡精度与速度参数调优复杂实际应用

📐 关键公式速查

本课涉及的核心数学公式汇总,方便快速参考:

🔬 深入理解:卷积的频域视角

从频域角度看,不同卷积核对应不同的频率响应:

均值核: 低通滤波器,截止频率与核大小成反比
高斯核: 最优低通滤波器(时频不确定性最小)
拉普拉斯核: 高通滤波器,检测高频变化
Sobel核: 带通滤波器,选择特定方向的高频

卷积的可分离性

2D高斯核可分离为两个1D核:G(x,y) = g(x)·g(y),计算复杂度从O(N·K²)降到O(N·2K)。

import cv2, numpy as np, time

# 可分离卷积加速
img = np.random.randint(0, 256, (512, 512), dtype=np.uint8)
kernel_2d = cv2.getGaussianKernel(15, 3.0)
kernel_2d = kernel_2d @ kernel_2d.T

t0 = time.time()
for _ in range(100):
    result_2d = cv2.filter2D(img, -1, kernel_2d)
t1 = time.time() - t0

t0 = time.time()
for _ in range(100):
    result_sep = cv2.GaussianBlur(img, (15, 15), 3.0)
t2 = time.time() - t0

print(f"2D卷积: {t1*10:.1f}ms, 可分离: {t2*10:.1f}ms")
print(f"结果差异: {np.mean(np.abs(result_2d.astype(float)-result_sep.astype(float))):.4f}")
print(f"加速比: {t1/t2:.1f}x")

🔬 深入理解:滤波器设计原则

设计滤波器时需要考虑以下因素:

因素低通滤波高通滤波带通滤波
目的去噪/平滑边缘检测纹理分析
核权重和=1(保均值)=0(去除DC)取决于设计
对称性中心对称反对称取决于设计
频域响应通低频阻高频通高频阻低频通特定频段

📈 自适应滤波

标准滤波器对整幅图像使用相同参数,自适应滤波根据局部特征调整:

import cv2, numpy as np

# 非局部均值去噪
img = np.zeros((200,200), dtype=np.uint8)
cv2.rectangle(img, (50,50), (150,150), 200, -1)
noise = np.random.normal(0, 25, img.shape)
noisy = np.clip(img + noise, 0, 255).astype(np.uint8)

# NLM去噪
nlm = cv2.fastNlMeansDenoising(noisy, None, h=10, templateWindowSize=7, searchWindowSize=21)
gauss = cv2.GaussianBlur(noisy, (5,5), 1.5)

nlm_mse = np.mean((nlm.astype(float)-img.astype(float))**2)
gauss_mse = np.mean((gauss.astype(float)-img.astype(float))**2)
print(f"NLM MSE: {nlm_mse:.2f}, 高斯 MSE: {gauss_mse:.2f}")
print("NLM保边性优于高斯,但计算慢100x+")

💻 滤波完整Pipeline代码

以下是一个完整的滤波处理管道,从数据准备到结果评估,包含所有关键步骤和参数调优建议:

import cv2
import numpy as np
import time

class 滤波Pipeline:
    """完整的滤波处理管道"""
    
    def __init__(self, params=None):
        self.params = params or {}
        self.results = {}
    
    def preprocess(self, image):
        """预处理步骤"""
        # 1. 去噪
        if self.params.get('denoise', True):
            image = cv2.GaussianBlur(image, (3, 3), 1.0)
        
        # 2. 对比度增强
        if self.params.get('enhance', False):
            if len(image.shape) == 2:
                image = cv2.equalizeHist(image)
            else:
                lab = cv2.cvtColor(image, cv2.COLOR_BGR2LAB)
                lab[:,:,0] = cv2.equalizeHist(lab[:,:,0])
                image = cv2.cvtColor(lab, cv2.COLOR_LAB2BGR)
        
        # 3. 尺寸调整
        max_size = self.params.get('max_size', 1024)
        h, w = image.shape[:2]
        if max(h, w) > max_size:
            scale = max_size / max(h, w)
            image = cv2.resize(image, (int(w*scale), int(h*scale)))
        
        return image
    
    def process(self, image):
        """核心处理步骤"""
        t0 = time.time()
        preprocessed = self.preprocess(image)
        t1 = time.time()
        
        # 主处理逻辑
        result = self._main_process(preprocessed)
        t2 = time.time()
        
        self.results['timing'] = {
            'preprocess': (t1-t0)*1000,
            'process': (t2-t1)*1000,
            'total': (t2-t0)*1000
        }
        return result
    
    def _main_process(self, image):
        """主处理逻辑(子类可重写)"""
        return image
    
    def evaluate(self, result, ground_truth=None):
        """评估处理结果"""
        metrics = {}
        if ground_truth is not None:
            if len(result.shape) == 2:
                mse = np.mean((result.astype(float) - ground_truth.astype(float))**2)
                metrics['mse'] = mse
                metrics['psnr'] = 10 * np.log10(255**2 / (mse + 1e-10))
        metrics['timing'] = self.results.get('timing', {})
        return metrics
    
    def visualize(self, image, result):
        """可视化对比"""
        if len(image.shape) == 2:
            image = cv2.cvtColor(image, cv2.COLOR_GRAY2BGR)
        if len(result.shape) == 2:
            result = cv2.cvtColor(result, cv2.COLOR_GRAY2BGR)
        h1, w1 = image.shape[:2]
        h2, w2 = result.shape[:2]
        h, w = max(h1,h2), w1+w2+10
        vis = np.zeros((h, w, 3), dtype=np.uint8)
        vis[:h1, :w1] = image
        vis[:h2, w1+10:] = result
        cv2.putText(vis, 'Input', (10, 25), cv2.FONT_HERSHEY_SIMPLEX, 0.7, (0,255,0), 2)
        cv2.putText(vis, 'Output', (w1+20, 25), cv2.FONT_HERSHEY_SIMPLEX, 0.7, (0,255,0), 2)
        return vis


# 使用示例
pipeline = 滤波Pipeline({'denoise': True, 'enhance': False, 'max_size': 512})
test_img = np.random.randint(0, 256, (256, 256, 3), dtype=np.uint8)
result = pipeline.process(test_img)
metrics = pipeline.evaluate(result)
print(f"处理时间: {metrics['timing'].get('total', 0):.1f}ms")

🔧 参数调优指南

不同参数对结果的影响:

建议从默认参数开始,根据结果逐步调整。先在少量数据上快速迭代,确定参数后再全量处理。

生产环境建议:1) 添加输入验证和异常处理 2) 记录每次处理的参数和结果 3) 设置性能监控和告警 4) 定期在测试集上评估质量 5) 保持模型/算法版本管理

📊 完整评估指标

import numpy as np

class Metrics:
    @staticmethod
    def mse(img1, img2):
        return np.mean((img1.astype(float) - img2.astype(float))**2)
    
    @staticmethod
    def psnr(img1, img2, max_val=255):
        mse_val = Metrics.mse(img1, img2)
        if mse_val == 0: return float('inf')
        return 10 * np.log10(max_val**2 / mse_val)
    
    @staticmethod
    def iou(mask1, mask2):
        inter = np.logical_and(mask1 > 0, mask2 > 0).sum()
        union = np.logical_or(mask1 > 0, mask2 > 0).sum()
        return inter / (union + 1e-10)
    
    @staticmethod
    def f1_score(pred, gt, threshold=0.5):
        pred_bin = (pred > threshold)
        gt_bin = (gt > 0)
        tp = np.logical_and(pred_bin, gt_bin).sum()
        precision = tp / (pred_bin.sum() + 1e-10)
        recall = tp / (gt_bin.sum() + 1e-10)
        return 2 * precision * recall / (precision + recall + 1e-10)

# 示例
img1 = np.random.randint(0, 256, (100,100), dtype=np.uint8)
img2 = img1 + np.random.randint(-10, 10, (100,100))
img2 = np.clip(img2, 0, 255).astype(np.uint8)

print(f"MSE: {Metrics.mse(img1,img2):.2f}")
print(f"PSNR: {Metrics.psnr(img1,img2):.2f} dB")
print(f"说明: PSNR>30dB为良好, >40dB为优秀")

✅ 实机验证

🏆

卷积大师

你已经掌握了滤波与卷积的核心知识,继续前进!