CN108469620A - 适用于辐射沙脊群浅水海域的水下地形测量方法 - Google Patents

Abstract

本发明涉及一种适用于辐射沙脊群浅水海域的水下地形测量方法，设置测深平台、PPK基准站和移动站；采集定位、导航数据，以及由测深仪测得的测深数据和由GNSS接收机测得的天线高数据；对获取的数据进行基本处理和数据精处理；采用特征点匹配方案或断面整体平移方案校正时延影响；以时间点对应进行测深数据和GNSS接收机天线高数据匹配，得到该海域地形椭球高数据；通过基准转换获得该海域地形水深数据。该方法应用于浅水滩宽的测量环境，避免对潮汐资料的依赖，具有较高的测量精度和测量效率。

Description

适用于辐射沙脊群浅水海域的水下地形测量方法

技术领域

本发明属于海底地形测量领域，特别是一种适用于辐射沙脊群浅水海域水下地形测量方法。

背景技术

辐射沙脊群是呈辐射状分布的出露于海面以上的沙洲与隐伏于海面以下的沙脊，以及其潮流通道的总称。辐射沙脊海域历来是水深测量的难点区域。究其原因，可以概括为以下几点：

(1)水浅滩宽，人船难至

辐射沙脊海域最显著的特点是滩宽、水浅，潮沟遍布。这样的环境既不便使用普通作业船施测，也不便采用人工跑滩的方法。很多时候，往往需要‘赶潮’作业，施工难度极大。另外，独特的地形，甚至会产生‘怪潮’等灾害现象，施工危险性也很高。

(2)水位改正困难

水位改正是一般水下地形测量作业流程的关键一环，直接关系到测量成果的精度。在辐射沙脊海域水深测量中，水位改正有二个方面困难，其一、潮位变化大，改正量大，水深比重更大；其二，验潮站布设、高程引测非常困难。

(3)测深精度低

主要有三个方面的影响。第一，测量平台小。浅水环境一般选用较小的测量船，一者防止吃水过深而搁浅，同时对于大比例尺测量密集的测线，小型船只也便于掉头及定位。如在江苏沿岸，一些小于5米的极浅水区，甚至会选用吃水 0.5m左右的泡沫筏子进行测量。但此类测量小平台稳定性差，受海况影响大，同时空间有限无法装配过多的硬件设备，直接影响测量信号的稳定性。第二，波浪影响大。波浪在浅水区变形及发生破碎，海面的不规则起伏，严重影响测量平台平稳运行，主要表现为平台高频的三维姿态改变。第三，障碍物多。近岸浅水区同样是海上养殖、捕捞的密集区。部分正位于测线上的水面、水下各类网具迫使船只随时改变航向甚至终止测量，测量连贯性差，测线容易不规则。

(4)多种新方法不适用

在学术研究层面，已有激光测深、遥感反演水深等手段应用于浅水区水深获得。但是，在本区域内，受制于水体浑浊等原因，这些手段测量精度仍比较低，还不能大范围的投入工程化应用。

(5)基准不一

海陆基准统一、海区内基准精化，涉及到基面的确定及高程联测。同以上第 2点，该项工作也需要建立在大量验潮站及长期的验潮资料的基础上，特别是对水文特复杂的本区域。但现实正好相反，江苏沿岸验潮站数量很少，长期验潮站则更少。各项工程的开展，通常只依托于一些临时验潮资料，资料周期短，精度低。

发明内容

本发明要解决的技术问题是提供一种适用于辐射沙脊群浅水海域水下地形测量方法，该方法应用于浅水滩宽的测量环境，避免对潮汐资料的依赖，具有较高的测量精度和测量效率。

为解决以上技术问题，本发明采用的以下技术方案：一种适用于辐射沙脊群浅水海域的水下地形测量方法，包括：

步骤1，系统安装；

设置测深平台

测深仪换能器悬挂在舷侧，测量实时导航采用信标机，信标机安装于换能器杆上部或船顶；

设置PPK基准站和移动站

基准站设置在已知点上，移动站GNSS接收机直接安装于测深杆正上方；

步骤2，采集数据；

采集定位、导航数据，以及由测深仪测得的测深数据和由GNSS接收机测得的天线高数据；

步骤3，数据处理

包括数据基本处理和数据精处理；

所述的数据基本处理，包括测深数据基本处理和PPK方法解算GNSS接收机天线高；

所述的数据精处理，包括校正时延影响、校正GNSS接收机天线高和精处理测深数据；

校正时延影响，采用特征点匹配方案或断面整体平移方案；

所述的特征点匹配方案是在试航阶段选取一些特征水域，选取路线，并对每条测线以不同的速度分别进行往、返测量；首先对所获取的往返断面数据进行测深数据处理，然后在所获取的往返断面数据中选择特征点对数据；通过式(5.1) 计算出该点对时延的量值Δt_k，L_k为同一特征点往返断面中的距离，往返速度分别为v1和v2；

Δt_k＝L_k/(v₁+v₂) (5.1)

逐一对上述特征点对进行时延计算；对所得的时延量取算术平均值即为最终延时量Δt，如式(5.2)；

所述的断面整体平移方案是通过对往返断面进行平移，并判断其最大相似性即可获得最佳匹配效果下的平移量，根据平移量以及船只运行速度，最终计算出时延，两个断面的相关系数可用式(5.3)表示：

式(5.3)表明，存在两个序列h^A和h^B，当两个序列完全一样时，则相关系数R为1；当两个断面不存在相似性时，相关系数R为0；

平移时，以其中一个断面h^A为基准，对另一个断面h^B进行平移，平移步长设置为固定距离d，d尽可能设置的小；连续移动的同时，根据式(5.3)即可以计算出一些列相关系数R，当R最大时，说明二者相关系数最大，此时对应的平移量d就是可以认为是时延引起的断面滞后量；若往返断面测量中的平均船速分别为和则系统延时Δt为：

采用断面整体平移匹配方案时，以高频率采集相应的断面数据；

所述校正GNSS接收机天线高是利用模糊算法，取4阶巴特沃斯滤波，截止频率取0.01HZ；

所述精处理测深数据采用巴特沃斯低通滤波方案或卡尔曼递归滤波方案；

所述巴特沃斯低通滤波方案是取4阶巴特沃斯滤波，截止频率取0.015HZ～0.02HZ之间；

所述卡尔曼递归滤波方案是将一维水深值视为一个离散控制过程的系统，简化建立线性随机微分方程，取d_(t)＝d_(t-1)+w_(t)，式中d(t)是t时刻的系统状态，意即用t-1时刻的系统状态表示t时刻的系统状态，w(t)为系统过程噪声；同时将t时刻系统的测量值D(t)表示为：D_(t)＝d_(t)+v_(t)，v(t)为测量误差；取w(t)、 v(t)的均方差为q和r；

在此前提下，对一维水深测量值，建立卡尔曼预测及估算方程如下：

预测过程：

d_(t|t-1)＝d_(t-1|t-1)认为水深状态前、后两点间连续且相等；

p_(t|t-1)＝p_(t-1|t-1)+q P为d对应的均方差。

估算过程：

d_(t|t)＝d_(t|t-1)+K_(t)(D_(t)-d_(t|t-1))最优估计水深值；

K_(t)＝p_(t|t-1)/(p_(t|t-1)+r)K为卡尔曼增益；

p_(t|t)＝(1-K_(t))p_(t|t-1)更新均方差，进入循环的下一步计算；

计算中，取q＝0.001～0.0015；r＝0.1～0.15之间；

步骤4，数据匹配；

以时间点对应进行测深数据和GNSS接收机天线高数据匹配，得到该海域地形椭球高数据；

步骤5，基准转换

根据椭球高—85高—平均海平面—深度基准面4级转换体系，获得该海域地形水深数据。

进一步地，步骤2中所述的采集定位、导航数据时，在水深测量定位过程中，所有定位数据均采用计算机自动记录存盘，并定时监测DGPS系统的质量参数，定位系统通过导航软件与测量仪器连接，实现同步定位，所有的定位信息均储存在文件中，文件至少包括事件、位置、测线号；在导航定位的同时，计算机通过设备驱动程序自动向测深仪同步输出注释信息，计算机同时采集测点号、测量时间和测点坐标及定位时的水深，并将它们记录在硬盘上。

进一步地，步骤2中，水深测量数据为数字记录和模拟显示。在导航定位的同时，计算机通过设备驱动程序给测深记录定点打标，将序号、日期、时间、经度、纬度、XY坐标和水深数据一同存入计算机中。现场测量结束后，经过水位改正，可以直接获得图载水深。

进一步地，步骤2中，PPK数据记录独立进行，现场不与测深数据同步；基准站、移动站记录设置相同，均采用静态测量模式；卫星高度角设置为10度，记录原始的rinex数据格式，采集频率均设置为1秒钟一次；PPK测量使用GNSS 双频接收机。

本发明提供PPK海测方案适用于辐射沙脊群海域浅水滩宽的测量环境，可替代传统验潮测量方案，能明显提高工作效率。PPK数据及测深数据分别单独处理后，再以时间点匹配的数据融合方案，原理上逻辑清晰，实践操作有效可行。建立在频谱定量分析基础上的PPK截频滤波方案，有效保留了相对较长周期振动信号，使其与测深处理数据对应。最终测量精度比潮改方案更高。

附图说明

图1显示由时延引起的‘错位’现象。

图2显示特征点对匹配方案时延改正效果图。

图3显示断面整体平移方案时延改正效果图。

图4显示6.8小时PPK高程数据。

图5显示10分钟PPK高程数据。

图6是PPK高程信号频谱图。

图7是PPK高程信号滤波效果图。

图8是实验中测线商用软件处理后测点信号图。

图9是巴特沃斯低通滤波效果图(4阶，截频0.015HZ)。

图10是卡尔曼滤波后效果图(取q＝0.001，r＝0.1)。

图11是巴特沃斯低通滤波效果图(4阶，截频0.015HZ)。

图12是卡尔曼滤波后效果图(取q＝0.001，r＝0.1)。

图13是剧变信号滤波效果(上图：巴特沃斯；下图：卡尔曼)。

图14是精处理前后数据等值线成图对比，

注：平滑线为巴特沃斯方案处理，卡尔曼方案处理效果与其近似。

图15是实测潮位对比图。

具体实施方式

本发明一种典型的实施方式提供一种适用于辐射沙脊群海域的无验潮海底地形测量方法，包括一下步骤：系统安装、采集数据、数据处理、数据匹配和基准转换。

本发明主要为测深部分和PPK观测部分，可用于实现本方明所述方法的主要设备包括：

(1)Sps351信标DGPS系统

GPS信标信号：12通道GPS+SBAS，双通道自动搜索信标信号；信标接收频率范围：283.5KHz-325KHz；有效作业距离：海上500km；定位精度：0.5m (1δ)；生产厂家：美国天宝公司。

(2)海星达IRTK2系统

多星多系统内核；通道数：220；BDS：B1、B2、B3；GPS：L1C/A、L1C、 L2C、L2E、L5；GLONASS：L1C/A、L1P、L2C/A、L2P、L3；定位输出频率：1Hz～50Hz；静态格式支持：GNS、Rinex双格式静态数据；RTK定位精度：平面：±(8+1×10^-6D)mm(D为被测点间距离)，高程：±(15+1× 10^-6D)mm；静态定位精度：平面：±(2.5+1×10^-6D)mm，高程：±(5+1 ×10^-6D)mm；生产厂家：中海达公司。

(3)HY1602型双频回声测深仪

生产厂家：海鹰公司；水深测量范围：0.5m—2000m；精度：1cm±0.1％水深，工作频率：24KHz(低频)和208KHz(高频)。

(4)ODOM DigiBarPro声速剖面仪

生产厂家：美国ODOM公司；声速测量范围：1400—1600m/s；声速分辨率：0.1m/s；声速精度：±0.3m/s。

(5)Aqualogger520PT验潮仪(自容)

生产厂家：英国AQUAtec公司；类型：压阻式；量程：10bar；精度：优于±0.2％满量程；采样：1～255秒或1～255分钟。

(6)测量平台

普通调查船，130吨左右，水深大于5米时使用；泡沫筏子，吃水约50cm，水深小于5米时使用。

步骤一，系统安装

1.设置测深平台

测深仪换能器悬挂在舷侧，测量实时导航采用信标机，信标机安装于换能器杆上部或船顶。

测深仪选择舷侧安装方式，即换能器被悬挂在舷侧，在这种安装方式中，固定换能器用的安装管的尺寸必须较好地保证换能器在水下足够深的位置，这样，舷侧安装的换能器在有浪的情况下就不至于露出水面而影响测量，并且要用钢缆把换能器分别向船前和船后拉紧，固定在船上结实的支架上。首选的换能器安装位置是在船的龙骨附近。这样将使船的纵、横摇角度产生最小的影响。

由于PPK基于后处理，不能用于实时导航。故测量实时导航仍采用信标机，信标机安装于换能器杆上部、或船顶(需位置归算)。

2.设置PPK基准站和移动站

基准站设置在已知点上，移动站GNSS接收机直接安装于测深仪的测深杆正上方。

PPK测量基准站需要设置在已知点上。由于基准站和移动站并不需要实时通信，故只要在一定的基线解算长度之内，基准站并不需要像RTK一样设置在码头或沿岸。可以选择测区更有利的观测位置(即便离岸线有一定的距离)，只需确保规范要求的净空条件及不受电磁信号干扰。本发明中观测基准站距离岸线约7 公里。

为准确量取移动站与测深探头间的准确垂直距离，移动站直接安装于测深杆正上方，同时保证观测净空条件。

步骤二，采集数据；

采集定位、导航数据，以及由测深仪测得的测深数据和由GNSS接收机测得的天线高数据。

1.采集定位、导航数据

定位使用信标机接收机系统。在水深测量定位过程中，所有定位数据均采用计算机自动记录存盘，并定时监测DGPS系统的质量参数(HDOP因子，视角等)，保证定位精度要求。定位系统通过导航软件与测量仪器连接，实现同步定位，所有的定位信息均储存在文件中，文件至少包括事件、位置、测线号。

在导航定位的同时，计算机通过设备驱动程序自动向测深仪同步输出注释信息，计算机同时采集测点号、测量时间和测点坐标及定位时的水深，并将它们记录在硬盘上。导航可使用南方测绘自由行海洋综合导航软件系统。

2.采集测深数据

水深测量使用双频/单频回声测深仪。正式测量工作开始前，对测深仪主要参数进行调试。测深数据为数字记录和模拟显示。在导航定位的同时，计算机通过设备驱动程序给测深记录定点打标，将序号、日期、时间、经度、纬度、XY 坐标和水深数据一同存入计算机中。现场测量结束后，经过水位改正，可以直接获得图载水深。另外，测区声速用声速剖面仪测定。

3.采集天线高数据，即PPK数据记录

PPK数据记录独立进行，现场不与测深数据同步；基准站、移动站记录设置相同，均采用静态测量模式；卫星高度角设置为10度，记录原始的rinex数据格式，采集频率均设置为1秒钟一次；PPK测量使用GNSS双频接收机。

PPK测量属于载波相位后处理差分，需要使用双频接收机，如采用中海带双星H32系统及三星IRTK2系统，实验表明，前者可接收到20颗左右卫星信号，后者能接收卫星数30颗以上，足以满足解算需求。

步骤三，数据处理；包括数据基本处理和数据精处理。

1.数据基本处理

所述的数据基本处理，包括测深数据基本处理和PPK方法解算GNSS接收机天线高。

(1)测深数据基本处理

为保证便捷及通用性，地形测量外业结束之后，一般采用测深设备自带的水深资料后处理软件，如自由行数据采集处理软件，对采集的原始数据进行初步处理，主要包括预处理和水深改正两部分：原始数据预处理主要是通过软件查找异常点，并将其剔除，同时按成图要求的密度筛选水深点；在预处理之后进行水深改正，包括换能器静态吃水改正、动态吃水改正和潮位改正三部分，前两部分一般在测量准备阶段测得并直接输入计算机改正。测深数据处理后处理之水位改正，主要即指潮位该改正，从而将采集的水深数据改正到统一高程基准。

(2)PPK数据解算GNSS接收机天线高

PPK数据解算可采用商业专业解算软件完成。可用于GPS动态后处理的软件很多，如：天宝TGO、TTC、TBC软件，莱卡LGO软件，南方GPS数据处理软件，中海达HGO等综合性后处理软件，以及南方GQuickPostion软件，中海达Hi-PPK 软件等单一功能处理软件等。这些软件在功能模块数量、应用范围、操作流程、参数设置、自动化程度等方面不尽相同。但经实践验证，均能有效的处理PPK 观测数据。

2.数据精处理

所述的数据精处理，包括校正时延影响、校正GNSS接收机天线高和精处理测深数据。

基于PPK方案的水深测量，主要的精度影响因子有四大类，即硬件性能因子、参数及解算因子、匹配时延因子和系统三维姿态影响因子。其中硬件性能因子主要取决于所使用的仪器设备，很难作后期的改善，只要仪器检测合格，其误差是可以被接受的。

GNSS定位本身的误差包括很多方面。与卫星有关的误差包括卫星钟差，卫星轨道误差，卫星天线相位偏差等；与信号传播路径有关的误差包括电离层延迟，对流层延迟，多路径影响等；与接收机有关的误差包括接收机钟差，接收机天线相位中心偏差，周跳及整周模糊度等；其他误差包括相对论影响，地球自转改正，潮汐影响等。消减这些误差的理论方法有很多种，比如使用双差消电离层观测量消除钟差和电离层延迟，使用模型估计对流层延迟，通过天线文件修正天线相位偏差，沿用Track等模块的潮汐估算模型计算潮汐影响等。事实上，任何一款处理软件均加载了这一类处理方案，而加载方案的优劣将直接决定软件最终的计算精度。对于测深部分而言，声速的设定，改正水位的测定等参数设置，直接会影响最终成果精度。正确的处理方式是严格按照规范及时测定这些参数并导入相关软件。

本发明特别针对校正时延影响、校正GNSS接收机天线高和精处理测深数据进行了改进。

(1)校正时延影响

由于测深工作是在动态的环境下进行，当测量船沿测深线航行测深时，必须进行定位。如果定位与测深不同步，将使测深值产生位移，使所测整个海区海底形状产生失真，称之为时延影响。无论是过去还是现在的海道测量方法，也无论定位采用的实时定位或后处理确定位置，由于定位和测深始终分属两个不同的系统(通常定位系统具有光学或电磁波特性，而测深系统具有声学特性，并且前者工作在海面之上，后者工作在海面之下)，所以始终存在不同程度的时延影响。

无论是现场使用信标机的实时定位还是基于PPK后期解算的位置，相应的时延影响值随着卫星数目等多种外围因素的不同而不同，是变化的，将其视为常数，则必然给后续测深计算引入时延偏差，造成定位、测深起算基准与深度不同步，从而也给水下地形测量成果带来较大影响。最为直接的反映即为对同一条测线实施往返测量时，同一目标出现‘错位’现象。

如图1所示，对同一测线进行往返断面测量时，在特征地形所处位置二者均出现错位现象(同一位置存在一定的位移量)。造成种现象的原因在于时延的存在，导致二者均存在滞后性。该位移量为二者速度对时间的积分量值。

校正时延影响，本发明采用特征点匹配方案或断面整体平移方案。这两种方法在单波束测量和多波束测量中均可有效的运用。

所述的特征点匹配方案是在试航阶段选取一些特征水域，选取路线，并对每条测线以不同的速度分别进行往、返测量；首先对所获取的往返断面数据进行测深数据处理，然后在所获取的往返断面数据中选择特征点对数据；通过式(5.1) 计算出该点对时延的量值Δt_k，L_k为同一特征点往返断面中的距离，往返速度分别为v1和v2；

Δt_k＝L_k/(v₁+v₂) (5.1)

逐一对上述特征点对进行时延计算；对所得的时延量取算术平均值即为最终延时量Δt，如式(5.2)；

上述方法在单一特征点对的基础上，通过增加特征点对数较好的解决了传统方法的稳定性，同时该方法简单易实现。特征点对匹配方案时延改正效果参考图2。

所述的断面整体平移方案是通过对往返断面进行平移，并判断其最大相似性即可获得最佳匹配效果下的平移量，根据平移量以及船只运行速度，最终计算出时延，两个断面的相关系数可用式(5.3)表示：

式(5.3)表明，存在两个序列h^A和h^B，当两个序列完全一样时，则相关系数R为1；当两个断面不存在相似性时，相关系数R为0；

平移时，以其中一个断面h^A为基准，对另一个断面h^B进行平移，平移步长设置为固定距离d，d尽可能设置的小；连续移动的同时，根据式(5.3)即可以计算出一些列相关系数R，当R最大时，说明二者相关系数最大，此时对应的平移量d就是可以认为是时延引起的断面滞后量；若往返断面测量中的平均船速分别为和则系统延时Δt为：

采用断面整体平移匹配方案时，以高频率采集相应的断面数据。只有这样才可以实现最大限度的相似性匹配，才可以准确的计算出平移量。由于该方法相当于采用了连续的特征点进行匹配，因此无论是精度还是稳定性，都较特征点方案有更进一步的提高。断面整体平移方案时延改正效果可参考图3。

(2)系统三维姿态影响

整个测量系统三维姿态的改变，对测深会产生最为复杂的影响。相比于传统测深方案，运用PPK方案，除了垂向上(heave)的变化可以在数据上相互抵消外，其余二维(roll、pitch)的变化会对测深数据产生更为复杂的影响，主要包括：

定位中心和测深中心偏移

GPS天线如直接架设在换能器连接杆上方，换能器垂直向下发射测量波束，那么在船体理想平衡状态下，定位中心和测深中心是一致的。但一旦船体发生横摇或纵摇，很明显，两个定位中心马上会不一致，两者的不一致程度与横摇、纵摇的角度以及水深(大角度时)直接相关。

水深量测误差

单波束水深测量的波束足印并非为一个点，而是有一定半径的一个圆，硬件采集足印内的返回信息合成为单一的水深返回值。对于小幅度的，这种设计保证了返回水深的稳定性。但当幅度较大时，所形成的波束倾斜导致水深虚假变深，而且由于上述圆形足印的问题，水深虚假变深程度并不一定完全符合三角定理。

GPS天线与换能器垂距改变

GPS天线与换能器通过换能器杆刚性连接。其平稳状态下的垂距固定不变且在测量前需严格测得。但一旦船体发生纵横摇，两者之间的垂距会相应的改变。不过和水深测量误差不同，只要严格测定三维姿态改变量值，这种垂距变化是可以通过三角定理推算的。

结合以上的分析，要对三维姿态变化引起的综合误差进行改正，最理想的方法，是在测量平台的重心位置加载三维姿态传感器(MRU)，并严格测定船、测深仪、GPS天线等设备的空间三维坐标，构建坐标系间的转换关系，同时根据MRU 测得的参数，推算姿态变化引起测量误差的改正量。但是，这一方法虽数学意义严格，但过程复杂且需加载如MRU等额外的硬件设备，不便于一般施工人员使用。

本发明针对这一情况，从信号频率变化的角度，建立了一种模糊分析策略，虽然不能像空间转换计算模型那样，对误差的分解有具体的认识、具体的改正，同时处理精度也可能会降低。但这种方案实施起来快速便捷，有利于运用到实际的测量工作中。

(1)校正GNSS接收机天线高

PPK观测的原始数据经专业软件解算后，得到逐时的天线高观测值。如图4 所示，为一次实验中历时6.8小时的天线高解算值。

放大后，明显可以发现信号中异常的高频振荡，图5所示，

按照上述分析理论及保证数据间匹配对应，对PPK信号中高频成分同样需要滤除，但截频应控制在一定的范围之内。为定量分析PPK高程信号的频谱分布，引入傅里叶变换，将信号由时域转化为频域，如图6。

由信号频谱图可知，除了低频有效信号外，明显振幅主要集中于0.01～ 0.025HZ左右，及0.1～0.4HZ之间。前者40～100秒左右周期的振动与涌浪周期相近，后者2.5～10秒周期的振动，从实践经验看，与风浪引起的船只高频晃动频率(引起垂向分量变化)相符。这两种振动是需要滤除部分的最主要成分。故此，在该海域的PPK测量数据，取截止频率0.01HZ进行滤波，符合上述高频部分滤除的总体思路。

经测试，同样利用巴斯沃特四阶低通滤波器，能够保证在通频带内有平稳的幅频特性。截频0.01HZ，滤波结果如图7所示。滤波后曲线总体平滑，但保留了局部细节及阶段性变化特点，符合设计预期。

本文研究的模糊处理方案，除了能有效的去除异常信号，保留有效信号外，最大的特点是便于在计算机环境下快速实现。

(2)测深数据精处理，采用巴特沃斯低通滤波方案或卡尔曼递归滤波方案。

在浅水区或高精度要求下，特别是密集取样后，成果信号中仍能明显发现噪声的存在。如图8所示，有振幅约0.2米的高频信号伴随整个测量过程，这类信号与实际地形不符，却符合由于三维姿态不稳定造成的数值变化。在精处理过程中，应该去除这一部分信号。本发明从数值滤波的角度，设计模糊化的处理方案，自动对这类噪声进行快速有效的精处理。

有两种方法可作为合理滤波器选择设计的指导。其一，如图可见，信号中的噪声成分具有高频的特点，真实地形信号应该是相对缓慢的变化，在时间尺度上则表现为低频。基于这一点，可选用合适的低通滤波器去除高频噪声成分；其二，地形是连续的，噪声量值也可估计，即可以建立受噪声干扰的系统模型。加上包含噪声干扰的系统观测量，可以运用卡尔曼滤波器这样的时变线性系统的递归滤波器进行真实信号的估计。以下在研究建立这两种滤波器的基础上，对软件粗处理后数据进行精处理及对比分析。

(1)巴特沃斯低通方案

从频率分析的角度，本文引入巴特沃斯低通滤波方案。这种数字滤波器通频带内的频率响应曲线最平滑的特点符合对水深值的滤波需求。将等时间间隔水深采样值视为随时间变化的一维变量。对该变量作FFT变换后进行频谱分析，可区分出高频信号的集中频率范围。对于巴特沃斯低通滤波而言，需设置阶数及截止频率两个参数。由于测量平台及作业海况的不同，参数的设定对不同的测线应有所调整，尤其是截止频率。在保证真实信号有效保留的前提下，经测试，本文数据精处理取4阶巴特沃斯滤波，截止频率取0.015HZ～0.02HZ之间，对应约50S～ 60S周期的振动。滤波效果如图9所示。

(2)卡尔曼递归方案

从数据递推估计的角度，本发明引入卡尔曼递归滤波方案。所述卡尔曼递归滤波方案是将一维水深值视为一个离散控制过程的系统，简化建立线性随机微分方程，取d_(t)＝d_(t-1)+w_(t)，式中d(t)是t时刻的系统状态，意即用t-1时刻的系统状态表示t时刻的系统状态，w(t)为系统过程噪声；同时将t时刻系统的测量值D(t)表示为：D_(t)＝d_(t)+v_(t)，v(t)为测量误差；取w(t)、v(t)的均方差为 q和r；

在此前提下，对一维水深测量值，建立卡尔曼预测及估算方程如下：

预测过程：

d_(t|t-1)＝d_(t-1|t-1)认为水深状态前、后两点间连续且相等；

p_(t|t-1)＝p_(t-1|t-1)+q P为d对应的均方差。

估算过程：

d_(t|t)＝d_(t|t-1)+K_(t)(D_(t)-d_(t|t-1))最优估计水深值；

K_(t)＝p_(t|t-1)/(p_(t|t-1)+r)K为卡尔曼增益；

p_(t|t)＝(1-K_(t))p_(t|t-1)更新均方差，进入循环的下一步计算；

计算中，取q＝0.001～0.0015；r＝0.1～0.15之间。滤波效果如图10所示。事实证明，这两种滤波方案均非常有效。

步骤4，数据匹配；

以时间点对应进行测深数据和GNSS接收机天线高数据匹配，得到该海域地形椭球高，其基本计算公式是：该海域地形椭球高＝用PPK方法计算的GNSS接收机天线高(T)-常数A-瞬时水深值D，常数A是测深仪与GPS天线间的距离。

步骤5，基准转换

PPK海测方案的实施过程中，涉及到多个垂向参考面。按一般的海上水深测量要求，最终成果的高程基准多用1985国家高程基准或理论深度基准，而GPS 天线直接得到是是椭球高，这就必然需要进行垂直基准的转换。

根据椭球高—85高—平均海平面—深度基准面4级转换体系，获得该海域地形水深数据。

(1)椭球高至85高转换

由椭球高转换至85高在陆地测量中，一般通过解算控制点计算转换参数拟合完成。但是，海上没有控制点，较为可行的处理方式是，通过全球或区域的似大地水准面模型(尽可能选择适用于当地的高精度模型)作两者之间的关系换算。

(2)85高至平均海平面转换

按大地测量学的传统观点，认为平均海平面与大地水准面吻合，但实际两者之间存在差异，其值可达3～4米，此差异称为海面地形。以江苏海域海域为例，这个差值最大仅为3厘米，因此，对于江苏海域的辐射沙脊群海域的测量，可将该值设为一个常数。

(3)平均海平面至深度基准面(理论最低潮面)转换

根据《海道测量规范》(GB12327-1998)，我国将理论最低潮面作为深度基准面，深度基准面的高度从当地平均海面起算，用L值表示(深度基准面在平均海面下的高度，单位cm)。

根据规范，L值的确定需要根据M₂、S₂、N₂、K₂、K₁、O₁、P₁、Q₁、M₄、MS₄、M₆、 S_a、S_sa等13个分潮的调和常数、节点因子根据一组公式求得最小值。前11个分潮的调和常数由30天水位观测资料，用潮汐调和分析法求得。S_a、S_sa分潮的调和常数则以一年的水位观测资料求得，对短期验潮站的S_a、S_sa分潮的调和常数，可采用邻近长期验潮站的S_a、S_sa分潮的调和常数。

要在大范围内密集的设置验潮站是个极费人力物力的事。本发明所述的PPK 海测方案的最大优点是不需要进行潮位观测，如为获得每个测点的深度基准面而仍然需要投置验潮站，则PPK海测方案的优点将无从体现。能否根据历时已知资料，建立区域内的深度基准面模型，是一个值得研究的问题。

海图深度基准面呈点状分布在各地沿岸的长期验潮站处，也就是说按照某种算法求得的海图深度基准面数值仅是对真正曲面形态基准面在特定点的采样。将这些离散的、跳变的基准‘点’转化为连续的基准‘面’，较为通用的数学方法是插值或拟合，两者的前提都要充分搜集区域内已知的海图深度基准点。本发明归纳两种对应的有效方法如下：

(1)Kriging插值方法

Kirigng的数学模型是建立在变异函数理论分析的基础之上,是对有限区域内的区域化变量取值进行无偏最优估计(best linear unbiased estimator,BLUE)的一种拟合方法。这种方法与传统插值方法的不同之处，在于估计元观测样本数值时，不仅考虑待插值点与邻近有观测数据点的空间位置，还考虑了各邻近点之间的位置关系，而且利用已有观测值空间分布的结构特点，使其估计比传统方法更精确，更符合实际，并可以有效避免系统误差产生的“屏蔽效应”。它在地质统计学及图像处理等方面得到了广泛应用。

具体地说，Kriging算法就是对每一采样值分别赋予一定的权系数，再进行加权平均来估计待估值。设待估点的真值进行估计，所用的实验数据是一组离散的某一指标值Z_a(a＝1,2,......,n)，待估点的Kriging估值可以表示为n个数值的线性组合，即Kriging方法就是求出权系数λ_a，使为的无偏估计值，并且估计方差最小。

为确保Kriging插值能有正的Kriging方差，必须根据经验半变异函数来调整一个模型(即一个连续函数或曲线图)。常用的模型包括球状模型、指数模型、高斯模型及其衍生模型等，然后根据来自经验半变异图的点来调整模型。根据前人研究，半变异函数Exponential-Cosine(type1)的拟合效果最佳，各项精度指标均最佳。

(2)多面函数拟合方案

多面函数拟合法的基本思想是：任何一个规则或不规则的连续曲面均可以用 n个规则曲面的叠加来拟合或逼近。该方法是从几何观点出发的一种曲面逼近法，解决根据数据点形成一个平差的数学曲面问题。该方法在大地水准面差距、重力异常、垂线偏差以及地壳变形等方面都取得了较好的拟合效果。

在某一观测区域，每个观测点都可以同该区域中各个已知点分别建立函数关系，该函数称为核函数，它的图形是一个规则数学曲面，将这些规则数学曲面按一定比例叠加起来，就可拟合出任何不规则的曲面，且能达到较好的拟合效果。多面函数在笛卡尔坐标系中的一般形式为：

式中：β_i为待定参数；F(x,y,x_i,y_i)为核函数；(x,y)为待求点坐标，其中心在 (x_i,y_i)处，为已知点坐标，核函数可以有多种表现形式，常用简单的核函数有以下几种：

1、锥面

2、双曲面

3、倒双曲面

4、三次曲面

上述各式中，为内插点到参考点之间的水平距离；C，δ是光滑系数，或称为平滑因子。

设已知海图深度基准面值有m个点，其m×1向量记为l，选取其中n个结点为j＝1,2,3,,n(n≤m)，β＝(β₁,β₂,β₃,,β_n)为n×1向量，当m＞n时，多面函数拟合算法对应的误差方程式为V＝Fβ-l，式中

根据最小二乘法原理可知β＝(F^TF)^-1F^Tl，则任意一点P(x_p,y_p)的海图深度基准面值l_p可表示为：

l_p＝F_pβ＝F_p(F^TF)^-1F^Tl

在选择已知点时，最好选择测区中变化显著的点，这些点能很好描述该区域内各点的分布特征，一般位于最高、最低处以及坡度变化处。多面函数拟合法主要是确定核函数和平滑因子δ的选取，需要不断试验和改进，以切合曲面局部的变化趋势，取得最佳的拟合效果。

研究表明，不同的核函数对平滑因子δ显示出不同的依赖性，对拟合结果的精度也有相应的影响。其中倒双曲面多面函数对平滑因子δ最为敏感，参数选取稍有偏离就可能产生较大误差；而三次曲面多面函数对平滑因子δ的依赖性很弱，可以在比较大的范围内保持一定的拟合精度；对于正双曲面多面函数，当可调参数在一定范围内时，拟合精度随参数的变化而逐渐变化。在各自都取得最佳参数值时，倒双曲面多面函数的拟合精度高于其它3种。

使用多面函数拟合法进行曲面拟合时，平滑因子δ的取值非常重要，在选取时有一定的规律可循。一般若结点间的平均距离增加，则选取的δ值相应也要提高；反之亦然。δ可以取已知数据点间距平均值的倍数，具体要对已知数据点进行反复试算来确定。

以上两种方法是构建深度基准模型的较为通用的方法。研究认为，这两种方法的有效运用有一个前提：即在研究区内搜集到足够多的‘已知点’信息，使点与点之间具有连续性，这是插分或拟合有效、准确的基础。但是，江苏辐射沙脊海域地形复杂，导致潮汐特征复杂，同时，验潮站数量极少(尤其是长期验潮站)，很难保证上述的连续性。故此，提出一种新的思路，可作为进一步研究的方向。

利用现有水动力数模软件，建立区域内网格化潮位预报模型。已知点信息可以作为特征点同化到模型中，提高模型预报精度。以网格点上的模型预报值反推各点上的潮汐调和常数，继而计算各点的理论最低潮位，由这些点构成区域内网格化的深度基准面模型。

本发明提供的思路得到的模型网格密度高、充分考虑了地形等边界条件影响下的水动力特征、对现有已知点的需求较小，便于在辐射沙脊海域实现。

下面提供一些对比结果对本发明方法的技术效果进行说明。

1.近岸浅水区水深数据精处理结果对比

近岸浅水环境对单波束水深测量的影响，有以下几个层次：

第一，外海潮波向近岸传播，在浅水区发生破碎，故近岸处海域，无论大范围天气如何，海面总难平静，海况差。

第二，海况差导致测量平台不稳定。特别的，浅水环境为防吃水过深搁浅及便于大比例尺测线船只掉头，通常选用较小的测量船，但小测量平台更易受海况影响(如在江苏沿岸，一些小于5米的极浅水区，甚至会选用相当轻泡沫筏子进行测量)。小平台的另一个缺点是作业空间小，加载不了过多的辅助性硬件设备。

第三，平台不稳定，导致测量设备不稳定，如测深仪换能器的三维姿态变化剧烈，直接导致测量信号不稳定，影响数据质量。

第四，近岸浅水区同样是海上养殖、捕捞的密集区，障碍物多。水上、水下各类网具除了容易产生假信号，也会迫使船只随时改变航向甚至终止测量，走走停停，测量连贯性差。

经过巴特沃斯及卡尔曼两种滤波方法处理后的信号，明显可以看到高频信号得到了压制，上文分析的约20cm左右的振幅得到削减，削减量可达2/3之多，滤波效果参考图11和图12。同时，两种方案虽略有差别，但均能有效跟踪低频信号，如实反映真实地形趋势。实验证明，即便是剧变地形，精处理信号在去噪同时，同样不会失真，参考图13。

2.交点差对比

测线交叉差是在开放数据处理后，对水深测量最主要的封闭精度指标，也是测量精度最直接的定量反映。

交点差统计按以下公式：

式中：M为重合点水深不符值中误差，单位为米；di为主测线与联络测线在重合点i处的深度不符值，单位为米；n为主测线与联络测线的重合点数。

以下是27个交点处，粗处理结果及两种方案精处理结果的详细对比表。为保证交点准确，测线数据均采用1米间隔密集取样数据。

表1测线交点差对比表(1米间隔取样数据)

对上表数据的分析有以下几点认识：①本实验所用测量数据经过了精确的水位改正(水位改正误差小)，主要误差为水深量测误差。数据精处理后，精度明显提高，反映为均值提高50％，均方值提高近70％。说明，两种精处理方案针对水深量测误差切实有效；②两种方案效果总体一致，但略有细微差异，主要还是取决于相关滤波参数的合理设定；③测区外部10米左右水深处用50吨船作业，对比近处泡沫筏子数据，处理前交点差总体更大，精处理效果相应的更显著。此与船只三维姿态变化(波束角变化)导致在水深处测量误差更大的特点相符。同时也说明两种精处理方案在一定意义上可发挥姿态传感器硬件的作用；④确实有部分交点差不降反升，但范围很小，不影响总体滤波效果。

3.地形图比对

测量精度的提升，同样在地形等值线图上得到了直观的反映。在全部参数一致的网格化处理后。明显可见，粗处理数据所成等值线由一系列小圈闭构成了较粗的等值‘线’，这是大量噪声包裹真实信号的典型特征，不是实际地形的反映。精处理数据所成等值线，平滑顺畅，几无周围散点，体现了地形的缓慢变化而不失美观，参考图14。

研究设计的巴特沃斯滤波及卡尔曼滤波精处理方案，对于近岸浅水区工程水深测量粗处理数据均能有效再提高其精度且数据不失真，总体改善量值在50％以上。最重要的，这两种方案均便于使用计算机语言编程后快速实施，能为一般施工人员接受使用，高效率高精度的服务于各类工程实践。

4.PPK方案与验潮方案对比

为全面了解PPK海测方案在整个辐射沙脊海域的适用性，PPK方案和验潮方案的比对分别针对长、短解算基线两种情况进行。

实验一测区位于大丰港以东，东沙以西海域，测区离岸约10公里，最终基线解算长度在18～21公里之间。近岸工程水深测量主要集中在离岸线20公里内，故探讨在这样的短基线情况下PPK方案的测量精度，有很大的实践意义。

对比分为两个方面。首先将获得PPK信号经上述研究方案处理后，与实测潮位的对比，能够反映PPK成果的总体相符程度，同时对部分滤波后数据与纯潮位曲线的异同有定量的认知。取测区内验潮站10分钟间隔连续验潮数据作对比，对比数据36组，计时长6小时。如图15所示。两组数据除总体几近一致，相符度好。反映的差异部分，均差4.5cm，最大差值10cm，最小差值0cm。以此表明， PPK信号完全反应了水位面的实时变化，有很好的替代作用。其次仍为交点差比对。由于测量比例尺大，主、辅测线密集，故在实验区域内统计得有效交点310 个(1秒取样最近点，不拟合)。PPK方案及验潮方案交点差统计结果如表2。

表2交点差统计对比表(单位：厘米)

	均差	标准差	最大值	最小值
					PPK	10.8	8.2	17	0
验潮	11.4	8.6	24	0

由此可见，PPK测量方案较验潮方案精度还略有提高，尤其是对个别差值较大点的约束。这和研究中采取的基于定量截频的部分滤波方案是直接相关的。

通过实验一的成果分析，总体有以下几点认识：①PPK海测方案适用于南黄海辐射沙脊群海域浅水滩宽的测量环境，可替代传统验潮测量方案，能明显提高工作效率。②PPK数据及测深数据分别单独处理后，再以时间点匹配的数据融合方案，原理上逻辑清晰，实践操作有效可行。③建立在频谱定量分析基础上的PPK截频滤波方案，有效保留了相对较长周期振动信号，使其与测深处理数据对应。最终测量精度比潮改方案更高。

实验三测区位于外磕脚临海基点附近海域，实际解算基线长度约80～90公里。主要验证PPK海测方案在长基线状态下运用的可行性。同时，这一离岸距离基本涵盖了辐射沙脊主要的潮流通道及浅水区域，有显著的实用意义。实验数据采集范围为4.5×4.5公里，按1:5万测图比例尺布设东西向主测线10条，南北向辅测线2条。为准确的作数据比对，测区中心投放了验潮仪做了同步潮位观测 (单站改正)。

全部实验采集数据均按前述方法进行了精细化处理。经数据垂向转换，统一高程基准，对比了采集到的20组交点差数据在潮改方案及PPK方案下的异同。

表3测线交点差对比表

上表数据的对比分析，反映了以下信息：①不同于短基线解算，长基线状态下，PPK方案的测量精度有所降低，偏差主要集中于10cm的级别，低于潮改方案(水位信息非常准确的前提)。但仍能达到测量规范要求。②有效的数据精处理方案对长基线解算，同样效果明显，能大幅提高成果精度。③PPK海测方案完全适用于辐射沙脊广阔海域的地形测量，高精高效，能够全面替代传统潮改方案。

Claims (4)

1.一种适用于辐射沙脊群浅水海域的水下地形测量方法，包括：

步骤1，系统安装；

设置测深平台

测深仪换能器悬挂在舷侧，测量实时导航采用信标机，信标机安装于换能器杆上部或船顶；

设置PPK基准站和移动站

基准站设置在已知点上，移动站GNSS接收机直接安装于测深杆正上方；

步骤2，采集数据；

采集定位、导航数据，以及由测深仪测得的测深数据和由GNSS接收机测得的天线高数据；

步骤3，数据处理

包括数据基本处理和数据精处理；

所述的数据基本处理，包括测深数据基本处理和PPK方法解算GNSS接收机天线高；

所述的数据精处理，包括校正时延影响、校正GNSS接收机天线高和精处理测深数据；

校正时延影响，采用特征点匹配方案或断面整体平移方案；

所述的特征点匹配方案是在试航阶段选取一些特征水域，选取路线，并对每条测线以不同的速度分别进行往、返测量；首先对所获取的往返断面数据进行测深数据处理，然后在所获取的往返断面数据中选择特征点对数据；通过式(5.1)计算出该点对时延的量值Δt_k，L_k为同一特征点往返断面中的距离，往返速度分别为v1和v2；

Δt_k＝L_k/(v₁+v₂) (5.1)

逐一对上述特征点对进行时延计算；对所得的时延量取算术平均值即为最终延时量Δt，如式(5.2)；

所述的断面整体平移方案是通过对往返断面进行平移，并判断其最大相似性即可获得最佳匹配效果下的平移量，根据平移量以及船只运行速度，最终计算出时延，两个断面的相关系数可用式(5.3)表示：

式(5.3)表明，存在两个序列h^A和h^B，当两个序列完全一样时，则相关系数R为1；当两个断面不存在相似性时，相关系数R为0；

平移时，以其中一个断面h^A为基准，对另一个断面h^B进行平移，平移步长设置为固定距离d，d尽可能设置的小；连续移动的同时，根据式(5.3)即可以计算出一些列相关系数R，当R最大时，说明二者相关系数最大，此时对应的平移量d就是可以认为是时延引起的断面滞后量；若往返断面测量中的平均船速分别为和则系统延时Δt为：

采用断面整体平移匹配方案时，以高频率采集相应的断面数据；

所述校正GNSS接收机天线高是取4阶巴特沃斯滤波，截止频率取0.01HZ；

所述精处理测深数据采用巴特沃斯低通滤波方案或卡尔曼递归滤波方案；

所述巴特沃斯低通滤波方案是取4阶巴特沃斯滤波，截止频率取0.015HZ～0.02HZ之间；

所述卡尔曼递归滤波方案是将一维水深值视为一个离散控制过程的系统，简化建立线性随机微分方程，取d_(t)＝d_(t-1)+w_(t)，式中d(t)是t时刻的系统状态，意即用t-1时刻的系统状态表示t时刻的系统状态，w(t)为系统过程噪声；同时将t时刻系统的测量值D(t)表示为：D_(t)＝d_(t)+v_(t)，v(t)为测量误差；取w(t)、v(t)的均方差为q和r；

在此前提下，对一维水深测量值，建立卡尔曼预测及估算方程如下：

预测过程：

d_(t|t-1)＝d_(t-1|t-1)认为水深状态前、后两点间连续且相等；

p_(t|t-1)＝p_(t-1|t-1)+q P为d对应的均方差；

估算过程：

d_(t|t)＝d_(t|t-1)+K_(t)(D_(t)-d_(t|t-1))最优估计水深值；

K_(t)＝p_(t|t-1)/(p_(t|t-1)+r)K为卡尔曼增益；

p_(t|t)＝(1-K_(t))p_(t|t-1)更新均方差，进入循环的下一步计算；

计算中，取q＝0.001～0.0015；r＝0.1～0.15之间；

步骤4，数据匹配；

以时间点对应进行测深数据和GNSS接收机天线高数据匹配，得到该海域地形椭球高数据；

步骤5，基准转换

根据椭球高—85高—平均海平面—深度基准面4级转换体系，获得该海域地形水深数据。

2.根据权利要求1所述的方法，其特征在于：步骤2中所述的采集定位、导航数据时，在水深测量定位过程中，所有定位数据均采用计算机自动记录存盘，并定时监测DGPS系统的质量参数，定位系统通过导航软件与测量仪器连接，实现同步定位，所有的定位信息均储存在文件中，文件至少包括事件、位置、测线号；在导航定位的同时，计算机通过设备驱动程序自动向测深仪同步输出注释信息，计算机同时采集测点号、测量时间和测点坐标及定位时的水深，并将它们记录在硬盘上。

3.根据权利要求1或2所述的方法，其特征在于：步骤2中，水深测量数据为数字记录和模拟显示，在导航定位的同时，计算机通过设备驱动程序给测深记录定点打标，将序号、日期、时间、经度、纬度、XY坐标和水深数据一同存入计算机中。现场测量结束后，经过水位改正，可以直接获得图载水深。

4.根据权利要求3所述的方法，其特征在于：步骤2中，PPK数据记录独立进行，现场不与测深数据同步；基准站、移动站记录设置相同，均采用静态测量模式；卫星高度角设置为10度，记录原始的rinex数据格式，采集频率均设置为1秒钟一次；PPK测量使用GNSS双频接收机。