CN104154911B - 一种具有旋转不变性的海底地形二维匹配辅助导航方法 - Google Patents

一种具有旋转不变性的海底地形二维匹配辅助导航方法 Download PDF

Info

Publication number
CN104154911B
CN104154911B CN201410394482.3A CN201410394482A CN104154911B CN 104154911 B CN104154911 B CN 104154911B CN 201410394482 A CN201410394482 A CN 201410394482A CN 104154911 B CN104154911 B CN 104154911B
Authority
CN
China
Prior art keywords
mrow
mfrac
msub
gray
munderover
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201410394482.3A
Other languages
English (en)
Other versions
CN104154911A (zh
Inventor
徐晓苏
岳增阳
李佩娟
张涛
汤郡郡
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Southeast University
Original Assignee
Southeast University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Southeast University filed Critical Southeast University
Priority to CN201410394482.3A priority Critical patent/CN104154911B/zh
Publication of CN104154911A publication Critical patent/CN104154911A/zh
Application granted granted Critical
Publication of CN104154911B publication Critical patent/CN104154911B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C21/00—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
    • G01C21/005—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 with correlation of navigation data from several sources, e.g. map or contour matching
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C21/00—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
    • G01C21/10—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration
    • G01C21/12—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning
    • G01C21/16—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning by integrating acceleration or speed, i.e. inertial navigation
    • G01C21/165—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning by integrating acceleration or speed, i.e. inertial navigation combined with non-inertial navigation instruments
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C21/00—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
    • G01C21/20—Instruments for performing navigational calculations

Landscapes

  • Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Automation & Control Theory (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Image Analysis (AREA)

Abstract

本发明公开一种具有旋转不变性的海底地形二维匹配辅助导航方法，包括以下步骤：1)从多波束测深系统获得实时的条带扫测深度值矩阵阵列，并将高程阵列值转化为灰度值得到模板灰度图；读取主导航系统误差范围内区域的地形二维高程阵列值转化为灰度值，以形成背景灰度图；2)采用旋转不变的图像矩函数作为特征量，将背景灰度图中搜索与模板灰度图特征量匹配的区域、获得多波束测深系统实时扫测时所处的位置。本发明与以往通过地形一维匹配进行的匹配算法(如TERCOM匹配算法和SITAN匹配算法)相比，克服了一维匹配具有的匹配精度低及受地形限制的缺点，增大了地形特征描述量的丰富度，具有较高精度和适应性。

Description

一种具有旋转不变性的海底地形二维匹配辅助导航方法
技术领域
本发明涉及一种具有旋转不变特性的海底地形二维匹配辅助导航方法，属于水下辅助导航技术领域。
背景技术
地形辅助导航技术的研究始于40多年前，并于上世纪80年代中期通过一架CV-580设备试验机首次进行了应用测试。目前，国内外已经研究并开发了多种地形辅助导航系统，其中几个较为成熟的地形辅助导航系统已在实战中装备武器并得到成功的检验。
目前，国外地形辅助导航技术的成熟应用，还只限于飞机以及低空航行器导航，随着水下潜器应用领域的拓宽，对相关技术的需求不断增加，促进了科学家对水下地形辅助导航技术的研究。国内在水下潜器领域还没有成功应用地形辅助导航技术的例子，国外的水下潜器是否装有地形辅助导航系统也不得而知，但从现有的资料可知，国外在水下地形辅助导航技术的研究方面已取得了阶段性的成果，由于涉及保密性，所以具体使用技术未知。
综上可以看出，国外到目前阶段已经进入水下导航设备及系统的生产和海试阶段，而我国仍停留在理论和算法研究阶段，其中，地形匹配辅助导航算法是研究重点，目前国内最新的相关研究算法有基于单波束测深的ICCP匹配算法，基于海床特征地貌的水下导航方法研究，基于等深线图像的地形匹配算法研究以及基于FFT的海底地形二维匹配技术研究等。但是，上述算法要么是基于单波束序列的一维匹配，在信息丰富度和算法适应性方面较弱；要么是没有考虑潜器的各种方向的运动适配问题，即抗旋转的匹配问题。
发明内容
发明目的：本发明的目的在于针对现有技术存在的缺陷提出一种具有旋转不变性的海底地形二维匹配辅助导航方法。
技术方案：本发明的具有旋转不变性的海底地形二维匹配辅助导航方法，包括以下步骤：
1)从多波束测深系统获得实时的条带扫测深度值矩阵阵列，并将高程阵列值转化为灰度值得到模板灰度图；读取主导航系统误差范围内区域的地形二维高程阵列值转化为灰度值，以形成背景灰度图；
2)采用旋转不变的图像矩函数作为特征量，将背景灰度图中搜索与模板灰度图特征量匹配的区域、获得多波束测深系统实时扫测时所处的位置。
进一步地，还包括步骤3)：运行步骤1)和步骤2)若干次，利用多次匹配的结果得到精确航迹并返回至主导航系统修正误差。
进一步地，所述的旋转不变图像矩函数为Zernike矩、Hu矩或小波矩。
进一步地，所述在背景灰度图中搜索与模板灰度图特征量匹配的区域的方法具体为，在背景灰度图中取出一系列子图，并在计算子图的图像矩函数特征量与计算得到的模板灰度图的图像矩函数特征量相比较，通过找到与模板灰度图最匹配的子图，然后以该子图在背景灰度图的位置作为模板灰度图的位置，以确定多波束测深系统实时扫描时所处的位置，其中在计算子图和模板灰度图图像矩函数特征量匹配的过程中采用圆窗口化处理后再计算特征量。
进一步地，所述步骤1)具体为：
从多波束测深系统获得实时的条带扫测深度值矩形阵列，根据主惯导误差范围从原始先验数据库得到背景深度值矩阵阵列；
把各个深度值转化成0-255之间的灰度值，通过下述公式归一化灰度转换，使深度值矩阵阵列形成待匹配的背景灰度图和模板灰度图；
其中，hi,j表示二维深度值阵列中第(i，j)个点的深度值，round(·)表示临近取整，|hi,j|表示第(i，j)个点深度值的绝对值，h′i,j表示该点转化后的灰度值，分别表示取二维深度值阵列中元素绝对值的最小值和最大值；
所述步骤2)具体为:
针对实时条带模板灰度图，以扫测条带宽度为边长取出初始位置正方形图作为模板图；针对背景灰度图，从左上角开始依次逐像素地取和模板正方形相同大小的图作为子图，直至遍历整个背景图；
对每一个正方形灰度子图和模板图进行如下圆窗口化操作:只取此正方形内切圆的像素点，其余像素点置零；
计算正方形灰度子图和模板图的Zernike多阶矩，具体如下：
由于Zernike矩是定义为极坐标系单位圆(r,θ)上的函数，Zernike多项式将不得不对每个像素位置求值，计算量显然很大。为了计算简便和快速，进行方圆变换；
经过方圆变换，有如下关系：
γ＝max{|x|,|y|}
其中，γ表示某点变换后在极坐标系里与原点的半径，x，y分别表示变换前某点在矩形的直角坐标系里的坐标，max{}表示取括号内元素的最大值；
如果|x|＝γ，则
如果|y|＝γ，则
其中，(γ,ξ)表示变换后对应的该像素点在圆形的极坐标系里的坐标；
若N代表图像的像素数，则坐标值的范围为：
1≤ξ≤8γ
像素(γ,ξ)的归一化极坐标(r,θ)为:
其中，(r,θ)表示归一化到单位圆内的极坐标，N表示图像像素总的点数；
把矩形图像二维分布f(x,y)变换为单位圆极坐标系里的图像分布为g(r,θ)；
根据Zernike矩的定义式：
其中，Zn,m表示图像的n阶m重Zernike矩，Bnmk表示图像的n阶m重径向矩，g(r,θ)即为极坐标系(r,θ)上的图像分布，r表示半径，θ表示角度；
可得，变换后圆形图像的Zernike矩的实部和虚部如下：
其中，分别表示p阶q重Zernike矩的实部和虚部，Rpq(·)为Zernike矩的p阶q重径向多项式；
根据Zernike矩的特性，旋转不变量有Zernike矩的各阶模值；
取其中的前四个模值|Z20|,|Z22|,|Z31|,|Z33|作为相似性度量的特征量，记为：
Z(i),(1≤i≤4)；
其中，Z(1)＝|Z20|,Z(2)＝|Z22|,Z(3)＝|Z31|,Z(4)＝|Z33|；
针对背景灰度图，采用遍历搜索法截取每个与实时模板图相同大小的正方形子图，对于每个子图计算四个模值|Z20|,|Z22|,|Z31|,|Z33|，得到每个子图的特征量；
对于模板和每个子图，使用归一化互相关算法作为相似性度量方法，得到最佳匹配子图，即模板的位置；公式如下：
其中，δ(x,y)表示模板图与坐标为(x,y)的子图之间的相关性系数，值越大，相似性越高；ZR(i),ZT(i)分别表示子图和模板图的各阶Zernike矩，分别表示子图和模板图的各阶Zernike矩的平均值；
找到δ(x,y)最大的位置子图，即为最佳匹配位置子图，以该子图范围作为获得多波束测深系统实时扫测时经过的位置。
进一步地，还包括步骤3)：运行步骤1)和步骤2)若干次，得到若干个最佳位置匹配子图中心连接形成的位置轨迹，然后将该位置轨迹返回主导航系统实时修正主导航系统的误差。
进一步地，所述的主导航系统为捷联惯性导航系统。
本发明与现有技术相比，其有益效果是：1、与以往通过地形一维匹配进行的匹配算法(如TERCOM匹配算法和SITAN匹配算法)相比，克服了一维匹配具有的匹配精度低及受地形限制的缺点，增大了地形特征描述量的丰富度，具有较高精度和适应性。2、创新性的把水深网格数据转化为图像像素数据，借鉴成熟的图像匹配领域的相关算法进行水下数字地图匹配，拓宽了地形匹配方法的范围，也使水下地形辅助导航更加高效，精确和快速。3、以图像旋转不变矩(本专利主要采用Zernike多阶矩)作为两个区域的特征量进行相似性度量，克服了传统二维阵列匹配中存在的旋转问题，在潜器未知运行方向下也能实现水下地形辅助匹配导航。4、在抗旋转匹配的搜索过程中，把正方形的模板和子图采用只取内切圆的窗口化处理，保证了在背景图上采用遍历搜索法也能搜索到某个正方形子图可以和任意旋转的模板图相匹配。5、采用归一化互相关算法作为相似性度量方法，能快速准确地实现模板图和子图的匹配。6、采用图像矩进行单次匹配得到的是潜器的位置范围，但是多次进行连续的动态匹配就可以指示出潜器的准确航迹，实时的辅助惯导修正误差。
附图说明
图1为仿真地形三维图；
图2为地形深度值阵列转化后的灰度图；
图3为潜器的两条运动方向扫测轨迹示意图；
图4为模板图和圆窗口化处理后的模板图；
图5为计算Zernike矩时的方-圆变换示意图；
图6为多次动态匹配形成的潜器运动轨迹示意图。
具体实施方式
下面对本发明技术方案进行详细说明，但是本发明的保护范围不局限于所述实施例。
本实施例的具有旋转不变特性的海底地形二维匹配辅助导航方法，采用以下方式实施：
步骤1：
从多波束测深系统获得实时的条带扫测深度值矩形阵列，根据主惯导误差范围从原始先验数据库得到背景深度值矩阵阵列。仿真数据采用国家海洋科学数据库中的东经127.528°～128.205°，北纬27.328°～28.005°海底地形高程数据为基础，插值形成100×100m的规则格网数据矩阵，矩阵大小为751×751，地形三维图如图1。
把各个深度值转化成0-255之间的灰度值，通过下述公式归一化灰度转换，使深度值矩阵阵列形成待匹配的背景灰度图和模板灰度图，如图2和图4所示。
其中，hi,j表示二维深度值阵列中第(i，j)个点的深度值，round(·)表示临近取整，|hi,j|表示第(i，j)个点深度值的绝对值，h′i,j表示该点转化后的灰度值，分别表示取二维深度值阵列中元素绝对值的最小值和最大值。
步骤2)：
在一次匹配过程中，对于模板灰度图，需要在背景灰度图中搜索与之相匹配的子图，具体过程如下：
考虑到计算图像不变矩(本方案以Zernike矩为例)时旋转不变问题，为了保证搜索到的子图和模板的初始匹配信息量在旋转条件下不改变，需要对正方形灰度子图和模板图进行如下圆窗口化操作:只取此正方形内切圆的像素点，其余像素点置零(即全黑)，如图4所示。
计算此正方形灰度图的图像不变矩(本方案以Zernike多阶矩为例)，具体如下：
由于Zernike矩是定义为极坐标系单位圆(r,θ)上的函数，Zernike多项式将不得不对每个像素位置求值，计算量显然很大。为了计算简便和快速，进行如下方-圆变换，如图5。
经过方-圆变换，有如下关系：
γ＝max{|x|,|y|}
其中，γ表示某点变换后在极坐标系里与原点的半径，x，y分别表示变换前某点在矩形的直角坐标系里的坐标，max{·}表示取括号内元素的最大值。
如果|x|＝γ，则
如果|y|＝γ，则
其中，(γ,ξ)表示变换后对应的该像素点在圆形的极坐标系里的坐标。
若N代表图像的像素数，则坐标值的范围为：
1≤ξ≤8γ
像素(γ,ξ)的归一化极坐标(r,θ)为:
其中，(r,θ)表示归一化到单位圆内的极坐标，N表示图像像素总的点数。
把矩形图像二维分布f(x,y)变换为单位圆极坐标系里的图像分布为g(r,θ)。
根据Zernike矩的定义式：
其中，Zn,m表示图像的n阶m重Zernike矩，Bnmk表示图像的n阶m重径向矩，g(r,θ)即为极坐标系(r,θ)上的图像分布，r表示半径，θ表示角度。
可得，变换后圆形图像的Zernike矩的实部和虚部如下：
其中，分别表示p阶q重Zernike矩的实部和虚部，Rpq(·)为Zernike矩的p阶q重径向多项式，表达式见下表。
径向多项式Rpq(r)的表达式
根据Zernike矩的特性，旋转不变量有Zernike矩的各阶模值。取其中的前四个模值|Z20|,|Z22|,|Z31|,|Z33|作为相似性度的特征量，记为：
Z(i),(1≤i≤4)
其中，Z(1)＝|Z20|,Z(2)＝|Z22|,Z(3)＝|Z31|,Z(4)＝|Z33|。
针对背景灰度图，采用遍历搜索法截取每个与实时模板图相同大小的正方形子图，对于每个子图计算四个模值|Z20|,|Z22|,|Z31|,|Z33|，得到每个子图的特征量。
对于模板和每个子图，使用归一化互相关算法作为相似性度量方法，得到最佳匹配子图，即模板的位置。公式如下：
其中，δ(x,y)表示模板图与坐标为(x,y)的子图之间的相关性系数，值越大，相似性越高。ZR(i),ZT(i)分别表示子图和模板图的各阶Zernike矩，分别表示子图和模板图的各阶Zernike矩的平均值。
找到δ(x,y)最大的位置图，即为最佳匹配位置图，此图为潜器的位置范围，由于潜器运行方向未知，所以进行步骤3。
步骤3：进行多次动态匹配(即把模板图逐像素移到下一位置，对于每一个模板图，重复一次步骤1)和步骤2))，找到一个多位置图轨迹，在背景图中就是潜航器的精确运行轨迹，此轨迹和背景图的正北方向夹角α即已知。由于一次匹配时间很短，故可以在一次匹配过程中把潜器看成匀速直线运动，根据一次匹配完成时间t，潜器运动速度V，最后一次匹配的位置图和轨迹的交叉点位置(X,Y)，已知的前进轨迹方向和正北方向夹角α，可以得到潜器的精确位置信息(X′,Y′)，公式如下，过程如示意图6。进而完成地形匹配辅助导航。
X′＝X+Vt sinα
Y′＝Y+Vt cosα
其中，(X′,Y′)表示潜器实际位置，(X,Y)表示最后一次匹配开始时的位置，α表示前进轨迹方向和正北方向夹角，t表示一次匹配完成时间，V表示潜器运动速度。并通过上述结果范围本实施例的主导航系统惯性捷联制导系统修正惯性捷联制导系统的累积误差。
如上所述，尽管参照特定的优选实施例已经表示和表述了本发明，但其不得解释为对本发明自身的限制。在不脱离所附权利要求定义的本发明的精神和范围前提下，可对其在形式上和细节上作出各种变化。

Claims (3)

1.一种具有旋转不变性的海底地形二维匹配辅助导航方法，其特征在于，包括以下步骤：
1)从多波束测深系统获得实时的条带扫测深度值矩阵阵列，并将高程阵列值转化为灰度值得到模板灰度图；读取主导航系统误差范围内区域的地形二维高程阵列值转化为灰度值，以形成背景灰度图；
所述步骤1)具体为：
从多波束测深系统获得实时的条带扫测深度值矩形阵列，根据主惯导误差范围从原始先验数据库得到背景深度值矩阵阵列；
把各个深度值转化成0-255之间的灰度值，通过下述公式归一化灰度转换，使深度值矩阵阵列形成待匹配的背景灰度图和模板灰度图；
<mrow> <msubsup> <mi>h</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> <mo>&amp;prime;</mo> </msubsup> <mo>=</mo> <mi>r</mi> <mi>o</mi> <mi>u</mi> <mi>n</mi> <mi>d</mi> <mrow> <mo>(</mo> <mfrac> <mrow> <mo>|</mo> <msub> <mi>h</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> </msub> <mo>|</mo> <mo>-</mo> <munder> <mi>min</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> </munder> <mo>{</mo> <mo>|</mo> <msub> <mi>h</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> </msub> <mo>|</mo> <mo>}</mo> </mrow> <mrow> <munder> <mi>max</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> </munder> <mo>{</mo> <mo>|</mo> <msub> <mi>h</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> </msub> <mo>|</mo> <mo>}</mo> <mo>-</mo> <munder> <mi>min</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> </munder> <mo>{</mo> <mo>|</mo> <msub> <mi>h</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> </msub> <mo>|</mo> <mo>}</mo> </mrow> </mfrac> <mo>&amp;times;</mo> <mn>255</mn> <mo>)</mo> </mrow> </mrow>
其中，hi,j表示二维深度值阵列中第(i，j)个点的深度值，round(·)表示临近取整，|hi,j|表示第(i，j)个点深度值的绝对值，h′i,j表示该点转化后的灰度值，
分别表示取二维深度值阵列中元素绝对值的最小值和最大值；
2)采用旋转不变的图像矩函数作为特征量，在背景灰度图中取出一系列子图，并在计算子图的图像矩函数特征量与计算得到的模板灰度图的图像矩函数特征量相比较，通过找到与模板灰度图最匹配的子图，然后以该子图在背景灰度图的位置作为模板灰度图的位置，以确定多波束测深系统实时扫描时所处的位置，其中在计算子图和模板灰度图图像矩函数特征量匹配的过程中采用圆窗口化处理后再计算特征量；所述的旋转不变图像矩函数为Zernike矩；
所述步骤2)具体为:
针对实时条带模板灰度图，以扫测条带宽度为边长取出初始位置正方形图作为模板图；针对背景灰度图，从左上角开始依次逐像素地取和模板正方形相同大小的图作为子图，直至遍历整个背景图；
对每一个正方形灰度子图和模板图进行如下圆窗口化操作:只取此正方形内切圆的像素点，其余像素点置零；
计算正方形灰度子图和模板图的Zernike多阶矩，具体如下：
为了计算简便和快速，进行方圆变换；
经过方圆变换，有如下关系：
γ＝max{|x|,|y|}
其中，γ表示某点变换后在极坐标系里与原点的半径，x，y分别表示变换前某点在矩形的直角坐标系里的坐标，max{}表示取括号内元素的最大值；
如果|x|＝γ，则
如果|y|＝γ，则
其中，(γ,ξ)表示变换后对应的该像素点在圆形的极坐标系里的坐标；
若N代表图像的像素数，则坐标值的范围为：
<mrow> <mo>-</mo> <mfrac> <mi>N</mi> <mn>2</mn> </mfrac> <mo>&amp;le;</mo> <mi>x</mi> <mo>&amp;le;</mo> <mfrac> <mi>N</mi> <mn>2</mn> </mfrac> <mo>;</mo> <mo>-</mo> <mfrac> <mi>N</mi> <mn>2</mn> </mfrac> <mo>&amp;le;</mo> <mi>y</mi> <mo>&amp;le;</mo> <mfrac> <mi>N</mi> <mn>2</mn> </mfrac> <mo>;</mo> <mn>0</mn> <mo>&amp;le;</mo> <mi>&amp;gamma;</mi> <mo>&amp;le;</mo> <mfrac> <mi>N</mi> <mn>2</mn> </mfrac> <mo>;</mo> <mn>1</mn> <mo>&amp;le;</mo> <mi>&amp;xi;</mi> <mo>&amp;le;</mo> <mn>8</mn> <mi>&amp;gamma;</mi> </mrow>
像素(γ,ξ)的归一化极坐标(r,θ)为:
<mrow> <mi>r</mi> <mo>=</mo> <mfrac> <mrow> <mn>2</mn> <mi>&amp;gamma;</mi> </mrow> <mi>N</mi> </mfrac> <mo>,</mo> <mi>&amp;theta;</mi> <mo>=</mo> <mfrac> <mrow> <mi>&amp;pi;</mi> <mi>&amp;xi;</mi> </mrow> <mrow> <mn>4</mn> <mi>&amp;gamma;</mi> </mrow> </mfrac> </mrow>
其中，(r,θ)表示归一化到单位圆内的极坐标，N表示图像像素总的点数；
把矩形图像二维分布f(x,y)变换为单位圆极坐标系里的图像分布为g(r,θ)；
根据Zernike矩的定义式：
<mrow> <msub> <mi>Z</mi> <mrow> <mi>n</mi> <mi>m</mi> </mrow> </msub> <mo>=</mo> <mfrac> <mrow> <mi>n</mi> <mo>+</mo> <mn>1</mn> </mrow> <mi>&amp;pi;</mi> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>k</mi> <mo>=</mo> <mi>m</mi> </mrow> <mi>n</mi> </munderover> <msub> <mi>B</mi> <mrow> <mi>n</mi> <mi>m</mi> <mi>k</mi> </mrow> </msub> <munderover> <mo>&amp;Integral;</mo> <mn>0</mn> <mrow> <mn>2</mn> <mi>&amp;pi;</mi> </mrow> </munderover> <munderover> <mo>&amp;Integral;</mo> <mn>0</mn> <mn>1</mn> </munderover> <msup> <mi>r</mi> <mi>k</mi> </msup> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mi>i</mi> <mi>m</mi> <mi>&amp;theta;</mi> </mrow> </msup> <mi>g</mi> <mrow> <mo>(</mo> <mi>r</mi> <mo>,</mo> <mi>&amp;theta;</mi> <mo>)</mo> </mrow> <mi>r</mi> <mi>d</mi> <mi>r</mi> <mi>d</mi> <mi>&amp;theta;</mi> <mo>,</mo> <mrow> <mo>(</mo> <mi>r</mi> <mo>&lt;</mo> <mn>1</mn> <mo>)</mo> </mrow> </mrow>
其中，Zn,m表示图像的n阶m重Zernike矩，Bnmk表示图像的n阶m重径向矩，g(r,θ)即为极坐标系(r,θ)上的图像分布，r表示半径，θ表示角度；
可得，变换后圆形图像的Zernike矩的实部和虚部如下：
<mrow> <msubsup> <mi>Z</mi> <mrow> <mi>p</mi> <mi>q</mi> </mrow> <mrow> <mo>(</mo> <mi>c</mi> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <mfrac> <mrow> <mn>2</mn> <mi>p</mi> <mo>+</mo> <mn>2</mn> </mrow> <msup> <mi>N</mi> <mn>2</mn> </msup> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>&amp;gamma;</mi> <mo>=</mo> <mn>1</mn> </mrow> <mrow> <mi>N</mi> <mo>/</mo> <mn>2</mn> </mrow> </munderover> <msub> <mi>R</mi> <mrow> <mi>p</mi> <mi>q</mi> </mrow> </msub> <mrow> <mo>(</mo> <mfrac> <mrow> <mn>2</mn> <mi>&amp;gamma;</mi> </mrow> <mi>N</mi> </mfrac> <mo>)</mo> </mrow> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>&amp;xi;</mi> <mo>=</mo> <mn>1</mn> </mrow> <mrow> <mn>8</mn> <mi>&amp;gamma;</mi> </mrow> </munderover> <mi>cos</mi> <mrow> <mo>(</mo> <mfrac> <mrow> <mi>&amp;pi;</mi> <mi>q</mi> <mi>&amp;xi;</mi> </mrow> <mrow> <mn>4</mn> <mi>&amp;gamma;</mi> </mrow> </mfrac> <mo>)</mo> </mrow> <mi>g</mi> <mrow> <mo>(</mo> <mi>r</mi> <mo>,</mo> <mi>&amp;theta;</mi> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>Z</mi> <mrow> <mi>p</mi> <mi>q</mi> </mrow> <mrow> <mo>(</mo> <mi>s</mi> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <mfrac> <mrow> <mo>-</mo> <mrow> <mo>(</mo> <mn>2</mn> <mi>p</mi> <mo>+</mo> <mn>2</mn> <mo>)</mo> </mrow> </mrow> <msup> <mi>N</mi> <mn>2</mn> </msup> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>&amp;gamma;</mi> <mo>=</mo> <mn>1</mn> </mrow> <mrow> <mi>N</mi> <mo>/</mo> <mn>2</mn> </mrow> </munderover> <msub> <mi>R</mi> <mrow> <mi>p</mi> <mi>q</mi> </mrow> </msub> <mrow> <mo>(</mo> <mfrac> <mrow> <mn>2</mn> <mi>&amp;gamma;</mi> </mrow> <mi>N</mi> </mfrac> <mo>)</mo> </mrow> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>&amp;xi;</mi> <mo>=</mo> <mn>1</mn> </mrow> <mrow> <mn>8</mn> <mi>&amp;gamma;</mi> </mrow> </munderover> <mi>sin</mi> <mrow> <mo>(</mo> <mfrac> <mrow> <mi>&amp;pi;</mi> <mi>q</mi> <mi>&amp;xi;</mi> </mrow> <mrow> <mn>4</mn> <mi>&amp;gamma;</mi> </mrow> </mfrac> <mo>)</mo> </mrow> <mi>g</mi> <mrow> <mo>(</mo> <mi>r</mi> <mo>,</mo> <mi>&amp;theta;</mi> <mo>)</mo> </mrow> </mrow>
其中，分别表示p阶q重Zernike矩的实部和虚部，Rpq(·)为Zernike矩的p阶q重径向多项式；
根据Zernike矩的特性，旋转不变量有Zernike矩的各阶模值；
取其中的前四个模值|Z20|,|Z22|,|Z31|,|Z33|作为相似性度量的特征量，记为：
Z(i),(1≤i≤4)；
其中，Z(1)＝|Z20|,Z(2)＝|Z22|,Z(3)＝|Z31|,Z(4)＝|Z33|；
针对背景灰度图，采用遍历搜索法截取每个与实时模板图相同大小的正方形子图，对于每个子图计算四个模值|Z20|,|Z22|,|Z31|,|Z33|，得到每个子图的特征量；
对于模板和每个子图，使用归一化互相关算法作为相似性度量方法，得到最佳匹配子图，即模板的位置；公式如下：
<mrow> <mi>&amp;delta;</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>)</mo> </mrow> <mo>=</mo> <mfrac> <mrow> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>4</mn> </munderover> <mrow> <mo>(</mo> <msub> <mi>Z</mi> <mi>R</mi> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>-</mo> <msub> <mover> <mi>Z</mi> <mo>&amp;OverBar;</mo> </mover> <mi>R</mi> </msub> <mo>)</mo> </mrow> <mrow> <mo>(</mo> <msub> <mi>Z</mi> <mi>T</mi> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>-</mo> <msub> <mover> <mi>Z</mi> <mo>&amp;OverBar;</mo> </mover> <mi>T</mi> </msub> <mo>)</mo> </mrow> </mrow> <msqrt> <mrow> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>4</mn> </munderover> <msup> <mrow> <mo>(</mo> <msub> <mi>Z</mi> <mi>R</mi> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>-</mo> <msub> <mover> <mi>Z</mi> <mo>&amp;OverBar;</mo> </mover> <mi>R</mi> </msub> <mo>)</mo> </mrow> <mn>2</mn> </msup> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>4</mn> </munderover> <msup> <mrow> <mo>(</mo> <msub> <mi>Z</mi> <mi>T</mi> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>-</mo> <msub> <mover> <mi>Z</mi> <mo>&amp;OverBar;</mo> </mover> <mi>T</mi> </msub> <mo>)</mo> </mrow> <mn>2</mn> </msup> </mrow> </msqrt> </mfrac> </mrow>
其中，δ(x,y)表示模板图与坐标为(x,y)的子图之间的相关性系数，值越大，相似性越高；ZR(i),ZT(i)分别表示子图和模板图的各阶Zernike矩，分别表示子图和模板图的各阶Zernike矩的平均值；
找到δ(x,y)最大的位置子图，即为最佳匹配位置子图，以该子图范围作为获得多波束测深系统实时扫测时经过的位置范围。
2.根据权利要求1所述的具有旋转不变性的海底地形二维匹配辅助导航方法，其特征在于，还包括步骤3)：运行步骤1)和步骤2)若干次，利用多次匹配的结果得到若干个最佳匹配位置子图中心连接形成的位置轨迹，然后将该位置轨迹返回主导航系统实时修正主导航系统的误差。
3.根据权利要求1或2所述的具有旋转不变性的海底地形二维匹配辅助导航方法，其特征在于，所述的主导航系统为捷联惯性导航系统。
CN201410394482.3A 2014-08-12 2014-08-12 一种具有旋转不变性的海底地形二维匹配辅助导航方法 Active CN104154911B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201410394482.3A CN104154911B (zh) 2014-08-12 2014-08-12 一种具有旋转不变性的海底地形二维匹配辅助导航方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201410394482.3A CN104154911B (zh) 2014-08-12 2014-08-12 一种具有旋转不变性的海底地形二维匹配辅助导航方法

Publications (2)

Publication Number Publication Date
CN104154911A CN104154911A (zh) 2014-11-19
CN104154911B true CN104154911B (zh) 2017-08-25

Family

ID=51880476

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201410394482.3A Active CN104154911B (zh) 2014-08-12 2014-08-12 一种具有旋转不变性的海底地形二维匹配辅助导航方法

Country Status (1)

Country Link
CN (1) CN104154911B (zh)

Families Citing this family (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105160665A (zh) * 2015-08-25 2015-12-16 东南大学 一种双圆子模板水下地形匹配方法
CN105447875A (zh) * 2015-12-09 2016-03-30 新疆畜牧科学院草业研究所 一种电子地形图自动几何校正方法
CN105787489B (zh) * 2016-03-04 2019-06-14 哈尔滨工程大学 一种基于水下地貌的匹配导航算法
CN107368528A (zh) * 2017-06-12 2017-11-21 北京航空航天大学 用于地形匹配的3D Zernike矩预计算技术
CN107643082B (zh) * 2017-09-05 2020-03-31 东南大学 基于多波束的多路径并行iccp水下地形匹配方法
CN108318034B (zh) * 2018-01-09 2020-05-22 浙江大学 一种基于声呐地图的auv回坞导航方法
CN109029434A (zh) * 2018-06-29 2018-12-18 电子科技大学 基于自适应尺度下推的桑迪亚惯性地形辅助导航方法
CN110378363B (zh) * 2019-05-16 2021-06-04 中国地质大学（武汉） 基于匹配控制子图的复合多边形相似性度量方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP1918677A2 (en) * 2006-10-31 2008-05-07 Honeywell International Inc. Systems and methods for a terrain contour matching navigation system
CN101339036A (zh) * 2008-08-20 2009-01-07 北京航空航天大学 地形辅助导航方法和设备
CN101887522A (zh) * 2010-07-13 2010-11-17 南京大学 基于多点同时快速行进的遥感影像道路提取方法
CN103047983A (zh) * 2013-01-14 2013-04-17 哈尔滨工程大学 水下机器人的面地形匹配导航方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP1918677A2 (en) * 2006-10-31 2008-05-07 Honeywell International Inc. Systems and methods for a terrain contour matching navigation system
CN101339036A (zh) * 2008-08-20 2009-01-07 北京航空航天大学 地形辅助导航方法和设备
CN101887522A (zh) * 2010-07-13 2010-11-17 南京大学 基于多点同时快速行进的遥感影像道路提取方法
CN103047983A (zh) * 2013-01-14 2013-04-17 哈尔滨工程大学 水下机器人的面地形匹配导航方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
基于3D Zernike矩的三维地形匹配算法及性能分析;叶斌 等;《宇航学报》;20070930;第28卷(第5期);第1241-1245页 *
基于Haar小波的地形匹配技术;王广君 等;《电波科学学报》;20020630;第17卷(第3期);第240-244页 *
水下地形辅助导航新方法仿真;徐遵义 等;《舰船科学技术》;20120630;第34卷(第6期);第53-56、106页 *

Also Published As

Publication number Publication date
CN104154911A (zh) 2014-11-19

Similar Documents

Publication Publication Date Title
CN104154911B (zh) 一种具有旋转不变性的海底地形二维匹配辅助导航方法
CN109269430B (zh) 基于深度提取模型的多株立木胸径被动测量方法
AU2007355942B2 (en) Arrangement and method for providing a three dimensional map representation of an area
CN104268935A (zh) 一种基于特征的机载激光点云与影像数据融合系统及方法
CN107907124A (zh) 基于场景重识的定位方法、电子设备、存储介质、系统
CN109471101A (zh) 一种弹载sar图像校正方法
Ozog et al. Building 3d mosaics from an autonomous underwater vehicle, doppler velocity log, and 2d imaging sonar
CN113177593B (zh) 一种水上交通环境中雷达点云与影像数据的融合方法
CN107218926A (zh) 一种基于无人机平台的远程扫描的数据处理方法
CN107766818A (zh) 一种启发式的水下结构化环境线特征提取方法
CN105160665A (zh) 一种双圆子模板水下地形匹配方法
CN103700135B (zh) 一种三维模型局部球面调和特征提取方法
Song et al. Application of acoustic image processing in underwater terrain aided navigation
CN112540371A (zh) 一种近底多波束坐标转换处理方法
Almanza-Medina et al. Deep learning architectures for navigation using forward looking sonar images
CN104166995A (zh) 一种基于马步测度的Harris-SIFT双目视觉定位方法
CN106408600B (zh) 一种用于太阳高分辨率图像中图像配准的方法
CN110927765B (zh) 激光雷达与卫星导航融合的目标在线定位方法
Zhuangsheng et al. Study on initial gravity map matching technique based on triangle constraint model
CN105205817B (zh) 一种基于声呐图像边缘角点直方图的水下地形匹配方法
CN109886988A (zh) 一种微波成像仪定位误差的度量方法、系统、装置及介质
CN113160191B (zh) 一种基于激光雷达的环境构图完整性判别方法和装置
CN115507752A (zh) 基于平行环境元素的单目视觉测距方法和系统
CN107907881A (zh) 一种长孔径星载sar地形估计方法
CN104700400A (zh) 一种基于特征点训练的高精度图像匹配方法

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant