ARTICLE DETAIL

资讯详情

深耕网站视觉设计与运营推广的一线实战洞察。

图像垂直条纹去除实战:频域陷波与列回归校正指南

图像垂直条纹去除实战:频域陷波与列回归校正指南 简介全局法图像垂直条纹去除是一份面向遥感图像处理、计算机视觉方向开发者的Matlab实用程序用于对高光谱图像或普通RGB图像中出现的规则垂直条纹噪声进行全局建模与消除。该思路适用于成像传感器像元响应不一致导致的条带伪影还针对倾斜条纹场景提示先旋转再处理并提供了无循环高速版与循环低速版两套实现经高光谱数据对比高速版运行效率可提升约20倍。资源压缩包共2个文件均为Matlab脚本.m体积仅2KB结构精简注释详尽便于直接运行、核对算法流程也可作为理解全局条纹去除原理的入门范例。当前已有886人学习下载适合希望快速验证算法效果、对比不同实现版本或在此基础上二次开发的工程人员与研究人员。1. 垂直条纹为什么难缠空域滤波失效背后的全局结构问题做图像算法的人几乎都在某个时刻被「图像垂直条纹」折磨过红外热像仪画面上一道道竖向亮暗纹工业线扫相机采集的图每隔几列就有一条遥感影像里经常整景都是条带噪声。这类噪声最烦人的地方在于你拿中值滤波、高斯滤波、双边滤波去试要么滤不干净要么把边缘和纹理一起抹花了属于典型的「越修越糟」。原因在于垂直条纹是覆盖整幅图像的周期结构而空域滤波只能看到局部窗口——窗口小了滤不掉窗口大了伤细节这是结构性矛盾。处理这类问题必须换到全局视角把图像变换到频率域在频谱上把条纹对应的能量精准剔除或者在全图范围内拟合列方向的全局趋势做修正。这份资源就是围绕「全局法」整理的一套可执行方案包含频域陷波滤波与列回归校正两条完整路线以及我踩过的各种坑。适合做红外图像处理、工业视觉检测、遥感影像预处理和低剂量 CT 重建的工程师也适合刚接触傅里叶变换但想直接把它落地到真实图像的初学者。2. 频域视角看条纹傅里叶变换如何定位周期噪声陷波器又该怎么设计2.1 垂直条纹在频谱上留下的记号一眼认出噪声峰要理解全局法先得知道垂直条纹在频域里长什么样。这里直接说结论一幅 M×N 的图像 I(x,y) 上叠加了垂直条纹条纹沿 y 方向延伸、灰度沿 x 方向呈周期变化可以写成 n(x) A·sin(2πu₀x φ) 的形式。对它做二维傅里叶变换后这个正弦分量的能量会集中到频域平面上的两个点 (±u₀, 0) 附近——注意第二个坐标是 0所以噪声峰全部落在频谱图的水平中轴线上也就是 v0 的那一行。在机器视觉和工业相机采集场景里图像傅里叶变换是分析周期噪声的第一选择原因就在于此空间域里横跨整幅图的条纹在频域里只是两个孤立的峰处理起来成本极低。用 fftshift 把零频移到中心后这两个峰就对称出现在原点左右两侧的水平轴线上。峰位到原点的距离和条纹周期直接相关图像宽度为 W 像素条纹周期为 T 像素那么峰位距原点的距离大约是 W/T 个频域像素。举个例子512 像素宽的图像上有周期为 16 像素的条纹峰位就在距原点 32 个像素的位置。实际图像往往不只一组条纹。传感器读出噪声、电源纹波、扫描机构抖动可能叠加出多个周期分量频谱水平轴线上会依次排开多个亮峰二次谐波、三次谐波也可能出现。这些峰的高度通常远高于周围频谱成分的平均水平这就是后面自动诊断算法的依据。需要特别说明的是空域滤波处理不了这类噪声的根本原因也在这里条纹是全域周期信号能量高度集中在特定频率你在空间域无论选多大窗口都只是在跟一个「无处不在」的信号较劲而频域只需动几个点。2.2 陷波器的设计为什么不能只置零一个像素明确了噪声峰的位置后常规做法是在频谱上做一个带阻滤波器把峰及其邻域的能量压掉叫陷波滤波Notch Filter。但这里有一个新手最容易犯的错以为把峰值那个像素置零就行了。实际做出来你会发现条纹确实淡了一点但整幅图像出现一圈圈水波纹一样的振铃而且条纹并没有彻底干净。原因在于离散傅里叶变换里周期噪声的能量并不是理想地集中在一个点上而是分布在以峰位为中心的一个小邻域内形状近似一个窄峰。只置零单个像素相当于用一根极细的针去戳一个面能量从旁边漏过去同时这个锐利的置零操作在反变换时又引入了新的高频振荡也就是振铃。正确做法是用一个平滑过渡的掩模把峰位周围的一小片区域压下去。我一般用高斯型陷波它的传递函数形式是H(u,v) 1 - exp(-(D²(u,v) - D₀²)² / (2·W²·D₀²))其中 D(u,v) 是当前像素到原点的距离D₀ 是噪声峰到原点的距离W 控制陷波的带宽。这个公式的好处是掩模从 1 到 0 是渐变过渡的不会在反变换时引入新的振铃。下面这段代码演示如何生成一个这样的陷波掩模import numpy as np import cv2 def gaussian_notch_mask(shape, center, radius2, bandwidth5): 生成高斯型陷波掩模 shape: (H, W) 频谱尺寸 center: (u0, v0) 噪声峰位置 radius: 需要压制的邻域半径 bandwidth: 陷波过渡带宽 H, W shape u np.arange(H) - center[1] # 行方向相对中心 v np.arange(W) - center[0] # 列方向相对中心 U, V np.meshgrid(v, u) D np.sqrt(U**2 V**2) D0 np.sqrt((center[0] - W//2)**2 (center[1] - H//2)**2) # 高斯带阻1 表示保留, 0 表示滤除 mask 1 - np.exp(-((D**2 - D0**2)**2) / (2 * bandwidth**2 * D0**2 1e-8)) # 只压制以 center 为中心的小范围, 远处恢复为 1 mask[D radius] 1 mask[center[1] - radius:center[1] radius, center[0] - radius:center[0] radius] mask[center[1] - radius:center[1] radius, center[0] - radius:center[0] radius] return mask这里几个参数需要讲清楚。radius 是压制范围经验上取 2~4 个频域像素就够因为能量扩散范围有限bandwidth 控制过渡的陡峭程度bandwidth 越小过渡越陡越接近二值掩模bandwidth 越大过渡越平缓但可能波及旁边真实的频谱分量。D0 的计算是为了让带阻的中心正好落在噪声峰上公式里减去的 W//2 和 H//2 是因为 fftshift 后原点在频谱中心。实际使用中我会对频谱峰值点周围的几个像素逐一测试看反变换后条纹能量下降多少、细节损失多少以此决定 radius 的取值这比纠结理论公式更快。3. 直接能跑的处理流水线频谱自动诊断、陷波生成与图像重建3.1 完整处理流程从读图到输出干净的图像理论铺垫完之后这里给出一套可以直接跑的完整方案。它的流程是读取灰度图 → 计算中心化傅里叶变换 → 在频谱水平轴线上自动寻找噪声峰 → 生成陷波掩模 → 频谱滤波 → 反变换回空间域。核心脚本如下import numpy as np import cv2 def remove_vertical_stripes(image, peak_thresh8, notch_radius3, border_width3): 频域陷波法去除垂直条纹 image: 输入灰度图, uint8 或 float, 建议先转 float peak_thresh: 噪声峰判定阈值, 相对频谱均值的倍数 notch_radius: 陷波压制半径, 单位是频域像素 border_width: 沿水平轴线检测峰的半宽 # 1. 转 float 并计算中心化频谱 img_float image.astype(np.float32) F np.fft.fft2(img_float) F_shift np.fft.fftshift(F) mag np.abs(F_shift) # 2. 在水平中轴线上找噪声峰 H, W mag.shape cy, cx H // 2, W // 2 # 取中轴线附近一小条带, 降低噪声干扰 strip mag[cy - border_width: cy border_width 1, :] profile strip.mean(axis0) # 沿竖直方向平均, 得到一维频谱能量 baseline profile.mean() # 排除原点附近的低频区, 避免把图像本身亮度误判为条纹 profile[:cx - 40] 0 profile[cx 41:] 0 profile[cx - 10: cx 11] 0 # 挖掉 DC 区域不要检测 # 找所有超过阈值的局部极大值点 threshold baseline * peak_thresh peaks [] for i in range(1, W - 1): if profile[i] threshold and profile[i] profile[i-1] and profile[i] profile[i1]: peaks.append(i) # 合并距离过近的峰, 只保留能量更高者 peaks _merge_close_peaks(peaks, profile, min_dist6) # 3. 构造陷波掩模并滤掉所有峰 mask np.ones_like(mag) for peak_x in peaks: mask _apply_notch(mask, (peak_x, cy), notch_radius) # 4. 滤波并反变换 F_shift_filtered F_shift * mask img_filtered np.fft.ifft2(np.fft.ifftshift(F_shift_filtered)).real return np.clip(img_filtered, 0, 255).astype(np.uint8), peaks这里_merge_close_peaks和_apply_notch是两个辅助函数前者把相互距离小于 6 个像素的峰合并只保留谱线能量更高的那个因为条纹周期非常接近时会出多个相邻峰全滤会导致带宽过大后者在指定位置生成高斯陷波并乘到掩模上。peak_thresh 是最难调的参数取 8 意味着「峰的能量至少是全频谱均值的 8 倍」才认为是条带噪声。这个值对大部分红外图像和遥感影像都适用但如果你处理的图像对比度特别高、细节纹理很丰富建议先看频谱剖面再决定。3.2 自动找峰阈值怎么定掩模参数怎么改上面代码里 peak_thresh 和 notch_radius 是两个最关键的手工参数。不少第一次用的人会直接把 peak_thresh 设得很低想「宁可错杀不可放过」结果把图像本身的横向纹理也滤掉了。我一般会让用户先跑一段频谱诊断脚本把水平轴线上的能量分布打印成数值列表肉眼看清楚噪声峰和正常频谱分量的高度差再定阈值。阈值合理性有个经验区间对典型的红外焦平面条纹噪声峰值通常是频谱均值的 15~30 倍peak_thresh 取 8~15 都行对比较微弱的光学条纹峰值可能只高 5~8 倍阈值就得降到 3~4。但注意阈值降到 3 以下时图像自身横向边缘比如地平线、建筑轮廓也会被当成峰这种情况我会改用第 4 章的列回归法而不是硬调阈值。notch_radius 建议从 3 开始检查反变换结果后按需增减 1。另外要注意处理顺序如果图像里叠加了多组不同周期的条纹一次把检测到的所有峰全部滤波往往比逐次处理效果更好因为各峰之间的陷波区域不会重叠。真正需要迭代处理的是那种峰位恰好落在另一个峰的过渡带里的情况这时先滤能量大的峰再重新检测剩余的小峰通常两轮就干净了。3.3 滤波完怎么验收三分钟看三个信号代码跑完别急着存图先做三件事确认效果。第一步看列均值曲线对去条纹后的图像逐列求平均一条平滑的曲线说明列间亮度连续如果再出现锯齿状抖动说明滤得不彻底。第二步看残差图把去条纹前后的图像相减残差图上应该只看到竖条纹的痕迹如果残差里出现明显的物体轮廓说明陷波带宽把真实细节也压掉了需要减小 notch_radius。第三步看频谱对去条纹后的图再做一次 FFT水平轴线上应该看不到突出的峰值。这三步可以用一个简单的评估函数串起来我通常把列均值差分绝对值之和作为量化指标这个指标在本文第 6 章还会展开讲。如果你在这个阶段发现条纹反而变得更明显多半是掩模中心位置算错了——确认一下你的代码里有没有用 fftshift以及峰值坐标是否和掩模生成时使用的坐标系一致这是最常翻车的地方。4. 频域解决不了的部分场景列回归与分块校正的思路4.1 非周期条纹为什么频谱里找不到明显的峰频域陷波不是万能的。红外焦平面阵列的非均匀性噪声、线扫相机的列暗电平不一致、某些 CMOS 传感器的列读出噪声这类条纹的特征是周期不固定、列与列之间的偏差是随机或缓变的频谱上不会出现干净利落的尖峰而是弥散在低频区域陷波滤波器无从下手。这类条纹在空间域的特征反而更明显逐列统计灰度均值后正常图像的列均值曲线是平滑变化的因为场景亮度在横向是渐变的而带条纹图像的列均值曲线上叠加了明显的逐列跳变像锯齿一样。于是有了另一条全局法思路——列回归校正拟合出列均值曲线中的低频趋势部分把它当作「场景真实亮度」用原始列均值减去这条趋势线剩下的残差就是列方向的噪声偏差最后把残差从每一列中减掉。这个思路在工业界叫列均一化或 flat-field correction 的简化版它对非周期性列噪声的效果比频域法可靠得多而且计算量小适合批量处理。4.2 列均值拟合与逐列修正代码与参数细节核心实现不复杂就是统计、拟合、修正三步。下面这段代码用多项式拟合列均值趋势并完成逐列修正import numpy as np def column_regression_correction(image, degree4, block1): 列回归去垂直条纹 image: 输入灰度图 float degree: 多项式拟合阶数, 推荐 3~5 block: 分块列数, 1 表示全局拟合; 图像亮度分布复杂时设为 128 H, W image.shape img_out image.copy() # 分块处理, 每 block 列作为一个拟合区间, 避免场景亮度突变干扰 for start in range(0, W, block): end min(start block, W) col_means img_out[:, start:end].mean(axis0) # 每列均值 col_positions np.arange(start, end) # 多项式拟合低频趋势 coeffs np.polyfit(col_positions, col_means, degree) trend np.polyval(coeffs, col_positions) # 残差 实际列均值 - 低频趋势, 即列噪声 residual col_means - trend # 逐列减去残差 img_out[:, start:end] - residual.reshape(1, -1) return np.clip(img_out, 0, 255) corr_img column_regression_correction(img.astype(np.float32), degree4, block256)多项式阶数 degree 是这里最需要小心的参数。阶数取 1~2 时拟合的是全局线性渐变适合亮度均匀的图像取 3~5 时可以表达场景中横向的明暗过渡比如红外图像中间亮四周暗但如果超过 6拟合曲线会开始跟着条纹本身的抖动走趋势线里混入噪声残差减小条纹就去除不干净了。分块 block 参数处理的是场景里有明显亮度分区的情况比如一行图里有天空有地面亮度相差很大全局一条多项式曲线无法同时拟合两段不同亮度分块后每段单独拟合就自然多了。我一般先用 block256 试如果分块边界处出现亮度跳变就减小到 128 或 64。这里还有一个容易忽略的点列均值用的是算数平均如果图像里有高亮饱和点或坏点会污染该列的均值估计导致这一列修正过头。做之前先用中位滤波或百分位截断把极端值压一下代价小、收益明显。这也解释了为什么有些场景下用列中位数代替均值会更稳——中位数对离群点不敏感但计算量大得多按需选择。4.3 频域法与列回归法怎么选一个对比表实际项目中两条路线可以互补选型主要看条纹形态。这里给一个对比表方便你根据现象快速决策对比维度频域陷波法列回归法适用噪声周期性条纹、固定周期条带非周期列噪声、增益不均频谱表现水平轴线有明显尖峰无明显尖峰, 低频弥散细节保留只压制特定频率, 细节保留好依赖拟合阶数, 阶数高会伤横向渐变参数数量阈值半径, 需要看频谱诊断阶数分块, 需要看列均值曲线典型场景遥感影像条带、扫描仪周期噪声红外非均匀性、线扫暗电平不一致风险点误伤横向纹理、振铃误把真实场景渐变当条纹抹掉一句话总结我的选型习惯先做频谱诊断水平轴线上有清晰尖峰的走频域陷波频谱上看不出峰的走列回归。遥感影像里的传感器条带噪声大多是周期性的低剂量 CT 重建出现的条状伪影也偏周期结构这两类用频域法更干净红外热像仪和工业线扫相机的列噪声则更接近随机偏差直接上列回归法更合适。5. 避坑与常见问题五个翻车现场与排查顺序5.1 忘了 fftshift掩模位置全错图像越滤越花现象滤波后的图像整个变暗或出现大范围条纹状残留和原图相比变得模糊不清检测到的峰值坐标看起来「完全对不上」。原因fft2 的输出中零频在左上角 (0,0)如果不做 fftshift 就去找峰、生成掩模所有坐标都偏移了半个图像尺寸。掩模没有对准任何真实噪声峰滤波相当于把低频能量和高频噪声同时削弱了。解决生成频谱后立刻 fftshift得到以中心为原点的坐标构造掩模时同样以中心为原点反变换前再用 ifftshift 把频谱还原。记住这个铁律fftshift 和 ifftshift 必须成对出现中间所有坐标运算都基于中心化之后的坐标否则等于白做。5.2 陷波只置零一个像素条纹没下去还出了振铃现象条纹颜色淡了一点但轮廓还在图像边缘处出现沿着轮廓扩散的水波纹高频细节看起来脏兮兮的。原因周期噪声的能量在离散频谱里分布在峰位周围至少 2~3 个像素的邻域内单点置零只去掉了很少一部分能量同时锐利的单点截断在反变换时产生吉布斯现象也就是振铃。解决改用高斯型或巴特沃斯型陷波掩模把峰位周围覆盖起来。notch_radius 从 3 起步逐步增加 1观察残差图。如果振铃严重而条纹没净还要检查是不是掩模只有 0 和 1 两个值——二值掩模很容易引发振铃换成渐变过渡的掩模立刻好转。5.3 陷波带宽太宽条纹清了图像也糊了现象条纹确实消失了但整幅图像像蒙了一层纱边缘不再锐利横向纹理比如云层、水面波纹明显变淡。原因带宽设太大时陷波器把条纹峰附近的正常频谱分量也压制了。尤其当图像自身有横向结构时它的频谱能量也分布在水平轴线附近和噪声峰重叠一起被滤掉了。解决把 notch_radius 降回 2~3或者把带宽参数调窄让陷波只覆盖峰的主瓣。用第 3.3 节的残差检查法确认残差图里如果出现物体轮廓说明带宽过宽。遇到图像自身横向纹理很强的情况建议换列回归法不要在频域里硬抠。5.4 拿一张图的频谱掩模去滤一组图每张都出问题现象同一批采集的图像A 图用这套参数去条纹效果很好B 图用了同一套掩模后出现新的条纹或者图像变糊。原因条纹周期不是恒定不变的。传感器温度变化、积分时间调整、电源波动都会让条纹频率发生轻微漂移峰位在频谱上会移动几个像素。固定掩模对不上新图的峰位自然滤不干净还可能压偏。解决对每帧图像重新做峰值检测至少做一次频谱诊断。如果帧率很高、逐帧诊断成本大就用多帧合并的方式生成一个公共掩模——把一组帧的频谱按像素取中位数再找峰这样得到的掩模对整组图都有适用性具体做法见下一章。5.5 对 uint8 图像原地做减法输出出现新的方块噪声现象列回归法跑完后图像里出现一行行规律排列的亮暗块或者原先的条纹没去掉又多了雪花点。原因uint8 是无符号 8 位整数范围 0~255做减法时低于 0 的值会溢出变成 255 附近的值形成亮斑最后如果用 np.clip 处理不当这些溢出就是新的噪声源。列回归修正时残差里有正有负直接在 uint8 上运算必然出问题。解决全流程用 float32 计算输入先转浮点所有减法在浮点域进行最后统一 clip 到 [0, 255] 再转回 uint8。这里没有捷径我第一次用 uint8 直接做列回归时就吃了这个亏输出图里多出一堆方块噪声排查了半天才发现是数据类型的问题。6. 进阶条纹能量指标与多帧自动掩模生成6.1 一个能写进验收报告的条纹能量指标去条纹做完了怎么跟上下游交代效果PSNR 和 SSIM 对整幅图敏感但条纹是细结构它们未必能灵敏地反映条纹能量的变化。我习惯用一个专门指标列均值曲线的相邻差分绝对值之和。它衡量的是「列与列之间亮度跳变的剧烈程度」条纹越重这个值越大去干净了这个值会显著下降。下面这段代码可以算出来def stripe_energy(img): 条纹能量: 列均值曲线相邻差分绝对值之和 值越大说明列间亮度跳变越剧烈, 条纹越重 col_means img.astype(np.float32).mean(axis0) diff np.abs(np.diff(col_means)) return float(diff.sum()) before stripe_energy(img) after stripe_energy(result) print(f去条纹前: {before:.2f}, 去条纹后: {after:.2f}, 下降: {(before-after)/before*100:.1f}%)这个指标对周期性和非周期性条纹都有效而且计算一次只要几毫秒。实践中红外条纹图像去完后这个值通常会下降 60%~85%如果下降不到 30%说明滤波没滤干净或者参数没找到点上。把它写进验收报告比截图更有说服力。需要注意的是这个指标要和残差检查配合使用——防止为了降低指标而把图像整体模糊化那种做法指标很漂亮但细节全没了。6.2 多帧自动掩模让一组图像共享同一个陷波器批量处理一组图像时逐帧诊断很烦。我的做法是先取 10~20 帧有代表性的图分别做 FFT取幅度谱的逐像素中位数合成一张「平均频谱」。中位数操作可以压掉单帧图像自身的纹理结构只留下所有帧共有的周期噪声峰。然后在这个合成频谱上做峰值检测生成统一的陷波掩模应用整组图像。这样参数只需要调一次后续全部自动。如果传感器工作状态改变比如换了积分时间重新采集一组帧再算一次掩模就行。从这个角度说频域陷波法的核心其实不在「滤波」本身而在「定位噪声峰的准确性」。定位准了剩下就是把峰值区域轻轻盖住。我处理过的项目里凡是效果翻车的基本都是定位环节偷懒了——没有做频谱诊断就套参数。从那以后我每次拿到一张条纹图都会先花三分钟跑一遍频谱诊断看一眼水平轴线上的峰分布再决定走频域还是列回归这个习惯帮我省下了大量反复调参的时间。希望帮到你。本文还有配套的精品资源点击获取
返回列表