做图像恢复的人基本都遇到过这种尴尬：一张看起来只是有点糊的图，直接在频域里除以传递函数，结果不但没变清晰，反而被密密麻麻的噪点淹没，连轮廓都认不出来。这个问题的根源，在于我们把噪声也一并做了除法。今天我以一套完整的Matlab实现为例，把逆滤波器、维纳滤波器这两种经典图像恢复方法讲透，不仅给代码，更把背后的退化模型、参数取舍和踩坑点说清楚。适合刚接触图像处理的学生，也适合需要实际解决运动模糊、散焦模糊图片问题的工程师。
1. 图像退化模型与恢复思路
图像恢复的第一步不是写代码，而是建立退化模型。只有把“图像是怎么变糊的”用数学写清楚，才知道恢复算法该在哪个环节入手。
1.1 退化模型：从空域卷积到频域乘积
一张清晰的图像在成像过程中，通常被两个因素破坏：一是点扩散函数造成的模糊，二是传感器或信道引入的噪声。点扩散函数就是光学系统对单个理想点源的响应，它决定了图像被“抹开”的程度。退化过程可以写成：
g(m,n) = f(m,n) * h(m,n) + n(m,n)
其中g是退化图像，f是原始清晰图像，h是点扩散函数（PSF），*代表卷积，n是加性噪声。
空域卷积算起来很麻烦，但转到频域就变成乘法：
G(u,v) = F(u,v) · H(u,v) + N(u,v)
这里H(u,v)就是 PSF 的傅里叶变换，也叫光学传递函数（OTF）。图像恢复本质上就是反卷积：已知G和H，尽可能估计出F。如果只看到时域卷积，很难理解为什么逆滤波会失败；一旦到了频域，“除以H”这个问题就暴露得很清楚。
我在实验中常用运动模糊来当例子，因为运动模糊的传递函数存在大量零点，最能体现逆滤波的病态性。运动模糊PSF用fspecial('motion', len, theta)直接生成，len是模糊像素长度，theta是运动方向角度。散焦模糊则常用高斯PSF模拟。无论哪种模糊，只要把PSF补零到与图像相同尺寸，再fft2得到H，就可以开始后续处理。
1.2 为什么直接逆滤波会失败：噪声放大与零点陷阱
如果忽略噪声，逆滤波是非常直接的思路：在频域把退化图像乘以1/H(u,v)，反变换回去就是原图的估计。但一旦加入噪声，事情就变了。对退化模型做一个除法：
F_hat = F + N / H
理想情况下，只要H不等于零，噪声项N/H被除以一个值后可能变大。尤其是当H的模值很小、接近零时，噪声那一项的放大倍数会极其恐怖。这就是为什么直接逆滤波恢复出来的图像常常全是高频噪点，而不是清晰图像。
打个比方：逆滤波像拿着一个放大镜看一张被涂花的照片，涂花的痕迹正好是那些高频模糊。如果你强行把模糊痕迹“除干净”，涂花过程中的灰尘颗粒也被同步放大了，结果照片变成了颗粒噪点。
更麻烦的是，理想运动模糊的 OTF 存在真正的零点。在那些频点上，无论你是模糊还是清晰，信息都已经丢失了，1/H直接无穷大。物理上，这些频点上的原图信息被光学系统完全滤除，无法恢复。所以逆滤波必须对H做约束，通常的做法是设定阈值，只保留H模值比较大的频点。
2. 逆滤波器：原理与实现
2.1 逆滤波的核心公式与边界条件
逆滤波的核心公式很简单：
F_hat(u,v) = G(u,v) / H(u,v)
实际操作时，我们要处理的不只是公式本身，还要处理分母为零的情况。严格说，数学上只要H(u,v) = 0，这个频点的F_hat就无法定义。但工程上，接近零的H也会有灾难性噪声放大，所以只靠避开严格零点不够。
常用做法是设定一个相对阈值。比如：
threshold = 0.01 × max(|H|)
只有|H(u,v)| > threshold的频点才参与除法，其余频点直接置零或保留为低频信息。这个阈值很关键：太小会带入噪声，太大又会导致图像变糊。我用过的经验值是0.01到0.05之间，具体要看噪声水平。
还需要注意边界条件。图像卷积时，如果PSF尺寸和图像尺寸不一致，直接fft2会产生尺寸不匹配，结果错位。这里用psf2otf最省心，它会自动把PSF放到正确位置并补零到指定尺寸。手动补零也可以，但要记住把PSF中心移动到左上角，否则相位会出错。
2.2 Matlab实现逆滤波的完整流程
下面是一段我实际在项目中用过的逆滤波函数。输入是退化图像和已知PSF，输出是恢复图像。
function img_inv = inverse_filter(blurred_img, PSF, threshold) % blurred_img: 退化图像，double 类型，范围 [0,1] % PSF: 点扩散函数 % threshold: 相对阈值，小于该阈值的频点不参与除法 if nargin < 3 threshold = 0.01; end % 1. 由 PSF 得到与图像同尺寸的 OTF H = psf2otf(PSF, size(blurred_img)); % 2. 构造逆滤波器 H_abs = abs(H); inv_filter = zeros(size(H)); valid = H_abs > threshold * max(H_abs(:)); inv_filter(valid) = 1 ./ H(valid); % 3. 频域滤波并反变换 F = fft2(blurred_img); img_inv = real(ifft2(F .* inv_filter)); end这段代码有三个细节值得多说几句。
第一，psf2otf是处理PSF尺寸的标准函数，它做两件事：补零到目标尺寸、按卷积的要求循环移位。如果自己补零，很容易忽略移位导致恢复图像整体偏移。第二，valid掩码不是单纯比的零，而是相对最大模值的百分比。这样不同PSF都能用一个统一阈值，不用针对每张图重新调绝对数值。第三，反变换后取real而不是abs。real保留正负号，因为图像在空域本来就有负值；abs会把所有负值翻正，导致亮度结构错乱。我在早期调试时在这里踩过坑，显示出来图像有明显“晕影”，后来才确认是取错幅值。
2.3 截断半径与振铃现象的取舍
如果PSF的OTF只有一个中心大包，高频区域快速衰减，那么一个固定半径的低通掩码也能凑合工作。但运动模糊的OTF沿方向有周期性零点，不能用简单圆形半径截断，只能逐个频点做阈值判断。
截断阈值选得太小，恢复图会被噪点淹没；选得太大，图像低频保留了，但高频细节也丢了，整体看起来像蒙了一层纱。更常见的副作用是“振铃”，也就是图像边缘附近出现一圈圈明暗交替的光晕。振铃的根源是频域截断和循环卷积假设，相当于在边缘产生了类似吉布斯效应的现象。
抑制振铃的几个土办法：一是恢复前对图像做边缘延拓或edgetaper；二是反变换前在频域加一个平滑窗，例如汉明窗，把高频硬截断变成缓过渡；三是在恢复结果上做中值滤波，去散粒噪声。第一种最有效，但本文重点是基础方法，振铃的完整处理放到后面章节展开。
3. 维纳滤波器：从理论到Matlab实现
维纳滤波是经典的最小均方误差准则下的最优线性滤波器。从工程角度看，它就是在逆滤波器的分母里加了一个“噪声保护项”，让恢复过程不至于把噪声放大到失控。
3.1 最小均方误差思路：恢复误差的期望最小化
维纳滤波的思想是找这样一个滤波器W(u,v)，使得恢复图像f_hat = w * g与真实图像f的均方误差最小。数学上写成：
E[ |f(m,n) - f_hat(m,n)|^2 ] → min
对频域求解，得到维纳滤波器的标准形式：
W(u,v) = H*(u,v) / ( |H(u,v)|^2 + S_n(u,v) / S_f(u,v) )
其中H*(u,v)是H的共轭，S_n是噪声功率谱，S_f是原图像功率谱。恢复时用F_hat = W · G。
对比逆滤波1/H，维纳滤波多了一个分母项S_n / S_f。在信号能量强、噪声能量弱的频段，这个比值小，维纳滤波接近逆滤波，尽量保留细节；在噪声占主导的频段，比值大，抑制放大。
可以把维纳滤波理解成一个“自适应均衡器”：它不是在每个频点上都盲目除以H，而是根据信噪比决定放大多少。信噪比高的频点多恢复一些，信噪比低的频点少恢复一些，这样整体误差最小。
3.2 常数近似：不知道原图功率谱怎么办
实际使用中最尴尬的问题是，S_f是未知的。我们手里只有退化图像，没有原图，无法精确计算信号功率谱。所以工程上常做简化，用一个常数K代替S_n / S_f：
W(u,v) = H*(u,v) / ( |H(u,v)|^2 + K )
这个K虽然不是严格最优，但实践证明只要量级合适，恢复效果已经很稳定。K的物理含义是噪声和信号的全局功率比，它控制的是“去卷积强度”和“噪声抑制强度”之间的天平。
如果噪声比较小，K可以取 0.001 到 0.005；如果噪声明显，则取 0.01 到 0.1。还有一种更客观的估计方式：先在退化图像中找一块平坦区域，计算局部方差作为噪声方差的近似sigma_n^2；再用全图方差近似信号方差sigma_f^2，于是K = sigma_n^2 / sigma_f^2。实测下来，这个估计值至少能给出一个合理起点，再手动微调一到两个数量级就够了。
3.3 Matlab实现维纳滤波的步骤与代码
维纳滤波的Matlab实现几乎是逆滤波的翻版，区别只在滤波器构造。下面是手动实现版本：
function img_wiener = wiener_filter(blurred_img, PSF, K) % K: 噪信功率比，常数或与图像同尺寸的数组均可 H = psf2otf(PSF, size(blurred_img)); wiener_filter = conj(H) ./ (abs(H).^2 + K); F = fft2(blurred_img); img_wiener = real(ifft2(F .* wiener_filter)); end这个版本本质是把逆滤波器的1/H换成了更温和的表达式。当K=0时，它退化为逆滤波；当K很大时，它趋近于conj(H)/K，表现为对退化图像做一个高通锐化，但不会放大噪声。
如果不想手写原理，Matlab也有现成函数deconvwnr(blurred_img, PSF, nsr)。注意第三个参数不是K，而是噪信功率比NSR。它可以是常数，也可以是和图像同尺寸的数组。如果已知噪声方差，可以传K = noise_var / signal_var。
我用这两种方式对比过，结果几乎一样。手写版本的好处是你可以把滤波器矩阵提取出来，分析不同频点被放大或抑制的程度；内置函数的好处是省事，适合快速验证。对学生来说，我建议先把手动版本跑通，再去用工具箱函数，不然出了问题很难定位。
4. 完整实验：生成退化图与两种恢复对比
接下来我把完整的实验流程走一遍，从测试图像生成到退化模拟，再分别用逆滤波和维纳滤波恢复。这部分内容可以直接抄到Matlab里运行，便于观察两种方法的差别。
4.1 测试图像与退化模拟脚本
为了避免依赖特定数据集，也为了演示效果更直观，我用合成图像做测试。一个包含方块、圆环和斜条的图案，频谱特征丰富，最适合观察恢复效果。
% 生成测试图像 img = zeros(256, 256); % 左上方块 img(60:120, 60:120) = 1; % 右下深灰方块 img(140:200, 140:200) = 0.7; % 同心圆环 [X, Y] = meshgrid(1:256, 1:256); dist = sqrt((X - 128).^2 + (Y - 128).^2); ring = (dist > 40) & (dist < 60); img(ring) = 0.4; % 斜线 for k = 0:10 img(100 + k * 2, 50 + k * 2 : 206 + k * 2) = 0.8; end % 显示 figure('Name', '原始图像'); imshow(img, []); % 退化：运动模糊 + 高斯噪声 PSF = fspecial('motion', 15, 30); H = psf2otf(PSF, size(img)); blurred_clean = real(ifft2(H .* fft2(img))); noise_std = 0.02; blurred_noisy = blurred_clean + noise_std * randn(size(img)); blurred_noisy = max(0, min(1, blurred_noisy)); figure('Name', '退化图像'); imshow(blurred_noisy, []);合成图像的好处是原始干净图像img完全可控，恢复质量好坏一眼就能看出来。fspecial('motion', 15, 30)生成一个长度15像素、方向30度的运动模糊PSF。加噪时，我用标准差0.02，范围在0-1的图像上算中度噪声，肉眼能看到明显颗粒，但不至于完全淹没结构。
4.2 逆滤波恢复与截断半径选择
对退化图直接做逆滤波，代码如下：
threshold = 0.01; img_inv = inverse_filter(blurred_noisy, PSF, threshold); figure('Name', '逆滤波恢复'); imshow(img_inv, []);如果threshold=0.01，恢复图像会是什么样？大概率是满屏噪点，但方块和圆环的大致轮廓还能隐约看到。为什么噪点这么明显？因为运动模糊OTF沿垂直运动方向有一些模值很小的频点，虽然被阈值滤掉了一部分，但剩下的频点中噪声仍然被放大了上百倍，高频噪点在视觉上主导了结果。
把阈值加大到0.05，噪声会明显减少，但图像看起来也更糊。原因是那些模值中等、信息量不小的频点也被一并丢掉了。逆滤波就是这样：要么忍着噪声，要么忍着模糊，很难两全。
如果加噪很小，比如noise_std=0.001，逆滤波效果会很惊艳，边缘能恢复到接近原始状态。这也是很多人初学时会困惑“逆滤波明明效果很好”的原因——那是在噪声几乎为零的理想环境下测的。真实场景中，噪声是躲不掉的。
4.3 维纳滤波恢复与K值调整
维纳滤波在这个退化场景下要稳得多：
K = 0.02; img_wiener = wiener_filter(blurred_noisy, PSF, K); figure('Name', '维纳滤波恢复'); imshow(img_wiener, []);K=0.02时，恢复结果仍然是平滑的，不会有密集噪点，方块边缘有轻微振铃，但整体轮廓清晰。相比于逆滤波，维纳滤波更像“带约束的锐化”，它在放大边缘的同时压制了噪声。
如果把K改成0.001，恢复边缘更锐利，但噪声跟着起来；改成0.1，图像非常干净，但边缘发肉，有点像被再次模糊过。这个趋势符合理论预期。实际调参时，我会观察两个指标：一是背景平缓区域是否出现细小噪点，二是方块边缘是否还是锐利的一字形。只要背景开始出现明显颗粒，说明K太小；只要边缘变成模糊斜坡，说明K太大。
4.4 恢复效果对比：什么时候用哪个
为了直观，我用对比表格总结一下两种方法的定位。
| 方法 | 抗噪声能力 | 边缘恢复能力 | 参数难度 | 适用场景 |
|---|---|---|---|---|
| 逆滤波 | 很差 | 很强 | 阈值需要小心 | 噪声极小、PSF零点分布简单的理想实验 |
| 维纳滤波 | 较好 | 中等偏强 | 一个K值即可调 | 实际图像、有噪声、PSF已知且稳定 |
逆滤波最大的价值在于告诉我们“病态反卷积”到底有多敏感。维纳滤波则用一个小小的正则项，把病态问题变成了良态问题。如果你的PSF不能准确估计，维纳滤波因为对H的除法不敏感，反而更不容易崩。很多工程场景中，PSF是从标定板估出来的，本身就带误差，这时候逆滤波几乎必然失败，维纳滤波则还有一定容错能力。
5. 常见问题与调参技巧实录
这一部分是我实际调试过程中反复遇到的问题，整理成速查表，方便你直接对照解决。
5.1 逆滤波恢复图全是噪点，轮廓都看不到
这是最典型的问题。原因多半是阈值设得太小，噪声被过度放大。可以先看H的模值分布，再用max(abs(H(:)))估算一下量级，把阈值调到最大模值的1%到5%。如果调到5%还是不行，说明噪声本身太大，不要继续用逆滤波了，赶紧切到维纳滤波。
另外检查PSF是否和退化过程匹配。我用过一个标定错误的PSF，逆滤波恢复出的图像全是条纹，当时还以为是阈值问题，后来发现是PSF角度偏了几度。PSF角度一旦不对，OTF相位就完全乱掉了，恢复等于把错误的反向模糊施加在图像上，效果自然惨烈。
5.2 维纳滤波图像发糊，K值怎么调
图像发糊说明K偏大。把K减小一个量级，比如从0.05变成0.005，通常边缘会明显锐化。但如果你看到背景出现噪声，说明降得太多了，回到两者之间取中间值。
我习惯先用“平坦区域估噪”的方法算基准K：
% 选一块平坦区域估计噪声标准差 patch = blurred_noisy(10:30, 10:30); noise_std_est = std(patch(:)); signal_var_est = var(blurred_noisy(:)) - noise_std_est^2; K = noise_std_est^2 / signal_var_est;这个K只能当初始值，最终还要根据恢复效果微调。通常最终值会在初始值附近一个量级内浮动。如果你的图像内容非常复杂，比如有密集纹理，信号方差会被高估，K偏小；如果图像大片都是平坦天空，信号方差被低估，K偏大。这时候手动微调是不可避免的。
5.3 图像边缘出现周期性条纹
周期性条纹几乎都是循环卷积假设导致的边界振铃。频域卷积默认图像是周期延拓的，但实际图像左右边缘不连续，反变换后会在边缘产生振铃。这种振铃不是滤波器的错，是模型边界假设的问题。
缓解办法有几种。最省事的是恢复前用edgetaper把图像边缘过渡平滑掉。也可以把图像先延拓一圈，恢复后再裁剪掉。如果只是轻微振铃，可以在频域对恢复滤波器乘以一个汉明窗，把高频段的硬截断变成软过渡，振铃会轻很多。
5.4 恢复结果有负值，显示不正常
ifft2的结果虽然是实数，但会包含负值。直接imshow(img, [])会自动把负值映射到黑色，这没问题。但如果你用uint8(img)转换，所有负值都会被截断成0，恢复图像会丢失一部分亮度层次，看起来暗区细节全没了。
正确做法是先归一化再显示：
img_display = mat2gray(img_inv); imshow(img_display);mat2gray会自动把最小值映射为0，最大值映射为1，这样显示效果和人眼感知更接近。如果在写文件时想保留动态范围，可以先做线性拉伸，再转uint8。
5.5 一点个人调参习惯
我实际操作中不会对着公式硬调参数。我会把K按对数间隔取一组值：0.0001、0.001、0.01、0.1，分别恢复，再把结果凑到一张图上对比。噪声大不大、边缘糊不糊，一目了然。很多时候，不是调一个“最优参数”就够了，而是在清晰度、噪声、振铃三者之间找一个最适合当下使用场景的平衡点。
如果你的设备是固定机位，比如工业相机一直拍同一产线，PSF基本不变，那你可以把PSF和K标定一次，存成配置文件，下次直接调用。这样既快又稳，比每次现场调参靠谱得多。逆滤波在这种固定场景下也够用，但一旦光照变化导致噪点增加，还是维纳滤波更省心。
