0 引　言

超短基线（ultra short baseline，USBL）系统是一类基于水声测量技术的定位系统，该类系统由于基阵尺寸小、操作简单、安装使用及维修方便，十分适用于母船对水下目标定位。目前，该类系统已在深水潜堤施工、潜水员作业、港湾水下目标探测、海峡电缆管道铺设、海洋勘探、海洋科学等领域广泛应用^[1-4]。

USBL定位系统，一般由应答器、声基阵、母船航姿及定位系统组成。其中，应答器固定安装在水下目标上，其余设备固定安装在母船上。超短基线定位系统通过实时测量水下目标相对于母船的相对位置矢量，再结合母船的实时位置信息，实现对水下目标的精确定位。在实际应用时，声基阵相对于母船的安装偏差是影响该类系统定位精度的一个重要因素^[5-6]。一般来说，1°的首向安装偏角将会产生1.7%斜距的定位误差^[7]，因此，在使用该类系统进行定位之前，必须先对声基阵安装偏差进行校准。

解决超短基线定位系统安装偏差校准问题的基本方法是最小二乘法^[8]，即通过两组一一对应的观测数据来拟合参数。由于该校准问题的模型是非线性的，因此通常运用循环迭代的思想来求解^[9-12]。唐秋华等^[9]提出了基于高斯-牛顿法的安装偏差校准算法，该方法的有效性在实际工程应用中不断得到验证^[10-11]。杨保国等^[12]比较了高斯-牛顿法、牛顿法、最速下降法等各类循环算法的性能，提出了将循环算法与非循环算法相结合的安装偏差求解算法。少数学者还根据一些特定应用背景提出了相应的校准方法，江林等^[13]根据工程船施工作业特殊性提出了基于太阳时角法的静态校准方法，Zheng等^[14]针对深海长程超短基线系统提出了基于声线跟踪理论的校准方法。

当声速测量值存在固定测量误差时，一般需要通过规划特定的测线轨迹来消除该测量误差的影响^[15-16]。这种处理方式相对费时费力，并且还需要母船循迹航行。针对上述应用背景条件，本文引入带约束的最小二乘法进行安装偏差校准，通过在模型上添加适当的系数矩阵，提出一种在声速测量误差条件下的校准方法，并通过仿真验证了本文方法的有效性。

1 超短基线的基本测量原理

设 $p_{AB}^u = {[{x_a},\; \; \; {y_a},\; \; \; {z_a}]^{\rm{T}}}$ 为声基阵中心A点至应答器中心B点的相对位置矢量（见图1），根据超短基线测量原理可知^[4]：

$\left\{ {\begin{array}{*{20}{c}} {{x_a} = \dfrac{c}{{2{\text{π}} df}}l{\phi _x},}\\ {{y_a} = \dfrac{c}{{2{\text{π}} df}}l{\phi _y},}\\ {{z_a} = \sqrt {{l^2} - {{\left(\dfrac{{cl{\phi _x}}}{{2{\text{π}} df}}\right)}^2} - {{\left(\dfrac{{cl{\phi _y}}}{{2{\text{π}} df}}\right)}^2}} } {\text。}\end{array}} \right. $

(1)

其中，坐标 ${x_a},\; \; \; {y_a},\; \; \; {z_a}$ 为相对位置矢量在超短基线声基阵坐标系（u系）下的投影，l为斜距（这里 $l = 0.5cT$ ，T为声波往返时间）， ${\varphi _x}$ 与 ${\varphi _y}$ 分别为声基阵端x轴、y轴上的2组水听器接收信号的相位差，关于声基阵坐标系的定义可参见文献[4]。另外，d为声基阵对角线长度，f为声波频率，c为平均声速。

2 超短基线安装偏差的校准方法 2.1 校准问题及其一般求解方法

超短基线定位系统的校准工作主要是完成对声基阵安装偏差的估计，其待估计的参数包括声基阵相对于母船基准位置的安装角度偏差和声基阵相对于母船基准位置的安装距离偏差。

超短基线声基阵安装偏差一般采用高斯-牛顿法进行校准。该方法需要2组测量数据，一组是参考基准值，另一组为测量值，这2组数据一般会拟定存在某种对应关系（即数学模型），将这2组数据按照最小二乘准则，完成对模型参数的最优估计。

就超短基线校准问题来说，参考基准值为母船基准位置至应答器的位置矢量 $p_{OB}^b$ ，测量值为USBL相对位置矢量 $p_{AB}^u$ ，两者之间理论上存在如下变换关系：

$p_{OB}^b = p_{OA}^b + {{C}}_u^bp_{AB}^u{\text{。}}$

(2)

其中， ${{C}}_u^b$ 为u系至船体坐标系（b系）的坐标变换矩阵， $p_{OA}^b$ 为母船基准位置至声基阵中心的位置矢量，关于船体坐标系的定义可参见文献[10]。

这里，参考基准值一般由下式求得：

$p_{OB}^b = C_n^b(p_{NB}^n - p_{NO}^n){\text{。}}$

(3)

其中： $p_{{\rm{N}}B}^{{n}}$ 为应答器在北东地坐标系（n系）下的位置矢量，其坐标值一般通过空间测距交汇法实现精确标定； $p_{{\rm{N}}O}^n$ 为母船在n系下的位置矢量，其坐标值一般根据母船差分定位系统（Differential Global Positioning System，DGPS）的导航参数转换得出。关于北东地坐标系的含义可参见文献[17]。

另一方面，参考基准值可根据测量值通过下式进行估计：

$\hat p_{OB}^b = \hat p_{OA}^b + \hat {{C}}_u^bp_{AB}^u {\text。}$

(4)

设估计偏差为 ${e_i} = {(p_{OB}^b)_i} - {(\hat p_{OB}^b)_i}(i=1,\;\cdots,\;N)$ ，若选取偏差平方和的算术平均值为优化目标，则可应用非线性最小二乘估计法求解 $\hat C_u^b$ 和 $\hat p_{OA}^b$ 。

2.2 基于带约束最小二乘的校准方法

1）带旋转矩阵约束的最小二乘估计法

由参考文献[18]可知，关于带旋转矩阵约束的最小二乘估计法存在如下结论：

引理1

假设m维空间的两组矢量集 $X = \{ {x_1},\cdots,{x_N}\} $ 和 $Y = \{ {y_1},\cdots,{y_N}\} $ ，满足关系 ${y_i} = t + c{{R}}{x_i}\;\;(i = 1,\cdots,N)$ ，其中R为旋转矩阵、t为平移量、c为尺度因子，且令

$\left\{ {\begin{array}{*{20}{l}} {{\mu _x} = \dfrac{1}{N}\displaystyle\sum\limits_{i = 1}^N {{x_i}} } {\text，}\\ {{\mu _y} = \dfrac{1}{N}\displaystyle\sum\limits_{i = 1}^N {{y_i}} } {\text，}\\ {\sigma _x^2 = \dfrac{1}{N}\displaystyle\sum\limits_{i = 1}^N {{{\left\| {{x_i} - {\mu _x}} \right\|}^2}} } {\text，}\\ {\sigma _y^2 = \dfrac{1}{N}\displaystyle\sum\limits_{i = 1}^N {{{\left\| {{y_i} - {\mu _y}} \right\|}^2}} }{\text，}\\ {{\displaystyle\sum _{xy}} = \dfrac{1}{N}\displaystyle\sum\limits_{i = 1}^N {({y_i} - {\mu _y})} {{({x_i} - {\mu _x})}^{\rm{T}}}} {\text，} \end{array}} \right.$

(5)

并且设 ${\displaystyle\sum _{xy}}$ 的奇异值分解为：

$\left\{ {\begin{array}{*{20}{l}} {{\displaystyle\sum _{xy}} = UD{V^{\rm{T}}}{\text，}\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;} \\ {D\; = {\rm{diag}}({d_i}),\;\;\;{d_1} \geqslant \cdots \geqslant {d_m} \geqslant 0} {\text，} \end{array}} \right.$

(6)

则关于变量R，t，c的均方误差函数为：

${e^2}({{R}},t,c) = \dfrac{1}{N}\sum\limits_{i = 1}^N {{{\left\| {{y_i} - (c{{R}}{x_i} + t)} \right\|}^2}}{\text，} $

(7)

最小值为：

${\varepsilon ^2} = \sigma _y^2 - {{tr{{(DS)}^2}} / {\sigma _x^2}}{\text，}$

(8)

并且当 $rank\left({\displaystyle\sum _{xy}}\right) \geqslant m - 1$ 时，变量R，t，c有最优解：

$\left\{ {\begin{array}{*{20}{c}} {{{R}} = US{V^{\rm{T}}}{\text，}\;\;\;\;\;\;\;\;\;} \\ {t = {\mu _y} - c{{R}}{\mu _x}{\text，}\;\;} \\ {c = {{tr{{(DS)}^2}} / {\sigma _x^2}}} {\text{。}} \end{array}} \right.$

(9)

其中：

$S = \left\{ {\begin{array}{*{20}{c}} {I\;,\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\det \left({\displaystyle\sum _{xy}}\right) \geqslant 0} {\text，}\\ {{\rm{diag}}(1,1, \cdots ,1, - 1)\;,\;\;\det \left({\displaystyle\sum _{xy}}\right) < 0}{\text。} \end{array}} \right.$

(10)

为方便起见，以下将该方法简称为带约束的最小二乘法。事实上，该方法已在水下航行器声学多普勒计程仪安装偏差的在线校准问题上得到了实际应用^[19-20]。

2）安装偏差的校准思路

根据引理1可知，要估计2组矢量之间的变换参数，先要知道2组矢量的观测值以及矢量之间的对应关系或者模型，然后才能应用该引理求解。

在超短基线系统校准问题中，具体求解思路可总结为：

首先，超短基线量测矢量组 $\{ p_{AB,i}^u\} (i = 1,\cdots,N)$ ，其中 $p_{AB,i}^u$ 表示 $p_{AB}^u$ 的第i个观测值，可通过超短基线直接测量得到；

其次，参考基准矢量组 $\bigg\{ p_{OB,i}^b\bigg\} (i = 1,\cdots,N)$ ，其中 $p_{OB,i}^b$ 表示 $p_{OB}^b$ 的第i个观测值，需要根据公式 $p_{OB}^b = $ $ C_n^b(p_{NB}^n - p_{NO}^n)$ 计算得到；

然后，将2组观测矢量之间的对应关系设定为 $p_{OB,i}^b = \hat p_{OA}^b + \hat C_u^bp_{AB,i}^u\;\;(i = 1,\cdots,N)$ ；

最后，定义 $X = \bigg\{ p_{AB,1}^u,\cdots,p_{AB,N}^u\bigg\}$ 和 $Y = \bigg\{ p_{OB,1}^b,\cdots, $ $ p_{OB,N}^b\bigg\}$ ，直接应用引理1求 $\hat C_u^b$ 和 $\hat p_{OA}^b$ 。

2.3 声速测量误差条件下的处理方法

当声速测量存在固定测量偏差（以下简称固定声速误差）时，2组观测值 $\bigg\{ \tilde p_{AB,i}^u\bigg\} (i = 1,\cdots,N)$ 和 $\bigg\{ \tilde p_{OB,i}^b\bigg\} (i = $ $ 1,\cdots,N)$ 将不满足模型 $\tilde p_{OB,i}^b = \hat p_{OA}^b + k\hat C_u^b\tilde p_{AB,i}^u\;(i = 1, \cdots , $ $ N)$ ，若直接应用带约束最小二乘法求解 $\hat C_u^b$ 和 $\hat p_{OA}^b$ ，由于模型存在一定的误差，因而参数估计将会有一定的估计偏差。

下面通过分析固定声速误差对超短基线相对位置矢量测量值的影响，给出在该条件下的校准方法。

1）固定声速误差的影响分析

设 ${p^*}= {[{x^*},\;\; {y^*},\;\; {z^*}]^{\rm{T}}}$ ， $\tilde p = {[\tilde x,\;\; \tilde y,\;\; \tilde z]^{\rm{T}}}$ ， $\Delta p = [\Delta x,\;\; \Delta y,\;\; $ $ \Delta z]^{\rm{T}}$ 分别为相对位置矢量 $p_{AB}^u$ 的真值、测量值和测量误差值， ${c^*}$ ， $\tilde c$ 和 $\Delta c$ 分别为声速 $c$ 的真值、测量值和测量误差值。因此， $p_{AB}^u$ 各分量的真值与测量值的比值（以下简称真测比）为：

$ \left\{ {\begin{array}{*{20}{l}} {\dfrac{{{x^*}}}{{\tilde x}} = \dfrac{{\tilde x - \left( {{{\left. {\dfrac{{\partial p_{AB}^u}}{{\partial c}}(1)} \right|}_{c = \tilde c}}} \right)\Delta c}}{{\tilde x}} = 1 - \dfrac{{2 \cdot \Delta c}}{{\tilde c}}} {\text，}\\ {\dfrac{{{y^*}}}{{\tilde y}} = \dfrac{{\tilde y - \left( {{{\left. {\dfrac{{\partial p_{AB}^u}}{{\partial c}}(2)} \right|}_{c = \tilde c}}} \right)\Delta c}}{{\tilde y}} = 1 - \dfrac{{2 \cdot \Delta c}}{{\tilde c}}} {\text，}\\ {\dfrac{{{z^*}}}{{\tilde z}} = \dfrac{{\tilde z - \left( {{{\left. {\dfrac{{\partial p_{AB}^u}}{{\partial c}}(3)} \right|}_{c = \tilde c}}} \right)\Delta c}}{{\tilde z}} = 1 - \left( {1 - \dfrac{1}{{2{{\cos }^2}{\theta _z}}}} \right)\dfrac{{2 \cdot \Delta c}}{{\tilde c}}} {\text{。}} \end{array} } \right.$

(11)

其中， $\dfrac{{\partial p_{AB}^u}}{{\partial c}}(i)$ 表示矢量 $\dfrac{{\partial p_{AB}^u}}{{\partial c}}$ 的第 $i$ 个分量 $(i = 1,2,3)$ ， ${\theta _z}$ 是声线与z轴的夹角，且 ${\theta _z} = \arctan \left( {\dfrac{{\sqrt {{{\tilde x}^2} + {{\tilde y}^2}} }}{{\tilde z}}} \right)$ 。

若将 ${{f}}(x) = \left| {\dfrac{{\tilde x - {x^ * }}}{{{x^ * }}}} \right|$ 定义为影响函数，则可得固定声速误差对相对位置矢量测量值的x轴和y轴分量的影响是一致的，且为固定值；而对z轴分量的影响是随 ${\theta _z}$ 而变化的，显然当 ${0^ \circ } \leqslant {\theta _z} \leqslant {45^ \circ }$ 时，固定声速误差对z轴分量的影响会比x轴和y轴分量的影响要小，且在 ${\theta _z}{\rm{ = }}{45^ \circ }$ 时取极值；当 ${45^ \circ } \leqslant {\theta _z} < {90^ \circ }$ 时，分析较为复杂，但能够证明，若声速测量误差在±1%以内，则当 ${45^ \circ } \leqslant {\theta _z} \leqslant {59.7^ \circ }$ 时，固定声速误差对z轴分量的影响不超过对x轴和y轴分量的影响；当 ${60.3^ \circ } \leqslant {\theta _z} < {90^ \circ }$ 时，固定声速误差对z轴分量的影响不小于对x轴、y轴分量的影响；而当 ${59.7^ \circ } < {\theta _z} < {60.3^ \circ }$ 时，该影响大小与具体的声速测量误差有关，须根据具体误差来判定。

2）处理方法及校准思路

由于固定声速误差对相对位置矢量测量值的x轴和y轴分量的影响是一致的，故考虑将z轴分量的测量值乘以系数 $m$ ，使得固定声速误差对相对位置矢量各分量的影响一致。

设 $\dfrac{{{z^*}}}{{m \cdot \tilde z}} = \dfrac{{{x^*}}}{{\tilde x}} = \dfrac{{{y^*}}}{{\tilde y}}$ ，则比例系数为：

$m = \dfrac{{1 - \left( {1 - \dfrac{1}{{2{{\cos }^2}{\theta _z}}}} \right)\dfrac{{2 \cdot \Delta c}}{{\tilde c}}}}{{1 - \dfrac{{2 \cdot \Delta c}}{{\tilde c}}}}{\text，}$

(12)

显然该系数是时变的，其实时估计值为：

$\hat m = \dfrac{{1 - \left( {1 - \dfrac{1}{{2{{\cos }^2}{\theta _z}}}} \right)\dfrac{{2 \cdot \Delta \hat c}}{{\tilde c}}}}{{1 - \dfrac{{2 \cdot \Delta \hat c}}{{\tilde c}}}}{\text。}$

(13)

其中， $\Delta \hat c$ 为固定声速误差的估计值。

从而，当声速测定存在固定偏差时，可以将两组观测值 $\bigg\{ \tilde p_{AB,i}^u\bigg\} (i = 1,\cdots,N)$ 和 $\bigg\{ \tilde p_{OB,i}^b\bigg\} (i = 1,\cdots,N)$ 的关系模型修正为 $\tilde p_{OB,i}^b = \hat p_{OA}^b + k\hat C_u^b\Lambda \tilde p_{AB,i}^u\;(i = 1, \cdots ,N)$ ，其中 $\hat \Lambda = {\rm{diag}}(1,\;1,\;\hat m)$ 。

在实际求解时，令 $p_i^u = \hat \Lambda \tilde p_{AB,i}^u$ ，则 $\bigg\{ p_i^u\bigg\} (i = 1,\cdots,N)$ 和 $\bigg\{ \tilde p_{OB,i}^b\bigg\} (i = 1,\cdots,N)$ 满足 $\tilde p_{OB,i}^b = \hat p_{OA}^b + k\hat C_u^bp_i^u\;(i = 1, \cdots , $ $ N)$ 。此时，可应用带约束最小二乘法求解 $\hat C_u^b\!\!$ 和 $\hat p_{OA}^b$ 。

因此，声速测量误差条件下的超短基线安装偏差校准步骤可归纳为：

1）获取2组观测值 $\bigg\{ \tilde p_{AB,i}^u\bigg\} (i = 1,\cdots,N)$ 和 $\bigg\{ \tilde p_{OB,i}^b\bigg\} (i = $ $ 1,\cdots,N)$ ，建立模型 $\tilde p_{OB,i}^b = \hat p_{OA}^b + {k_1}\hat C_u^b\tilde p_{AB,i}^u\;(i = 1, \cdots , N)$ ，应用带约束最小二乘法求解尺度因子 ${k_1}$ ；

2）根据尺度因子 ${k_1}$ 和声速测量值 $\tilde c$ ，利用公式 $\Delta \hat c = \dfrac{{\tilde c}}{2}\left( {1 - {k_1}} \right)$ 估算声速误差 $\Delta \hat c$ ；

3）将 $\Delta \hat c$ 代入式（7），估算比例系数 $\hat m$ ；

4）根据 $\hat m$ 计算 $\bigg\{ p_i^u\bigg\} (i = 1,\cdots,N)$ ；

5）针对矢量集 $\bigg\{ p_i^u\bigg\} (i = 1,\cdots,N)$ 和矢量集 $\bigg\{ \tilde p_{OB,i}^b\bigg\} (i = $ $ 1,\cdots,N)$ ，建立模型 $\tilde p_{OB,i}^b = \hat p_{OA}^b + {k_2}\hat C_u^bp_i^u\;(i = 1, \cdots ,N)$ ，应用带约束最小二乘法求解 $\hat C_u^b$ 和 $\hat p_{OA}^b$ 。

3 仿真分析

设母船在水面航行，初始点位置为n系的坐标原点，应答器三维坐标位置为（0 m，20 m，20 m），母船绕应答器航行的轨迹如图2所示。假设DGPS和USBL的数据采样频率均为1 Hz，数据采样总数为400对。

假设安装角度偏差为（3°，2°，5°），安装位置偏差为（–1 m，0 m，1 m）；设DGPS的随机测量误差服从高斯分布N（0，0.1²），USBL相对位置矢量的随机测量误差服从N（0，0.1²）和N（0，r²）的高斯混合分布，其中后者以0.05的概率出现，这里误差标准差的单位为m；设定声速真值为1500 m/s。若将r≥0.3作为存在测量粗差的判据（即认为在r=0.1时无测量粗差），则可以通过调节r来改变测量粗差的大小。

假设高斯-牛顿法的角度偏差、角度修正量、距离偏差、距离修正量的迭代初值均置零，设定迭代次数为8次，设4对样本点分别为第100对、第150对、第250对和第300对。

为方便起见，下面将存在固定声速误差情形下的校准方法称为抗固定声速误差法。

1)仿真实验1

本实验为抗粗差能力测试，假设不存在固定声速测量误差。

在抗粗差能力方面，经过多次仿真后发现，当相对位置矢量无测量粗差（对应r=0.1）时，高斯-牛顿法、带约束最小二乘法和抗固定声速误差法的校准精度相当；当相对位置矢量存在测量粗差（对应r=1），特别是在上述4对样本点中含有粗差时，高斯-牛顿法的校准结果将出现较大偏差，如表1所示。然而，粗差的出现对带约束最小二乘法和抗固定声速误差法的校准结果的影响相对较小。可见，带约束最小二乘法和抗固定声速误差法都具有一定的抗粗差能力。

2)仿真实验2

本实验为抗固定声速测量误差能力测试，假设USBL相对位置矢量无测量粗差。

在抗固定声速测量误差能力方面，考虑到固定声速误差会使校准结果产生系统误差，为研究这种系统误差，特意针对不同的声速误差，分别进行100次仿真模拟，然后将校准结果求平均，得到不同声速测量误差条件下的校准结果。

首先，根据仿真结果可知，在不存在固定声速测量误差的条件下，3种校准方法的校准精度相当；当出现10 m/s（或–10 m/s）的固定声速误差时，高斯-牛顿法和带约束最小二乘法的校准结果均受到一定的影响，然而抗固定声速误差法的校准结果基本不受声速测量误差的影响，如表2所示。

其次，若将角度校准的误差定义为角度校准结果向量与安装角度偏差向量的矢量差向量的模，将距离校准的误差定义为距离校准结果向量与安装距离偏差向量的矢量差向量的模，将角度校准的误差和距离校准的误差统称为校准误差，则根据仿真结果可知，当出现一定的固定声速误差时（见图3），固定声速误差对高斯-牛顿法和带约束最小二乘法的影响相对较大，校准误差随着声速测量误差的增大而逐渐增大，而固定声速误差对抗固定声速误差法的影响相对较小。

由此可见，抗固定声速误差法具有一定的抗固定声速测量误差的能力。

4 结　语

1)针对超短基线安装偏差的校准问题，提出应用带约束的最小二乘法来进行校准。仿真结果表明，该校准方法适用于求解超短基线安装偏差校准问题，并且具有一定的抗粗差能力。

2)在上述校准方法的基础上。仿真结果表明，该方法具有一定的抗粗差能力和抗固定声速测量误差能力，因此，该方法适用于存在固定声速测量误差的情况。

本文的校准方法是以声线为直线传播为前提的，适用于声速变化范围较小的应用场景（比如浅水区）。若用于声速变化范围较大的应用场景（如深水区甚至深海），且声线的初始掠射角偏小，则会导致声线严重弯曲，本文方法不再适用。这里有一种应对思路可供参考：尽量缩小母船的绕行半径，进而缩小与应答器的水平距离，以增大初始掠射角，降低声线弯曲程度，使声线接近直线，以便运用本文方法求解。