等值线平滑的噪声地图矢量化优化方法

2022-04-21李志强刘兰华潘永琛

噪声与振动控制 2022年2期

李志强，刘兰华，潘永琛

（中国铁道科学研究院集团有限公司节能环保劳卫研究所，北京100081）

噪声地图是将噪声数据渲染于地图上，并以不同的色彩代表噪声的大小，直观地展示噪声的空间分布、传播衰减及影响范围的专题地图。噪声地图可用于辅助定位未知噪声源，确定各类噪声源贡献量，指导相关部门采取合理的环境噪声控制措施，已成为当前国内外环境噪声监管中重要技术手段[1]。按照噪声数据来源通常可将噪声地图分为监测型、计算型、监测-预测复合型3类[2]。这3种类型噪声地图在噪声渲染时都可通过对一定数量和密度的离散点处的实测或预测噪声值作插值处理，如克里金插值、反距离加权插值、样条函数插值等，形成绘图范围内的栅格型噪声云图。

在地理信息系统（GIS）软件中利用《环境影响评价技术导则声环境》中铁路噪声预测模型进行二次开发并绘制某长度为120 km、铁路两侧250 m 范围内噪声地图时，设计像元尺寸为5 m×5 m，得到的单幅栅格型噪声云图数据量高达1.5 GB。当噪声地图比例尺较小、栅格像元尺寸远小于地图展示范围时，噪声云图的展示效果基本不受像元尺寸影响。当观察某一较小范围内噪声分布并将噪声地图比例尺逐渐放大时，栅格型数据的缺陷会逐渐显现，表现在不规则绘图边界处锯齿化严重、噪声变化显著的区域有明显的马赛克现象等，给噪声地图展示的连续性和美观性带来较大影响[3]。像元尺寸对图像展示效果的影响问题在其他领域的栅格型数据可视化中也普遍存在[4－5]。

解决栅格型数据锯齿形边界和马赛克现象最简单的方式是减小像元尺寸，使其满足将地图放大到最大比例尺条件下的应用需求，但该方法会导致栅格数据量大幅增加，如：空间分辨率增大n倍，像元边长缩小1/n，栅格文件占用的存储空间会增大为原来的n2倍。栅格图像文件过大，对图像数据的生成、存储、迁移、发布、批量处理和展示等均不利，这种图像格式不适合作为大尺度的噪声云图绘图格式。

为提升噪声地图中数据展示效果，并降低噪声云图数据量，本文对栅格图像的双线性插值法[6]和双三次卷积插值法[7]的重采样平滑优化效果及栅格数据直接矢量化效果进行了对比研究，进而提出了一种基于等值线平滑的栅格数据矢量化的优化方法，不仅使噪声云图展示效果不再受地图缩放的比例尺影响，还使噪声云图的数据量也大幅降低。

1 栅格数据重采样方法

重采样是栅格图像数据平滑的主要方法，常见的算法有双线性插值法和双三次卷积插值法等。

1.1 双线性插值法

双线性插值法是线性插值法的二元扩展，通过待插值点处最邻近的4个像元的像素值的权重关系对目标点处进行插值，其原理如图1 所示。假设栅格数据像元边长为单位1，已知中心位于坐标(i，j)、(i+1，j)、(i，j+1)、(i+1，j+1)的4 个像元的像素值分别为f(i，j)、f(i+1，j)、f(i，j+1)、f(i+1，j+1)，为获得图1中某一点(i+p，j+q)的值f(i+p，j+q)，需要对x坐标和y坐标分别进行线性插值（结果与x、y插值顺序无关）。首先在y坐标为j和j+1的位置分别向x方向进行线性插值得到f(i+p，j)和f(i+p，j+1)：

图1 双线性插值示意图

之后在x坐标为i+p的位置向y方向进行线性插值得到f(i+p，j+q)：

将式(1)、式(2)待入式(3)得到f(i+p，j+q)为

根据式(4)，利用待插值点最邻近4 个像元的像素值，即可得到待插值点处的值。进行双线性插值法计算时仅需要考虑插值点最近的4 个像元，计算公式简单，可作为快速插值算法对栅格数据进行平滑，在放大图像时快速抑制图像马赛克现象效果较好，但由于该方法自带低通滤波性质，因此会引起图像高频分量受损，在像素值变化剧烈及轮廓边界处图像会趋于模糊化。

1.2 双三次卷积插值法

双三次卷积插值法与双线性插值法计算过程类似，但其根据对待插值点最邻近的16个像元的像素值计算加权平均值，并基于利用三次多项式逼近最佳插值函数sinπx/πx的方式解决线性法插值的图像模糊问题。如图2 所示，已知图中待插值点(x，y)处最邻近的16 个像元的像素值f(i,j)（i=1，2，3，4；j=1，2，3，4），满足这16 个像元像素值的插值函数为g(x,y)，g(x,y)函数表达式可表示为：

图2 双三次卷积插值示意图

式中：W(x)、W(y)分别为x方向和y方向的插值卷积核函数。根据最佳插值函数特征，三次卷积核函数表达式可写为：

s为待插值点与邻近像元中心在x或y方向的距离（|x-i|或|y-j|），Ai～Di为各多项式待定系数。

三次卷积核函数为对称函数，三次卷积核函数为对称函数，|s|在[0，2]区间内连续且其1 阶导数连续，将上述边界条件代入W(s)，求解各参数，可得W(s)的通用表达式为：

式中：a为满足三次卷积核函数边界条件的任意值。

为选择合理的a值，在图3中对不同a值的三次卷积核函数曲线与sinπx/πx曲线进行对比，根据对比结果可知，当a=-1.0时三次卷积核函数曲线与sin πx/πx曲线一致性最高，因此将a=-1.0 代入式(7)进一步确定三次卷积核函数为：

图3 不同a值的卷积核函数曲线与最佳插值函数曲线对比

与双线性插值方法相比，双三次卷积插值方法能较好地保持图像的高频部分，在有效平滑栅格数据的同时避免了双线性插值方法所引起的图像模糊问题的产生，但由于涉及最邻近16 个像素值的加权，计算过程相对复杂，计算量较大，放大图像时重采样耗时也更长，当图像数据量较大或计算机性能较差时图像重采样过程可能存在明显可觉察的时延现象。

1.3 重采样效果分析

图4为我国某高速铁路沿线局部区域的栅格型噪声云图的展示效果，在该视图范围内像元尺寸为5 m×5 m的噪声云图在绘图区边界处存在明显的锯齿形，在建筑物声影区边缘附近也可观察到像元马赛克现象。下面分别采用双线性插值法和双三次卷积插值法优化噪声云图，并对比分析不同方法的图像显示优化效果。

图5 为采用双线性插值法优化后的噪声云图，对比图4 可见，在绘图区边界的锯齿形虽未完全消除，但已被模糊处理，锯齿角变得柔和；建筑物声影区边缘附近的马赛克现象在一定程度上被抑制，但放大局部区域后，仍能观察到一定的马赛克残留痕迹。

图4 噪声地图的初始噪声云图

图5 采用双线性插值优化后的噪声云图

图6为采用双三次卷积插值法优化后的噪声云图，对比图4 可见，在绘图区边界的锯齿形仍然显著，未受到插值的影响；建筑物声影区边缘附近的马赛克现象消失，噪声云图变化趋于连续，将局部区域放大后也没有再出现明显马赛克痕迹。

图6 采用双三次卷积插值优化后的噪声云图

双线性插值法和双三次卷积插值法虽然都能在一定程度上平滑栅格型噪声云图，提升栅格图像放大后的视觉效果，但各自仍存在较为明显的不足之处，且插值后噪声云图数据量并未减小，因此仍需进一步研究改进的新方法。

2 栅格图像矢量化方法

矢量图像是将基本图形要素如点、线、面等以数学方式表示的图像，其图像可无限放大而不失真，在表征相同效果的图像时，矢量图像占用的存储空间远小于栅格图像，因此矢量化是减小噪声云图数据量的重要手段。

2.1 直接矢量化

所使用的GIS软件可直接对栅格型噪声云图按照像素值进行矢量化转化，转化后的部分区域视图如图7 所示。新生成的矢量型噪声云图缩小为370 MB，数据量约减小为原来的1/4。对比图4 的栅格型噪声云图，直接矢量化后得到的噪声云图在绘图区边界处的部分锯齿形栅格转变为直线，对改善边界锯齿形有一定效果，但矢量化后得到的各面要素边界均由折线组成，不同噪声值的相邻面要素之间边角分明，对噪声传播及空间分布特征的展示效果大打折扣，没有达到优化视图效果的目的。

图7 直接矢量化后的噪声云图

2.2 基于等值线平滑的矢量化

对栅格型噪声云图直接进行矢量化所引起的边角锐化现象的显著性与像元大小及像素值空间变化率有一定关系，虽然通过插值方式减小像元尺寸能使矢量云图中面要素边界看起来平滑，但将局部区域进一步放大时该问题还会继续出现。为解决栅格图像直接矢量化所造成的问题，本文研究了一种基于等值线平滑的矢量化方法，下面分步介绍该矢量化方法。

第一步，利用栅格型噪声云图或密度较高的离散点噪声值，初步绘制一定噪声间隔的等值线要素集（本文以图4所示的完整栅格数据、1 dB(A)间隔为例作介绍）。由于栅格数据及离散点数据具有空间离散性，初步形成的等值线由一条条折线段组成，两相邻折线段之间有尖锐的拐点或折点，这与实际的噪声传播及分布规律不符。

第二步，为模拟相对均匀、连续的噪声变化特征，利用插值方法对折线型等值线进行平滑处理，得到过渡流畅、轮廓平滑的曲线型等值线要素集。常见的等值线平滑方法有二次多项式平滑法、距离加权最小二乘平滑法、负指数加权平滑法[8]、三次样条函数平滑法[9]、贝塞尔插值平滑法和指数核多项式近似平滑法[10]等。图8 为平滑前的等值线和利用指数核多项式近似法平滑后的等值线对比。

图8 等值线平滑前后对比

第三步，利用绘图区边界、建筑物轮廓线以及平滑后的等值线，将绘图区域切割、划分为由多个面要素组成的面要素集。此时即得到噪声分布区域不同噪声值所在的矢量图层，但各面要素中尚没有噪声值数据。

第四步，对各面要素按照边界处包含的等值线赋予噪声值属性。例如，构成某面要素的两条等值线的值分别为70 dB(A)和71 dB(A)，则赋予该面要素起始值为70 dB(A)，结束值为71 dB(A)。

第五步，按照面要素噪声起始值或结束值属性对各面要素关联色彩，相同噪声值的面要素渲染同种颜色，得到矢量化的噪声云图。

矢量化后噪声云图如图9 所示，优化后的矢量噪声云图不仅在绘图区边界处非常平滑，在绘图区内部的噪声分布随空间变化的层次感也更加清晰，噪声传播规律及分布特征的表现效果不弱于经双三次卷积插值优化后的栅格型噪声云图，且将局部区域放大后不再出现马赛克现象，缩放过程中噪声云图也能无时延地即刻显示。矢量化后的噪声云图仅有51 MB，数据量减小为原来栅格云图的1/30。

图9 优化后的矢量噪声云图

上述基于等值线平滑的矢量化方法绘制精度灵活可调，如当噪声地图中噪声精度要求为0.5 dB(A)时，仅需将第一步操作中的等值线间隔设置为0.5 dB(A)即可实现。进一步减小等值线间隔虽然会增大矢量云图数据量，但仍远小于栅格云图数据量。具体选用的精度应根据绘图范围内噪声的大小范围和实际需求合理确定。根据《环境噪声监测技术规范噪声测量值修正》（HJ706—2014）[10]，环境噪声监测数据需修约到整数，1 dB(A)的间隔设置可满足我国当前的环境噪声监管需求。

3 不同优化方法对比分析

表1对优化栅格型噪声云图不同方法的特点进行了对比分析。对栅格型数据直接进行插值平滑的优化方法不能降低文件数据量，虽然在一定程度上能提升图像显示质量，但平滑质量较高的算法计算量也更大，可能带来显示延迟的问题，适合用于绘图范围和像元尺寸均不大的栅格图像的优化。

表1 不同优化方法效果对比

直接对栅格型图像进行矢量化转化能减小图像数据量，部分改善绘图区边界的锯齿形，但无法提升图像的平滑度，图像展示效果总体变差。利用等值线平滑的矢量化优化方法不仅能有效消除绘图区边界的锯齿形和内部的马赛克现象，使图像展示效果明显增强，还大幅降低图像数据量，图像的展示具有即时性，图像可随视图范围的缩放同步变化而不失真。该方法的绘图精度灵活可调，可满足不同比例尺和不同精度要求下的噪声地图绘制需求。