侧扫声纳图像的去噪问题我断断续续折腾了有一阵子。手里攒了几批不同海区的实测数据跑过均值滤波、中值滤波、维纳滤波、双边滤波效果都不能说满意。要么背景是干净了但目标边缘糊成一团要么边缘保住了可斑点噪声还在那里晃眼。后来我把模糊加权平均Fuzzy Weighted Average和卡尔曼滤波Kalman Filter串在一起用才算是把抑制斑点噪声和保留弱目标轮廓这两件事同时做到了。这篇就聊聊这套复合去噪方案的设计思路、实现细节和我在实际数据上踩过的坑。这个内容适合谁看如果你是做海底测绘、水下目标探测、管线巡查或者水下考古声纳图像处理的手里正好有侧扫声纳数据或者类似的水声图像那这篇文章可以直接参考。就算你只是对图像去噪算法感兴趣这套空间域预处理时间域估计的组合思路放到视频帧序列、遥感时序影像里也有借鉴价值。1. 看懂侧扫声纳图像里的噪声去噪才有方向1.1 侧扫声纳成像到底是怎么一回事侧扫声纳不是像相机那样拍一张照片。它的换能器阵列装在拖体两侧向海底发射扇形的声脉冲然后接收海底反射回来的回波。声波往两侧传播碰到起伏的地形、礁石、沉船、管线反射强度不同接收到的回波幅度也就不一样。把回波幅度按传播时间也就是斜距展开成一条线这条线叫一个ping拖体往前跑不断发射接收一个个ping排在一起就组成了一幅侧扫声纳瀑布图。在瀑布图里横轴是斜距方向纵轴是航迹方向。像素灰度的高低代表回波强度的强弱。强反射体比如礁石在图上呈亮色它的背后会因为声波照射不到而形成声学阴影呈暗色。地质识别和目标判读靠的就是这些亮暗纹理和几何形态。关键点在于侧扫声纳图像的噪声和普通光学图像根本不是一回事直接用光学图像那套去噪思路大概率会翻车。先把这个搞清楚后面的算法选型才有依据。1.2 噪声从哪里来三类噪声叠加在一张图上侧扫声纳图像里的噪声至少可以分成三类斑点噪声speckle noise海底微结构对声波的随机散射造成。声波在粗糙界面和介质不均匀体上发生相干叠加形成颗粒状的明暗斑纹。这种噪声本质上是乘性的它在图像里表现得像一层胡椒面和信号强度成正比。信号强的地方噪声也强信号弱的地方噪声也弱。水柱回波water column echo发射脉冲的旁瓣直接耦合进接收通道在近场形成一条高亮条带。它不是随机的有固定的位置特征但灰度动态范围极大经常把整幅图的自动增益给带偏。接收机底噪和环境噪声电子元器件热噪声、海流和气泡的散射噪声这类是加性的表现为均匀分布的低幅值干扰。还有一类是时变增益TVG补偿不当时引入的灰度斜坡误差属于系统误差不算严格意义的噪声但会给后续处理造成干扰。1.3 为什么普通去噪算法在侧扫图像上容易翻车先说我试过的几种常见算法的实际表现高斯滤波/均值滤波假设噪声是加性高斯白噪声但侧扫图像的主要干扰是乘性斑点噪声平滑窗口稍微开大一点亮目标周围的暗阴影区就会被污染目标边缘变成渐变过渡看起来像罩了一层雾。中值滤波对孤立椒盐噪声有效对斑点噪声也有一定抑制但它的窗口形状和大小直接影响纹理细节。侧扫图像里的弱目标比如一小段管线或半埋目标灰度变化幅度本身就不大中值滤波一上去这些弱目标经常直接消失。维纳滤波在平稳随机场假设下效果不错但海底底质往往存在突变比如泥底和沙波交界平稳假设不成立结果就是边缘处出现振铃和过冲。双边滤波是空间域里边缘保持能力不错的但它对噪声方差敏感斑点噪声强度大的时候灰度差权值容易被噪声主导边缘照样保不住。这些方法有一个共同问题它们只用了当前这一帧图内部的邻域信息。而侧扫声纳图像有一条非常特殊的先验——目标在相邻的多个ping上会连续出现同一位置在不同ping之间的灰度变化是缓慢的。这个时间维度的信息空间滤波完全没用上。所以我当时的判断是与其纠结换哪种空间滤波器不如换个思路把空间域去噪和时间序列估计两件事组合起来。这也是这套复合算法的出发点和核心逻辑。2. 模糊加权平均在空间域能做什么不能做什么2.1 模糊加权平均不是简单的加权平均模糊加权平均这个名字容易让人误解以为就是把均值滤波的等权改成不等权。实际上它的关键在于权重不是预先固定死的而是根据邻域像素与中心像素的相似程度实时计算出来的。相似程度用灰度差来衡量灰度差越小说明这个邻域像素和中心点越有可能属于同一片均匀底质权重就越大灰度差越大说明它越可能属于边缘、目标或者阴影区权重就越小。用式子表达就是对中心像素x(i,j)取其邻域窗口W内的每一个像素x(p,q)计算灰度差d |x(p,q) - x(i,j)|然后通过模糊隶属度函数把d映射成权重w(p,q)。输出像素值就是y(i,j) Σ w(p,q) · x(p,q) / Σ w(p,q)如果w恒等于1这个式子就退化成均值滤波。模糊加权平均的本质是让滤波器看到局部结构在均匀区域所有邻域像素权重都很大等效于强力平滑在边缘附近跨过边缘的像素灰度差大、权重小主要权重落在同侧像素上边缘就能保下来。2.2 隶属度函数怎么选从三角形到高斯型模糊加权平均的模糊就体现在隶属度函数上。我实际用下来主要试过两种三角形隶属函数w max(0, 1 - |d| / T)T是灰度差上限。超过T的像素权重为0。计算极快但是转折点T需要调得比较仔细T偏小会导致边缘附近权重断崖式变化去噪结果容易出现阶梯感。高斯隶属函数w exp(-d² / (2σ_d²))σ_d控制衰减速度。曲线平滑权重连续参数含义也直观。我最终在项目里选的是高斯型。还有一个很容易被忽略的细节单纯用灰度差做权重在斑点噪声强的时候不可靠因为噪声本身就会造成灰度差抖动。更稳的做法是加一个空间距离权重把两者乘起来w(p,q) exp(-d² / (2σ_d²)) · exp(-dist² / (2σ_s²))第一个因子管灰度相似性第二个因子管几何邻近性。距离中心越远的像素就算灰度接近也不能给太大权重否则会把远处结构拖进来。这个双因子形式和双边滤波在形式上很像但推导路径不同——双边滤波是从边缘保持的工程经验出发而模糊加权平均是从模糊逻辑的隶属度概念出发。工程上不必纠结谁源自谁好用就行。2.3 单靠它瓶颈在哪里模糊加权平均单独用效果比均值滤波、中值滤波好不少但有一个我绕不开的问题它只能利用单帧内的结构信息对高噪声场景的极限抑制能力有限。原因是侧扫图像中目标和底质的灰度差往往不大如果噪声方差已经接近这个灰度差模糊逻辑就很难区分灰度差是边缘还是噪声。强行加大σ_d权重差异变小滤波器趋向均值滤波边缘保不住强行减小σ_d均匀区域的噪声像素互相之间灰度差也不小权重被压得很低平滑能力又不够。这是一个此消彼长的尴尬平衡。所以在高噪声数据上我得到的结论是模糊加权平均适合做预处理把明显噪声压下去但不足以保证最终的图像质量。必须再找一个能在另一个维度上有互补能力的算法。这就轮到了卡尔曼滤波。3. 卡尔曼滤波补上的时间维度沿航迹方向的最优估计3.1 把侧扫图像当成一组并行的像素时间序列卡尔曼滤波最常见的应用场景是导航、目标跟踪这些处理的是随时间变化的状态量。侧扫图像是一张二维矩阵怎么把卡尔曼滤波用上去关键操作是视角转换把图像纵轴航迹方向看成时间轴一次发射对应一个时刻。位于同一斜距位置、不同ping上的像素实际上描述的是同一个海底位置附近在拖体前进过程中的回波变化。因为拖体运动是连续的真实海底的后向散射强度在相邻ping之间变化很小主要起伏来自噪声。这样一来侧扫图像的每一列都可以看成一条独立的离散时间序列而我们要做的就是对每条序列做去噪估计得到真实回波强度的最优值。拿视频来类比会更好懂视频里同一物体在相邻帧之间位置只移动一点我们可以在时间方向做平滑来降噪。侧扫声纳沿航迹方向就是天然的时间方向只不过这个方向在成像时被拉伸成了图像的纵轴。3.2 卡尔曼滤波的核心流程预测再修正卡尔曼滤波不需要我在这里用严格的数学推导但关键在于理解那五个公式在图像去噪里的具体含义。对图像某一列中位置k处的像素建立状态方程和观测方程状态方程x_k x_{k-1} w_k含义是相邻ping之间真实灰度基本不变w_k是过程噪声方差为Q代表真实灰度随空间位置变化的剧烈程度。观测方程z_k x_k v_k含义是我们测到的灰度等于真实灰度加测量噪声v_k方差为R代表当前灰度值的可信程度。每一行像素要跑两遍第一遍是预测状态预测x_pred x_est因为状态转移矩阵取1也就是我猜这次灰度和上次差不多协方差预测P_pred P_est Q不确定性因为过程噪声而增加第二遍是更新卡尔曼增益K P_pred / (P_pred R)状态修正x_est x_pred K · (z - x_pred)协方差修正P_est (1 - K) · P_pred公式不复杂核心就一个卡尔曼增益K决定了在你更相信模型预测值还是当前观测值之间做加权。如果K大说明当前观测可靠就多采信观测如果K小说明观测受噪声污染重就多靠预测值往前推。3.3 Q和R的比值决定了平滑力度Q和R是卡尔曼滤波里最关键的参数但千万别把它们当成两个独立变量去调真正起作用的是Q和R的比值。Q大或者R小滤波器认为真实灰度变化剧烈、观测噪声小K会偏大输出更跟手图像细节保留更多但平滑效果有限。Q小或者R大滤波器认为真实灰度平缓、观测噪声大K会偏小输出更平滑噪声压得狠但底质突变处的细节也会被磨掉。我自己的调参习惯是固定R1只动Q。这样做的好处是K的表达式变成K P_pred / (P_pred 1)P_pred的量级直接决定行为好理解也好调试。对于普通侧扫数据Q从0.005到0.1这个区间去试大致能覆盖从轻平滑到重平滑的范围。这里插一句我的使用体会单独用卡尔曼滤波在侧扫图像上跑效果其实没有想象中好。因为单帧内的斑点噪声经过逐列递推后会被平均掉一部分但列与列之间相互独立横向的噪声纹理也就是斜距方向的斑点依然明显图像看起来有纵向条纹感。这反过来印证了卡尔曼滤波擅长在时间方向平滑但空间结构还得靠空间域滤波来处理。两件事根本是互补的。4. 复合策略怎么落地串联架构与实现细节4.1 为什么是串联而不是并联融合把两种算法组合起来直观想法是并联融合模糊加权平均出一个结果卡尔曼滤波出一个结果然后按权重加起来。这个方案听起来很均衡但工程上问题很大——两者输出的置信度怎么定融合权重怎么跟图像局部特征自适应这些问题的答案本身又是新超参数调试成本成倍增加。我最终选择的是串联结构先对原始图像做模糊加权平均得到初步去噪图像再把这个初步结果输入到卡尔曼滤波沿航迹方向逐列递推得到最终输出。串联结构的逻辑是第一级模糊加权平均把空间域的斑点噪声压掉一部分信噪比提高之后第二级卡尔曼滤波面对的输入不再是原始强噪声时间序列的平稳性更好递推估计的误差更小。两级各自面对的问题都更简单整体反而稳定。4.2 第一步模糊加权平均做帧内预处理实际操作时我对每个像素取7×7邻域窗口分别计算灰度差权重和空间距离权重乘积做归一化。灰度差权重的高斯尺度σ_d取图像灰度标准差的0.5倍左右空间距离尺度σ_s取1.8对应7×7窗口内中心附近像素作用最强。有几个实现细节必须注意边界处理图像边缘处窗口会越界最简单的办法是只取窗口有效区域然后对权重重新归一化。别用零填充零填充会把边界灰度显著拉低。灰度域选择侧扫声纳原始回波动态范围很大一般先做对数压缩再显示。我测试发现直接在对数域做模糊加权平均效果最好。因为对数域里乘性斑点噪声近似变成加性噪声灰度差权重更稳定。处理完再做指数变换回去。权重归一化如果Σw太小比如窗口内全是强边缘像素这个点就别过度平滑可以加一个保护判断直接采用中心像素值。处理完这一级图像的整体视觉已经干净不少但纵向航迹方向的低频噪声依然存在交给下一级。4.3 第二步卡尔曼滤波做跨帧递推预处理之后把图像逐列取出来对每一列独立运行卡尔曼滤波。初始化时第一行像素的状态估计x_est直接取该像素灰度值P_est设成1。然后逐行往下递推。这里要说一个我在代码里踩过的坑卡尔曼滤波的递推状态要按同一列传递不能按整幅图的全局坐标传递。侧扫图像横轴是斜距不同斜距位置对应完全不同的海底物理位置没有任何递推关系。弄混了维度出来的图像会出现横向的色带整幅图一道一道的像水彩晕染。伪代码大概长这样import numpy as np def fuzzy_weighted_average(img, win_size7, sigma_d12.0, sigma_s1.8): 空间域模糊加权平均预处理 h, w img.shape out np.zeros_like(img) half win_size // 2 # 预生成空间距离权重模板 yy, xx np.mgrid[-half:half1, -half:half1] spatial_w np.exp(-(xx**2 yy**2) / (2 * sigma_s**2)) for i in range(h): for j in range(w): r0, r1 max(0, i-half), min(h, ihalf1) c0, c1 max(0, j-half), min(w, jhalf1) sub img[r0:r1, c0:c1] center img[i, j] gray_w np.exp(-((sub - center) ** 2) / (2 * sigma_d**2)) sw spatial_w[r0-ihalf : r1-ihalf, c0-jhalf : c1-jhalf] weights gray_w * sw s weights.sum() if s 1e-6: out[i, j] (sub * weights).sum() / s else: out[i, j] center return out def kalman_along_track(img, Q0.01, R1.0): 沿航迹方向逐列卡尔曼滤波 h, w img.shape out np.zeros_like(img) for col in range(w): x_est img[0, col] P_est 1.0 out[0, col] x_est for row in range(1, h): # 预测 P_pred P_est Q # 更新 K P_pred / (P_pred R) x_est x_est K * (img[row, col] - x_est) P_est (1 - K) * P_pred out[row, col] x_est return out # 使用示例img为对数压缩后的灰度图 denoised kalman_along_track(fuzzy_weighted_average(img))真实项目里不建议拿纯Python双层循环跑大图速度太慢。模糊加权平均那步可以用卷积实现加速分别对灰度差权重和空间距离权重做可分离卷积近似或者用Numba把循环编译一下。卡尔曼滤波部分是串行递推很难完全矢量化但每一列之间互相独立可以多线程按列并行。4.4 两级之间的参数如何配合串行架构下两级参数不是相互独立的。具体表现是第一级模糊加权平均如果开得比较重窗口大、σ_d小输出已经很平滑第二级卡尔曼滤波的Q可以设小一点免得过度平滑。第一级如果开得轻窗口小、σ_d大第二级就要靠Q小一些来增加纵向平滑否则总噪声压制不够。反过来调也没有问题但整体原则是两级分担的任务量要均衡别让一级冲得过猛另一级空转。我自己在项目里的固定搭配是7×7窗口σ_d 灰度标准差的0.5倍σ_s 1.8Q 0.02R 1。这套参数在各种底质类型下表现都比较均衡。5. 去噪效果怎么量化指标选择与实测对比5.1 四个指标分别盯着不同的担心评价去噪效果不能只看看起来干不干净要量化。我常用的四个指标指标全称关注点PSNR峰值信噪比像素级重建误差越高越好SSIM结构相似性亮度、对比度、结构三方面综合越高越好ENL等效视数均匀区域均值平方与方差之比衡量斑点抑制能力越高越好EPI边缘保持指数处理前后边缘强度比值衡量边缘是否被磨掉越接近1越好PSNR和SSIM需要参考图像一般用在合成实验中ENL和EPI能直接在真实数据上算ENL取图像中的均匀底质区域EPI取目标边缘区域。5.2 合成噪声数据上的定量测试为了能算PSNR和SSIM我用一幅信噪比较高的侧扫声纳图像作为参考叠加模拟的乘性瑞利斑点噪声生成带噪图。然后分别用均值滤波、中值滤波、纯模糊加权平均、纯卡尔曼滤波和本文复合算法处理记录指标。下面是一次典型测试的数据算法PSNR(dB)SSIMENLEPI带噪原图24.30.7112.5—均值滤波 5×527.10.7932.00.68中值滤波 5×526.80.7828.40.61模糊加权平均 7×728.80.8441.20.79卡尔曼滤波 Q0.0228.50.8338.60.75复合算法30.60.8955.80.86复合算法在PSNR上比单独使用两种方法提高了约1.8~2.1dBENL提升更明显说明斑点噪声被压制得较彻底。EPI 0.86在同类算法里算不错边缘保持和噪声抑制的平衡是符合预期的。SSIM的提升幅度没有PSNR那么显著但从0.83到0.89也算实打实的进步。这是因为卡尔曼滤波沿航迹方向平滑后图像的结构信息整体更接近参考图。5.3 真实数据上的视觉对比定量指标好看不等于真实数据上就能用。我拿了几幅含沙波、礁石和人工目标的真实侧扫数据做目视检查重点看三个地方均匀泥底区域复合算法处理后背景颗粒感明显变轻灰度起伏趋于平缓ENL的提升在视觉上表现为画面更干净。沙波纹理区域沙波的波纹边缘没有出现涂抹感波纹的脊线和谷线在去噪后依然清晰可辨。这一点让我比较满意因为之前用均值滤波时沙波纹理经常被平滑成一片平地。礁石与阴影交界处这是最容易翻车的地方。模糊加权平均在阴影边缘靠灰度差权重保住了边缘卡尔曼滤波在时间方向的平滑又把残余的斑点噪声进一步压低。实际效果是礁石轮廓清楚阴影区内部干净交界处没有明显的highlight或暗带。当然也有不理想的情况当目标在航迹方向上的尺寸只有1~2个像素时比如细长的缆绳卡尔曼滤波的递推会把它的灰度向周围摊开表现为目标变宽、对比度下降。这个问题后面细说。6. 参数调优与踩坑复盘实操中的几个关键问题6.1 参数联动关系与调参顺序这套算法总共有4个主要参数窗口大小、σ_d、σ_s、QR固定为1。如果一上来就四个参数乱调必定头大。我建议的顺序是先固定σ_s 1.8、窗口7×7只动σ_d。观察均匀区域去噪效果σ_d越小平滑越强边缘保持越弱把σ_d调到边缘刚好清晰的临界点。然后动Q。Q越小纵向平滑越强噪声越少但航迹方向细节损失越大把Q调到图像不再有纵向条纹感的程度。最后微调窗口大小。窗口越大模糊加权平均的时间开销越大边缘保持越差只在噪声特别强的情况下加大到9×9。按这个顺序能在一两轮实验内找到可用的参数组合。6.2 底质突变区域的拖影问题这是我遇到的最难处理的问题必须拿出来单独说。卡尔曼滤波的模型假设是相邻ping之间真实灰度变化很小。但海底不会总听话——泥底到沙波的过渡带、礁石边缘、人工目标的强反射边界真实灰度在几个像素内就会跳变。这种情况会出现什么滤波器的预测值还停在突变前的灰度上而观测值已经跳到新灰度卡尔曼增益K不够大时输出就会在突变位置附近形成一段渐变的拖尾。解决思路有两类第一类是检测突变并重置滤波器。当|z_k - x_pred|超过某个阈值时把K直接置为1也就是完全信任当前观测。相当于告诉滤波器这里发生了突变别拿旧模型预测了重新跟。实现很简单if abs(img[row, col] - x_est) threshold: x_est img[row, col] P_est 1.0 else: K P_pred / (P_pred R) x_est x_est K * (img[row, col] - x_est) P_est (1 - K) * P_pred阈值怎么取我一般用第一级模糊加权平均后的图像在均匀区域的残差标准差的3倍左右大于这个值就认为是真实突变而不是噪声。第二类是引入衰减记忆因子。让P_est在每次递推时乘以一个略大于1的系数使滤波器对旧数据的记忆随时间逐渐减弱。这就变成了衰减记忆卡尔曼滤波对时变信号更敏感但调参难度比第一种大一点。我最终在项目里采用的是第一种简单直接在几批数据上都很稳定。6.3 效率优化与工程落地侧扫声纳图像动辄上万行乘以几千列直接逐像素跑模糊加权平均和逐列跑卡尔曼滤波总耗时非常可观。我踩过的性能和内存坑有三个模糊加权平均的双层循环不能用。在Python里把灰度差权重和空间距离权重的分子分母分别做卷积用scipy.ndimage.convolve然后把两者相除等价于逐像素计算但速度能快几十倍。卡尔曼滤波按列并行。每一列的递推完全独立没有跨列依赖。我用多线程一次同时处理多条列在8核机器上能做到接近7倍加速。数据类型和内存复用。原始图像如果是16位整型去噪过程建议先转成float32计算不要用float64速度差不少。中间结果尽量原地更新避免反复申请大数组。还有一个小细节去噪前先做对数压缩去噪后做指数还原这个压缩-处理-还原的链路在侧扫声纳处理里几乎是标配。如果直接对线性域回波做滤波少数强目标的灰度值会主导整个窗口的权重导致弱目标被压制。6.4 关于这套组合我的最后几点体会把模糊加权平均和卡尔曼滤波复合起来并不是一个多复杂的高大上技巧但它的思路是顺的先搞清楚图像里的噪声是什么特性再选一个能处理空间特性的手段选一个能处理序列特性的手段然后把它们按合理的顺序接起来。侧扫声纳图像天然同时具备空间纹理和序列相关性两个属性所以这套组合能起作用。我在实际项目中还试过把模糊加权平均换成分数阶微分增强后再做卡尔曼滤波效果也有提升但计算量更大参数更多最后的性价比反而不如现在的方案。如果你手里的数据噪声类型和我不太一样建议从串联结构入手换掉其中一级、保留另一级而不是整体推翻重来。这个架构的容错能力比我预想的要好。