Abstract
Global positioning system signal does not penetrate into the water volume, and autonomous underwater vehicle have to use the inertial navigation system that will cause an inevitable cumulative error. Terrain reference position can zero out the inertial navigation system error and have been widely used. To improve the positioning accuracy, the planning path of the underwater robot is required to pass through the suitable matching areas in the path planning stage. So, a gridding map is needed to quantify the suitability of the terrain and to partition the suitable matching area and the unsuitable matching area optimally. In this article, we will focus on the quantitative method of terrain suitability and the grid parameter solution method for optimal partition of suitable matching area. Finally, the validity of the algorithm is obtained using the ship borne data. The results show that in the low suitability blocks, the influence of measurement error on terrain reference position accuracy is higher than that in high suitability blocks. And the average of terrain reference position deviation obtained of the path passed through the high suitability area is 53.3% which is lower than that of the path passed through the low suitability area, and the standard deviation of the position deviation is reduced by 21.15%.
Keywords
Introduction
Terrain reference navigation (TRN) is an ideal aided navigation method that has no time accumulation error. At present, underwater TRN technology uses multi-beam as the measurement equipment, which improves both acquisition speed and matching accuracy. 1,2 However, there are some disadvantages associated with TRN, such as large errors in terrain measurement, fewer terrain features, and the lower acquisition speed of the measurement data, which always cause large position deviation. 1 At present, the research on improving the TRN/terrain reference position (TRP) accuracy mainly focuses on the filtering algorithm, including particle filter (PF), 1,3 –6 robust filtering algorithm, 7,8 –10 point mass filter (PMF), 5,11,12 and so on. The latest research progress on TRN can be referred to the review article. 13 –17 Another way is to improve the positioning accuracy from the perspective of path planning, allowing autonomous underwater vehicle (AUV) to go through the terrain with high adaptability. As shown in Figure 1(a), the error of the reference navigation system exceeds the threshold at the position C and the TRP is needed to correct the navigation error when the AUV navigation is from A to B. If AUV can independently quantify the adaptability of the digital elevation map (DEM) and identification of suitable areas and unsuitable area. Then AUV will be able to obtain a suitability map (Figure 1(b)) to plan a path through a high suitable area (blue route in Figure 1(b)), so that the positioning error at the end of the navigation point will be in a small range. Otherwise, AUV may not reach the target point B through the low suitable area (green route in Figure 1(b)). To make AUV use the topographic map path planning and correct navigation error, the first problem that needs to be solved is to quantify the suitability of the terrain and the problem of the division of the suitable area and the unsuitable area. Now, the TRN path planning method mainly focuses on the path search algorithm. For example, the waypoints selected are based on the binary-tree search algorithm, 18 and the terrain matching navigation method based on the A*. 19 –21 It is worth noting that the key to TRN path planning include two prerequisites: (1) suitability quantization and (2) division of suitable area and unsuitable area. Therefore, this study is mainly to solve these two problems. Because large number of path planning algorithms are based on grid map. 19,21,22 So, the suitable area and the unsuitable area partition are based on topographic gridding. In this article, a method of quantization of terrain suitability and optimal partition of terrain suitability is mainly studied. Based on the DEM gridding, the suitable area and the unsuitable area can be partitioned off optimally.

Terrain suitability quantification and suitable area partition can realize AUV navigation mode decision and TRN path planning autonomously. AUV: autonomous underwater vehicle; TRN: terrain reference navigation.
There are many quantization methods for terrain suitability 21,22 that is based on terrain height standard deviations. But terrain height standard deviation is not the only factor affecting the suitability. In the following, we will also deduce another important influence factor of terrain suitability, direction, and a more reliable method of suitability quantification is obtained and part the suitable and unsuitable area with gridding method. The size of the mesh of gridding map will affect the suitability of the subblocks. The DEM has been gridded optimally when the suitable area and the unsuitable area are parted into different subblocks to the maximum extent at a certain value of subblock size. In this article, we primarily study parameterized representation methods of terrain suitability and then propose the optimal partition methods for the terrain suitability area based on the grid premise.
The following research contents of the article include the following parts. In second section, the influence factors and the quantitative evaluation method of TRP accuracy are analyzed. At last, the quantization parameter signal-to-noise ratio (SNR) of the local terrain area is defined. In the third section, the suitability matrix SNR
Suitability and TRP Accuracy
Influencing factors on TRP accuracy
Suppose AUV has obtained measurement terrain map (MTM)
In the equation,
In the equation,
In the equation,
Measurement sequence
Consider the fact that TRP is a two-dimensional search process and that TRP accuracy is related to measurement error and local topography gradients; Increases in the terrain gradient can improve positioning accuracy; Because the terrain gradients have different values in different directions, TRP accuracy will be different in different directions.
The last conclusion is critical. We can infer that even if the terrain changes drastically, if the gradient in a certain direction is small, it may also cause the TRP results to deteriorate as analyzed in “Terrain suitability quantification” section.
Accuracy evaluation of TRP
The terrain surface is an anisotropic surface. According to estimation theory, we can obtain more information for the positioning point.
2
A reliability evaluation of TRP can be obtained by calculating the information content at TRP position point (6). According to the conclusion of Nygren,
2
we assume that the variance of elevation error of terrain nodes are independent and identically distributed
In this equation,

Eight discrete directions of the location terrain information content.
Discrete the results of equation (6) to get equation (7). The inverse of
In the equation, ei
represents the unit direction vector and will have eight values; d represents the distance of two terrain nodes when calculating the degradation of the local terrain. So, equation (7) reflects the accuracy of the TRP.
2
The greater the The gradient variation item The terrain gradient will change with the directions of The accuracy of TRP is related to the measurement error. The smaller the measurement error is, the higher the positioning accuracy.
From the above analysis, we get the influence factors and the quantitative evaluation of the positioning accuracy. Then, we get the quantitative function of the local terrain adaptation by further analysis.
Terrain suitability quantification
According to the conclusion of Nygren,
2
the Cramer–Rao bound of TRP

The
We refer to SNR as the suitability function of TRP, its value is the minimum value of

The likelihood of MTM1 (a) and MTM2 (b). MTM: measurement terrain map. MTM: measurement terrain map.
Equation (8) can calculate the suitability value of a terrain block with a size of
Optimal partition of terrain suitable area
Suitability quantization of DEM under arbitrary gridding
We are concerned with the problem of gridding the DEM with optimal partition of the suitable area and unsuitable area. The next step is to find the optimal partition method for the suitable region and unsuitable region in DEM. Figure 5 shows the division process of a DEM. Figure 5(a) to (d) show the distribution of suitable and unsuitable area during the process of mesh size changes.

An example of DEM partitioning process. (a) Divide the DEM into four blocks, (b) divide the DEM into nine blocks, (c) divide the DEM into 16 blocks, and (d) divide the DEM into 25 blocks. DEM: digital elevation map.
With changes in the mesh size, the suitable area and unsuitable area has been parted not different blocks or in the same block. The problem we care about is whether there is a mesh size that can divide the suitable area and the unsuitable area into different mesh as much as possible. Next, we will study the partition method. Our research object is the DEM as shown in Figure 6. The previous content has studied the method of terrain suitability quantification. On this basis, the following sections will introduce the partition method of the suitable area and the unsuitable area. First, we quantify the suitability of DEM according to the previous method. Then, we further discuss the solution of partition of the suitable area and the unsuitable area optimally. The suitability of DEM can be quantified as following steps:
1. Calculate the SNR value of node
2. The size of the mesh is defined by the number of nodes on the edge of mesh. The interval of the number of nodes on the edge of the mesh is
3. Obtain the quantitative expression of terrain suitability of the block

The DEM used in this article. DEM: digital elevation map.

Discrete information content for eight directions of DEM. DEM: digital elevation map.
Assumed that when the block size is
Partition of suitable and unsuitable area optimally
Now we can quantify the suitability of each block. Our goal is to divide the suitable area and the unsuitable area into different terrain blocks as far as possible. According to the previous analysis, the larger the
In the equation, K and L represent the rows number and columns number of map blocks; a
1 represents the average of SNR of map block under the blocks number of

A good partition indicate that the average of suitability of all high suitability blocks and the average of suitability of all low suitability blocks close to right and left end point.

Relationship curve between the number of block side nodes and the block evaluation index
According to the above requirements, we establish the constraint conditions for an optimal block. Based on the above analysis, a mathematical description of optimal block is obtained. In the process of optimal partitioning, we must to determine the range of the block sizes. The value interval is based on the measurement range of multi-beam sonar, supposing that the maximum and minimum values of the number of block side nodes are
According to the above analysis, we can obtain the optimal block of the DEM. Next, we calculated and analyzed an actual submarine topography as shown in Figure 6. The size of DEM is
Calculate the optimal block
As shown in Figure 10, the block results are obtained according to the 87 boundary nodes. Figure 10(a) shows the SNR of each subblock, while Figure 10(b) shows the DEM results.

DEM and its gridding results. DEM: digital elevation map. (a) Gridding of adaptation map; (b) Optimal gridding result of DEM.
TRP comparison and analysis of suitability block
The validity of the suitability region segmentation is shown below. The DEM map will provide the priori information and the MTM will be used to simulate the real-time surveying terrain. We evaluate the effectiveness of the optimal partition algorithm by comparing the TRP accuracy in high suitability blocks and low suitability blocks. In the DEM (Figure 6) and the MTM (Figure 11), they have surveyed in the same sea area as shown in Figure 12. The DEM survey was taken in 2012 and the MTM measurement time was October 2016. Figure 13 shows the connection of measurement equipment while surveying MTM. The measurement device and their parameters are shown in Appendix Table 1A.

The DEM map and the MTM map used in the simulation. DEM: digital elevation map; MTM: measurement terrain map.

Location of DEM. DEM: digital elevation map.

Topographic map surveying equipment and their connections.
Comparison of commonly quantitative parameters
We compare the quantitative method proposed in this article with the two commonly used terrain adaptation quantization methods (topographic entropy and topographic standard deviation).
10,18,19
The expressions of topographic entropy and topographic standard deviation are as follows: 1. Topographic entropy H
In the equation, the number of nodes that measure the terrain is
2. The standard deviation of terrain S
We use the MTM1 (Figure 3(a)) and MTM2 (Figure 3(b)) as a contrast experiment. From Figure 3(c) and (d) we can see that the gradient of MTM1 has a distinct directional distribution, but the direction of MTM2 is not obvious. The suitability of MTM1 and MTM2 are compared using the suitability parameters SNR, Terrain Entropy, and Terrain Standard Deviation, respectively. The test result shown in Figure 13, and the deviation average and standard deviation of TRP deviation shows that the suitability of MTM2 is better than that of MTM1. The SNR of MTM1 is smaller than MTM2, and it shows that the suitability of MTM2 is better than that of MTM1, but Terrain Entropy and Terrain Standard Deviation of MTM1 and MTM2 are almost the same. The major reason is that the Terrain Entropy and Terrain Standard Deviation can’t show the characteristics of anisotropy of terrain. As one of the most important features of terrain suitability is its directionality, and the most important feature of the suitability quantization parameter SNR used in this article is the ability to show the orientation of suitability.
Terms and abbreviations.
TRP accuracy of path in different suitability blocks
As shown in Figure 15(a), we choose suitable blocks and unsuitable blocks from the suitability map (Figure 15(a)). The unsuitable blocks and the suitable blocks are marked in white rectangle (marked A and B in Figure 15(a)) and red rectangle (marked C in Figure 15(a)), respectively. The suitability average of unsuitable blocks is 0.007, and the suitability average of suitable blocks is 0.019. The planning path (red line) and waypoints (green dots) in the suitable blocks and unsuitable blocks are shown in Figure 15(a). The simulated real-time measurement terrain is obtained by interpolation in the Figure 14 along the planning path as shown in Figure 15(a), and the simulated real-time measurement terrain are shown in Figure 15(b) and marked as a, b, c, respectively. To compare the position accuracy of high suitable blocks and unsuitable blocks under the condition of different measurement error. The design of comparative experiments is as follows:
The terrain blocks in red rectangles (marked as A and B) shown in Figure 15(a) are unsuitable terrain blocks. The terrain blocks in white rectangle (marked as C) shown in Figure 15(a) are suitable terrain blocks.
To simulate different measurement error conditions, the Gauss white noise added to the simulated measurement terrain map are
Repeat 10 times positioning simulation under the noise of

Comparison of position accuracy and quantization parameters of MTM1and MTM2. TRP: terrain reference position; MTM: measurement terrain map; SNR: signal noise ratio.

The data used in the simulation experiment. (a) The suitability map of DEM, high suitable blocks (C), low suitable blocks (A, B), and waypoint (•). (b) The counter map of DEM and planning measurement path (__), MTM. MTM: measurement terrain map; DEM: digital elevation map.
Finally, we get the statistical results of TRP positioning accuracy as shown in Figures 16 and 17.

Statistical accuracy of positioning accuracy in suitable blocks (SNR = 0.019) and unsuitable blocks (SNR = 0.007), when the measurement error is 0.4. SNR: signal-to-noise ratio.

Statistical accuracy of positioning accuracy in suitable blocks (SNR = 0.019) and unsuitable blocks (SNR = 0.007), when the measurement error is 0.9. SNR: signal-to-noise ratio.
From the results of Figures 16 and 17, we can see that when the suitability of the terrain block is SNR = 0.019 and the measurement error is 0.4, the average of TRP deviation is 3.69 m and the standard deviation of TRP deviation is 1.93 m. When the measurement error is 0.9, the average of the TRP deviation is 4.18 m and the standard deviation of TRP deviation is 1.96 m. The average of TRP deviation is increased by 13.3%, and the standard deviation of TRP deviation is increased by 1.6%. When the suitability of the topographic block is SNR = 0.007 and the measurement error is 0.4, the average of TRP deviation is 8.13 m and the standard deviation of TRP deviation is 5.62 m. When the measurement error is 0.9, the average of the TRP deviation average of the TRP is 9 m and the standard deviation of TRP deviation is 7.25 m. The average of TRP deviation is increased by 9.7%, and the standard deviation of TRP deviation is increased by 29%. The experiment shows that the TRP accuracy in suitable blocks will be better than that of unsuitable blocks. In addition, the standard deviation of TRP deviation in the suitable blocks will increase faster than that of the unsuitable blocks when the measurement error increased. That is to say, in the low suitable blocks, the TRP accuracy will be affected by the measurement error more easily.
Accuracy comparison between different suitability paths
The simulation process is shown in Figure 18. At first, the DEM is gridded as shown in Figure 10, and the high suitability area and the low suitability area are separated using the method proposed in the previous part. Then, we plan two paths (path A and path B), path A passes through the high suitability blocks and path B passes through the low suitability blocks. Next, we use the planning path (path A and path B) to interpolate in MTM (Figure 11), simulate the real-time measurement process of multi-beam, and get the simulated terrain along the path. Then, the simulated path (path A and path B), the simulated terrain, and DEM are put into the TRP system.

Flow chart of TRP simulation system. TRP: terrain reference positioning.
The planning path in the DEM is shown in Figure 19. Figure 19(a) and (b) show the location in the DEM and SNR map of path A, and the green dots are waypoints. As the same, Figure 19(c) and (d) shows path B. In addition, the block covered by path A will have bigger SNR values than those of path B.

Two tested planning path and waypoints for path A and path B. (a) Path A and the waypoints in the DEM, (b) path A pass through high suitability blocks, (c) path B and the waypoints in the DEM, and (d) path B pass through high suitability blocks. DEM: digital elevation map.
Figure 20 shows the TRP result of path A through the high suitability block. Figure 20(a) shows the distribution of simulation of real-time measurement map and TRP points (•), while the SNR values of the covered block are shown in Figure 20(b). The SNR values are greater than 0.02. Figure 20(c) shows the TRP deviation at each waypoint. From the figure, we can see that the TRP deviation is limited below

TRP results for the route through the high SNR block. (a) Path A and the simulated real-time measurement terrain at waypoints, (b) the suitability value of the blocks that path A pass through, and (c) the position deviation of the waypoints in path A. TRP: terrain reference positioning; SNR: signal-to-noise ratio.
Figure 21 shows the TRP results for path B pass through the low suitability block. Figure 21(a) shows the distribution of MTM and TRP points (•), while the SNR values of the covered block are shown in Figure 21(b). The SNR values are below 0.02. Figure 21(c) shows the TRP deviation at each waypoint. From the figures, we can see that the TRP stability of path B is lower than path A (Figure 19). In blocks 5 and 6, the TRP deviation increases sharply.

TRP results with the route through the high SNR block. (a) Path B and the simulated real-time measurement terrain at waypoints, (b) the suitability value of the blocks that path B pass through, (c) the position deviation of the waypoints in path B. TRP: terrain reference positioning.
Figure 22 shows a comparison of TRP results for the two routes. As we can see, the average and standard deviation for TRP deviation along path B (the low SNR path) is larger than that for path A (the high SNR path). Path A, through a terrain block with high SNR, and path B, through a terrain block with low SNR. From the TRP results, we can see that path A produces more accurate TRP positions, and the result is more stable.

TRP position comparing path A and path B. TRP: terrain reference positioning.
Conclusions
The main contribution of this article is to propose the terrain suitability parameter SNR and global map gridding method based on partition of suitable area and unsuitable area optimally. The difference with the commonly used quantitative parameters is that SNR not only takes the measurement error of terrain, the change in rate of terrain gradient into account, but also takes the directional characteristics of adaptation into account. The simulation experiment shown that the influence of measurement error on TRP accuracy is higher than that in high suitability area, and higher positioning accuracy and stability can be obtained in the high suitability blocks. Notably, the optimal partition criterion in this article is that the maximum distance between the average of high suitability blocks and the suitability of the low suitability blocks. This optimal criterion is not the only one and there may be a better criterion. In addition, the next goal is to establish a more perfect underwater environment information map, including underwater threats, task location, priorities, and so on. In this way, AUV can plan an optimal path independently according to the task requirements and environmental threats, navigation errors, terrain suitability, and other factors.
Footnotes
Declaration of conflicting interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.
Funding
The author(s) disclosed receipt of following financial support for the research, authorship, and/or publication of this article: This study is supported by the Central University Freedom Exploration Fund (Projects GK2010260235) and National Key Research and Development Plan (No. 2017YFC0305700).
Appendix 1
Major measure device parameters.
Compass
Topography measurement device
Device name
Harbin engineering university self-developed fiber-optic gyro
Geo Swath Plus
Specifications
Temperature bias stability: 0.05° h−1
Maximum water depth: 50 m
Random walk: 0.0025° h−1
Depth resolution: 1.5 mm
Angular resolution: 0.2 arc s
Max Swath update rate: 30 s−1
