CN103473812B - 一种基于tps的等深线追踪算法 - Google Patents

Abstract

本发明公开了一种基于TPS的等深线追踪算法，包括下述步骤：S1、导入具有三维坐标的离散点作为控制点；S2、利用TPS对离散点进行曲面拟合，得到一个通过每个离散点的光滑曲面公式；S3、根据离散点的密度和范围构建初始平面格网；S4、计算每个格网节点所在拟合曲面上的高程值，将平面格网转换成曲面格网；S5、根据给定的高程值在曲面格网上追踪等深线。本发明基于TPS的曲面拟合算法能很好地根据采集点模拟水下地形，通过建立一定密度的规则格网可以精细地追踪等深线。又因曲面是光滑的，因此生成的等深线也具有一定的平滑性，节点分布比较紧密和均匀，最后利用张力样条或多次样条拟合就可得到光滑的等深线，有效避免曲线拟合过程中产生凸起或陡变现象。

Description

一种基于TPS的等深线追踪算法

技术领域

本发明涉及等深线绘制的技术领域，特别涉及一种基于TPS的等深线追踪算法。

背景技术

一直以来，等深线的绘制，是通过手扶数字化完成，然而这样的作用方式劳动强度大，周期长。在自动数字化成图的趋势中，出现了比较简单的等深线生成方式，就是通过建立三角网。在这样的过程当中水深数据的质量直接影响到等深线的生成效果，因为在等深线自动连线前缺少了数值内插和拟合的过程。然而对于等深线要求在表达信息准确的前提下，尽量要求等深线的线型要流畅顺滑，也就是对数据拟合有较高的要求。现在常用生成等深线的算法是线性三角网法，连续样点数据间的连线形成三角形，覆盖整个研究区域，所有三角形的边都不相交。随后在三角网上追踪等深线待经过的节点，这些节点就会成为等深线连线的控制节点，这些控制节点相连得到原始的折线形的等深线，最后对这样的原始等深线进行拟合。也就是说在整个过程中只进行了一次拟合，加上三角网上控制节点的数量是随机的，故最终生成的等深线的线型会常出现线型走向突然变化，而且拐弯的幅度大，也就是凹凸变化明显，突兀处无法控制的情况。因此急需一种可以精确追踪等深线的方法。

发明内容

本发明的目的在于克服现有技术的缺点与不足，提供一种基于TPS的等深线追踪算法，可以有效减少等深线走向突变的情况，最终得到顺畅平滑的等深线样式。

本发明的目的通过下述技术方案实现：

本发明一种基于TPS的等深线追踪算法，包括下述步骤：

S1、导入具有三维坐标的离散点作为控制点；

S2、利用TPS对离散点进行曲面拟合，得到一个通过每个离散点的光滑曲面公式；

S3、根据离散点的密度和范围构建初始平面格网；

S4、计算每个格网节点所在拟合曲面上的高程值，将平面格网转换成曲面格网；

S5、根据给定的高程值在曲面格网上追踪等深线。

优选的，步骤S2中，TPS曲面方程如下：

z＝a₁+a₂x+a₃y+Δ

$\mathrm{\Δ}=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{w}_i\·s^2\·\log \left(s\right)$

$s=\sqrt{{\left(x_j{-x}_i\right)}^2+{(y_j-y_i)}^2}$

式中z为点(x,y)处的高程值，a₁、a₂、a₃和w_i为方程系数，Δ为修正值，s为当前控制点与其余控制点的平面距离，n为控制点总点数；

将曲面方程用矩阵的形式表达如下：

由

$\left[\begin{array}{cc}K& P\\ P^T& O\end{array}\right]\·\left[\begin{array}{c}\stackrel{\→}{w}\\ \stackrel{\→}{a}\end{array}\right]=\left[\begin{array}{c}\stackrel{\→}{z}\\ \stackrel{\→}{o}\end{array}\right]$

得

$\left[\begin{array}{c}\stackrel{\→}{w}\\ \stackrel{\→}{a}\end{array}\right]={\left[\begin{array}{cc}K& P\\ P^T& O\end{array}\right]}^T\·\left[\begin{array}{c}\stackrel{\→}{z}\\ \stackrel{\→}{o}\end{array}\right]$

由上式即可求得方程系数，式中和为待求向量，形式如下：

$\stackrel{\→}{w}={\left[\begin{array}{cccc}{w}_1& w_2& \·\·\·& w_n\end{array}\right]}^T,\stackrel{\→}{a}={\left[\begin{array}{ccc}{a}_1& a_2& a_3\end{array}\right]}^T;$

K为n×n矩阵，每个元素由下式计算：

$K_{\mathrm{ij}}=s_{\mathrm{ij}}^2\·\log \left(s_{\mathrm{ij}}\right)+\frac{1}{n^2}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{s}_{\mathrm{ij}}$

其余矩阵和向量如下：

$P=\left[\begin{array}{ccc}1& x_1& y_1\\ 1& x_2& y_2\\ \·& \·& \·\\ \·& \·& \·\\ \·& \·& \·\\ 1& x_n& y_n\end{array}\right],$ $O=\left[\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ 0& 0& 0\end{array}\right],$ $\stackrel{\→}{z}=\left[\begin{array}{c}{z}_1\\ z_2\\ \·\\ \·\\ \·\\ z_n\end{array}\right],$ $\stackrel{\→}{o}=\left[\begin{array}{c}0\\ 0\\ 0\end{array}\right]$

解算过程如下：先根据每个控制点的平面坐标计算矩阵K，再组成矩阵P和向量，最后通过矩阵求逆和相关运算即可求得TPS曲面方程的系数。

优选的，步骤S2中，从TPS曲面方程的解算过程中可知，当控制点越多运算量越大，因此可根据控制点的分布，对其进行分块处理，每块建立一个子格网，具体步骤如下：

S21、将所有控制点的分布区域分成若干个正方块，确保每块里包含一定数量的控制点，若数量不足则合并到邻接的块；

S22、解算每一块的TPS曲面方程，根据方形块的范围建立子格网，利用方程系数计算每个格网节点的高程值；

S23、合并相邻的格网，最终形成一个覆盖整个区域的格网。

优选的，步骤S23中，相邻格网衔接处的节点有两个高程值，对该两个高程值取其平均值。

优选的，步骤S3中，初始平面格网是所有格网节点没有高程值或高程值被设置为一个固定值，由于格网中的格子一般取为方形，要追踪等深线的区域往往不是矩形，因此利用控制点的分布情况来限制格网范围。

优选的，利用控制点的分布情况来限制格网范围的具体方法为：生成控制点的外包多边形或构建三角网，在多边形里面或三角网里面的格网节点设为待求点，在外面的节点设为无效点。

优选的，步骤S4中，通过求解TPS曲面方程的系数后，即可根据公式

z＝a₁+a₂x+a₃y+Δ

计算任意点(x,y)的高程值z，为保证所求高程值的精度，要求点(x,y)要求在该曲面的控制点范围内，即控制点构成的外围多边形内。

优选的，步骤S4，面格网也即初始格网，当求解出曲面方程后，即可根据平面坐标(x,y)计算对应的高程值z，由于曲面方程覆盖整个格网所在的范围，每个格网节点的平面坐标均可由格网的起点坐标和节点所在的行列号计算得到，因此只要循环格网中每个节点，逐一计算对应的高程，格网就沿着曲面分布，也就成为曲面格网。

优选的，步骤S5中，在曲面格网上追踪等深线的过程如下：设等深线的高程值为H，相邻两个格网节点的高程分别为H_i和H_j，若H的大小在H_i和H_j之间，则表明这两个节点的连线与该等深线有相交点，通过高差的比例关系求出相交点的坐标并保存起来，通过遍历整个格网得到该等深线与格网的所有交点，由这些交点连成的线即为等深线。

本发明相对于现有技术具有如下的优点及效果：

1、本发明基于TPS的曲面拟合算法能很好地根据采集点模拟水下地形，通过建立一定密度的规则格网可以十分精细地追踪等深线。又因曲面是光滑的，因此生成的等深线也具有一定的平滑性，最后利用张力样条或B次样条拟合就可得到光滑的等深线。也就是经过两次的拟合，第一次对整个面进行了拟合计算，在较平缓的拟合面中求得原始等深线，然后对线进行二次拟合计算，故可以大大减少等深线走向突变的情况，最终得到顺畅平滑的等深线样式。

2、本发明利用控制点的分布情况来限制格网范围，如生成控制点的外包多边形或构建三角网，在多边形里面或三角网里面的格网节点设为待求点，在外面的节点设为无效点。通过这样种方法来约束等深线的追踪范围，且在利用TPS曲面方程计算格网节点的过程中，减少了不必要节点的计算。

附图说明

图1是本发明TPS模型的结构示意图。

图2是本发明的流程图。

具体实施方式

下面结合实施例及附图对本发明作进一步详细的描述，但本发明的实施方式不限于此。

实施例

TPS是一种插值方法，它寻找一个通过所有的控制点的弯曲最小的光滑曲面。它的物理意义是：假设在原形状中有N个点An，这N个点在形变之后新坐标之下对应新的N个点Bn。用一个薄钢板的形变来模拟2D形变，确保这N个点能够正确匹配，那么怎样的形变，可以使钢板的弯曲能量最小？TPS插值是这个问题的数值解法。在进行了TPS运算后再进行等深线的追踪处理。

薄板样条函数（Thin Plate Spline，简称TPS）是由Duchon提出的一种关于几何设计的插值方法，通过模拟金属薄板的弯曲得到一个经过所有控制点且弯曲最小的光滑曲面。

等深线生成算法一般是基于三角网或规则格网进行追踪得到的，其精度依赖于模型的精度和内插算法。三角网对离散点的分布比较敏感，分布不均匀会导致稀疏处出现比较严重的失真，模型构建的精度和光滑程度不高，这些都直接影响等高线的效果。

利用TPS方法拟合得到的曲面，严格通过拟合点，曲面有严密的公式表达，可以精确计算曲面上每一点的高程值，而且计算过程中充分考虑所有拟合点的空间关系。计算出曲面参数后，构建一个覆盖曲面的格网模型，如图1所示，格网中的每个节点都可利用曲面公式计算，从而保证格网在光滑的曲面上。由于等深线是在该模型上进行追踪的，因此这些等深线肯定也在光滑曲面上，具有很好的光滑性。

如图2所示，本法发明基于TPS的等深线追踪算法，其步骤如下：

1、导入具有三维坐标的离散点作为控制点。

2、利用TPS对离散点进行曲面拟合，得到一个通过每个离散点的光滑曲面公式。若离散点数量较多，可以采取分块拟合的思路来提高计算效率。

3、根据离散点的密度和范围构建初始平面格网。

4、计算每个格网节点所在拟合曲面上的高程值，将平面格网转换成曲面格网。

5、根据给定的高程值在曲面格网上追踪等深线。

步骤2中，TPS曲面方程如下：

z＝a₁+a₂x+a₃y+Δ

$\mathrm{\Δ}=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{w}_i\·s^2\·\log \left(s\right)$

$s=\sqrt{{\left(x_j{-x}_i\right)}^2+{(y_j-y_i)}^2}$

式中z为点(x,y)处的高程值，a₁、a₂、a₃和w_i为方程系数，Δ为修正值，s为当前控制点与其余控制点的平面距离，n为控制点总点数；

将曲面方程用矩阵的形式表达如下：

由

$\left[\begin{array}{cc}K& P\\ P^T& O\end{array}\right]\·\left[\begin{array}{c}\stackrel{\→}{w}\\ \stackrel{\→}{a}\end{array}\right]=\left[\begin{array}{c}\stackrel{\→}{z}\\ \stackrel{\→}{o}\end{array}\right]$

得

$\left[\begin{array}{c}\stackrel{\→}{w}\\ \stackrel{\→}{a}\end{array}\right]={\left[\begin{array}{cc}K& P\\ P^T& O\end{array}\right]}^T\·\left[\begin{array}{c}\stackrel{\→}{z}\\ \stackrel{\→}{o}\end{array}\right]$

由上式即可求得方程系数，式中和为待求向量，形式如下：

$\stackrel{\→}{w}={\left[\begin{array}{cccc}{w}_1& w_2& \·\·\·& w_n\end{array}\right]}^T,\stackrel{\→}{a}={\left[\begin{array}{ccc}{a}_1& a_2& a_3\end{array}\right]}^T.$

K为n×n矩阵，每个元素由下式计算：

$K_{\mathrm{ij}}=s_{\mathrm{ij}}^2\·\log \left(s_{\mathrm{ij}}\right)+\frac{1}{n^2}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{s}_{\mathrm{ij}}$

其余矩阵和向量如下：

$P=\left[\begin{array}{ccc}1& x_1& y_1\\ 1& x_2& y_2\\ \·& \·& \·\\ \·& \·& \·\\ \·& \·& \·\\ 1& x_n& y_n\end{array}\right],$ $O=\left[\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ 0& 0& 0\end{array}\right],$ $\stackrel{\→}{z}=\left[\begin{array}{c}{z}_1\\ z_2\\ \·\\ \·\\ \·\\ z_n\end{array}\right],$ $\stackrel{\→}{o}=\left[\begin{array}{c}0\\ 0\\ 0\end{array}\right]$

解算过程如下：先根据每个控制点的平面坐标计算矩阵K，再组成矩阵P和向量，最后通过矩阵求逆和相关运算即可求得TPS曲面方程的系数。

另外，当离散点数量较多时，可以采取分块拟合的方法来提高计算效率，从TPS曲面方程的解算过程中可知，当控制点越多运算量越大，因此可根据控制点的分布，对其进行分块处理，每块建立一个子格网，步骤如下：

首先，将所有控制点的分布区域分成若干个正方块，确保每块里包含一定数量的控制点，若数量不足则合并到邻接的块。

然后，解算每一块的TPS曲面方程，根据方形块的范围建立子格网，利用方程系数计算每个格网节点的高程值。

最后，合并相邻的格网，最终形成一个覆盖整个区域的格网，相邻格网衔接处的节点有两个高程值，可取平均值。

在本实施例中，初始格网是所有格网节点没有高程值或高程值被设置为一个固定值，由于格网中的格子一般取为方形，要追踪等深线的区域往往不是矩形，因此可以利用控制点的分布情况来限制格网范围，如生成控制点的外包多边形或构建三角网，在多边形里面或三角网里面的格网节点设为待求点，在外面的节点设为无效点。通过这样种方法来约束等深线的追踪范围，且在利用TPS曲面方程计算格网节点的过程中，减少了不必要节点的计算。

通过求解TPS曲面方程的系数后，即可根据公式：

z＝a₁+a₂x+a₃y+Δ

计算任意点(x,y)的高程值z，为保证所求高程值的精度，要求点(x,y)要求最好在该曲面的控制点范围内（控制点构成的外围多边形内）。

当求解出曲面方程后，即可根据平面坐标(x,y)计算对应的高程值z，由于曲面方程覆盖整个格网所在的范围，每个格网节点的平面坐标均可由格网的起点坐标和节点所在的行列号计算得到，因此只要循环格网中每个节点，逐一计算对应的高程，格网就沿着曲面分布，也就成为曲面格网。

在曲面格网上追踪等深线的过程如下：设等深线的高程值为H，相邻两个格网节点的高程分别为H_i和H_j，若H的大小在H_i和H_j之间，则表明这两个节点的连线与该等深线有相交点，通过高差的比例关系求出相交点的坐标并保存起来，通过遍历整个格网得到该等深线与格网的所有交点，由这些交点连成的线即为等深线。

本实施例的具体应用如下，但并不限于下述两种场景的应用：

1、对于复杂地形，根据离散的高程点生成等高线，由于高程点的采集带有随机性，密度和分布不均匀，导致所建的模型失真，内插精度不高。因为TPS是过点拟合，曲面上每一处的高程都能根据所有拟合点的空间关系进行计算得到，充分利用测量点的精度，而且拟合的曲面较传统方法更贴近真实地形，从而能生成精度更高的等高线。

2、在航道测量中，根据一定的密度采集水深点的深度值，据此生成等深线。由于水下地形不可见，很难采集到特征线或特征点，而且大部分地形都比较平缓，从而大大增加追踪等深线的难度。基于TPS的曲面拟合算法能很好地根据采集点模拟水下地形，通过建立一定密度的规则格网可以十分精细地追踪等深线。又因曲面是光滑的，因此生成的等深线也具有一定的平滑性，最后利用张力样条或B次样条拟合就可得到光滑的等深线。

上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。

Claims (8)

1.一种基于TPS的等深线追踪算法，其特征在于，包括下述步骤：

S1、导入具有三维坐标的离散点作为控制点；

S2、利用TPS对离散点进行曲面拟合，得到一个通过每个离散点的光滑曲面公式；

TPS曲面方程如下：

z＝a₁+a₂x+a₃y+Δ；

式中z为点(x,y)处的高程值，a₁、a₂、a₃和w_i为方程系数，Δ为修正值，s为当前控制点与其余控制点的平面距离，n为控制点总点数；

将曲面方程用矩阵的形式表达如下：

由

得

由上式即可求得方程系数，式中和为待求向量，形式如下：

K为n×n矩阵，每个元素由下式计算：

其余矩阵和向量如下：

解算过程如下：先根据每个控制点的平面坐标计算矩阵K，再组成矩阵P和向量最后通过矩阵求逆和相关运算即可求得TPS曲面方程的系数；

S3、根据离散点的密度和范围构建初始平面格网；

S4、计算每个格网节点所在拟合曲面上的高程值，将平面格网转换成曲面格网；

S5、根据给定的高程值在曲面格网上追踪等深线。

2.根据权利要求1所述的一种基于TPS的等深线追踪算法，其特征在于，步骤S2中，从TPS曲面方程的解算过程中可知，当控制点越多运算量越大，因此可根据控制点的分布，对控制点进行分块处理，每块建立一个子格网，具体步骤如下：

S21、将所有控制点的分布区域分成若干个正方块，确保每块里包含一定数量的控制点，若数量不足则合并到邻接的块；

S22、解算每一块的TPS曲面方程，根据方形块的范围建立子格网，利用方程系数计算每个格网节点的高程值；

S23、合并相邻的格网，最终形成一个覆盖整个区域的格网。

3.根据权利要求2所述的一种基于TPS的等深线追踪算法，其特征在于，步骤S23中，相邻格网衔接处的节点有两个高程值，取其平均值作为该节点的高程值。

4.根据权利要求1所述的一种基于TPS的等深线追踪算法，其特征在于，步骤S3中，初始平面格网是所有格网节点没有高程值或高程值被设置为一个固定值，由于格网中的格子一般取为方形，要追踪等深线的区域往往不是矩形，因此利用控制点的分布情况来限制格网范围。

5.根据权利要求4所述的一种基于TPS的等深线追踪算法，其特征在于，利用控制点的分布情况来限制格网范围的具体方法为：生成控制点的外包多边形或构建三角网，在多边形里面或三角网里面的格网节点设为待求点，在外面的节点设为无效点。

6.根据权利要求1所述的一种基于TPS的等深线追踪算法，其特征在于，步骤S4中，通过求解TPS曲面方程的系数后，即可根据公式

z＝a₁+a₂x+a₃y+Δ；

计算任意点(x,y)的高程值z，为保证所求高程值的精度，要求点(x,y)在该曲面的控制点范围内，即控制点构成的外围多边形内。

7.根据权利要求4所述的一种基于TPS的等深线追踪算法，其特征在于，步骤S4，面格网也即初始格网，当求解出曲面方程后，即可根据平面坐标(x,y) 计算对应的高程值z，由于曲面方程覆盖整个格网所在的范围，每个格网节点的平面坐标均可由格网的起点坐标和节点所在的行列号计算得到，因此只要循环格网中每个节点，逐一计算对应的高程，格网就沿着曲面分布，也就成为曲面格网。

8.根据权利要求1所述的一种基于TPS的等深线追踪算法，其特征在于，步骤S5中，在曲面格网上追踪等深线的过程如下：设等深线的高程值为H，相邻两个格网节点的高程分别为H_i和H_j，若H的大小在H_i和H_j之间，则表明这两个节点的连线与该等深线有相交点，通过高差的比例关系求出相交点的坐标并保存起来，通过遍历整个格网得到该等深线与格网的所有交点，由这些交点连成的线即为等深线。