Abstract
A high-resolution digital elevation model (DEM) is an important element that determines the performance of terrain referenced navigation (TRN). However, the higher the resolution of the DEM, the bigger the memory size needed for storing it. It is difficult to secure such large memory spaces in small, low-priced unmanned aerial vehicles. In this study, a high-precision terrain regression model to fit the DEM is generated using the extreme learning machine technique based on the multilayer radial basis function. The TRN results using the proposed method are compared with existing studies on various DEM fitting methods. This study verifies that the proposed method obtains improved fitting accuracy and TRN performance over existing DEM fitting methods such as bilinear interpolation, SVM for regression, and bi-spline neural network, without the DEM storage space.
1. Introduction
Terrain referenced navigation (TRN) is a navigation technology that estimates the position of an aircraft by comparing the altitude measured by an altimeter with the digital elevation model (DEM). Currently, the global navigation satellite system (GNSS) is commonly used for precise navigation [1] but an alternative is needed for environments in which the GNSS is unavailable due to interference such as hostile jamming. Thus, studies on TRN are actively being carried out as an alternative that will ensure precise navigation performance when GNSS is denied [2–5].
The DEM is a map that stores the elevation values of terrain to represent a terrain’s surface in three dimensions, and a high-resolution DEM is essential for improving the performance of TRN. However, large amounts of memory space are needed to store a high-resolution DEM and the loading time of the DEM in a separate storage space interferes with real-time computation. In addition, various interpolation methods are being used to obtain terrain elevation at the desired location using a discontinuous DEM. This method also generates additional computational errors, which can be larger as the resolution of the DEM decreases. To alleviate this problem, studies are continuously being conducted on devising a continuous terrain elevation model through a regression method based on the machine learning technique [6–10]. However, studies on DEM fitting using machine learning techniques up until now were conducted on DEM data in small areas of several hundred square meters or were unable to acquire higher accuracy than interpolation methods. Furthermore, excluding the method proposed by Liu et al. [6], there have been no studies that applied actual TRN using trained regression models. Also, in the case of Liu et al. [6], applying low-resolution DEM learning results in narrow areas would be inappropriate for actual high-performance TRN. In addition to the necessity of very high fitting accuracy for application in TRN, there is also the constraint of guaranteeing real-time operation. Taking into consideration such conditions, we applied the extreme learning machine (ELM) technique to train a DEM. The ELM has a simple structure and a fast learning speed, while it is still able to obtain high levels of accuracy and it is therefore currently being researched in various fields [11–14]. In addition, studies are actively being conducted on advanced algorithms to improve the performance of the conventional ELM or to resolve the complex problems beyond the abilities of the conventional ELM [15–22]. Kan et al. [9] and Hu and Yuan [10] estimated the DEM using a conventional ELM, but it was conducted on a very narrow area, and it was unable to obtain higher accuracy than the Delaunay triangulation-based interpolation method.
Even with the advantage of not needing large memory space, unless it has high levels of fitting accuracy that can replace the DEM, it cannot be applied to TRN. Moreover, when fitting a DEM as a complex model to enhance accuracy, it can preclude the real-time operation for navigation. Considering such restrictions, we used the generalized ELM autoencoder [15, 20] and multilayer radial basis function (RBF) [16, 18] to design a model that enhances the fitting accuracy of the ELM and guarantees real-time computing. We compared systems that used the proposed model and bilinear interpolation, which demonstrated the highest TRN performance in past studies, and our results verified that the proposed technique is the most suitable for TRN.
In the next section, we introduce existing studies on DEM fitting methods, namely, bilinear interpolation, the bi-spline neural network (B-spline NN), and the support vector machine for regression (SVMR). Then we introduce DEM learning techniques that use the conventional ELM and propose the multilayer RBF-based ELM (ML-RBF-ELM) to overcome the limitations of the conventional ELM. Furthermore, we compare the existing various DEM fitting methods and use of the proposed ML-RBF-ELM regression method and bilinear interpolation. Finally, we verify the design of the proposed technique by comparing it with the TRN performance.
2. Conventional DEM Fitting Method
There are different types of DEM database: the shuttle radar topography mission (SRTM) elevation database, digital terrain elevation data (DTED), etc. In this study, we used 10 m resolution DTED level 3 and 30 m resolution DTED level 2. There are various interpolation methods and machine learning regression techniques that can be used to obtain terrain elevation at the desired position using discontinuous DEM data. This section introduces previous studies that are aimed at obtaining continuous terrain elevation information.
2.1. Bilinear Interpolation
Interpolation methods frequently used to obtain continuous elevation information include Delaunay triangulation [9] and bilinear and bicubic interpolations. This study introduces bilinear interpolation, which is generally used in TRN systems. Bilinear interpolation is a method for linearly determining the linear distance of the four nearby DEM grid points, , , , and , when estimating the elevation value at any point within the four points as shown in Figure 1.
If the terrain elevation value estimated at the point is , the elevation values stored in DTED, , , , and , can be used to compute as follows.
2.2. B-Spline Neural Network
As discussed by Liu et al. [6], the B-spline NN-based regression method was proposed to fit the 30 m resolution DEM. The method divides all the data into multiple bins to supplement the weaknesses of polynomial regression, and it fits the divided data with separate polynomial functions. The two adjacent polynomials’ first and second differentials were made the same for smooth connection without making noncontinuous points in bins where division occurs. Liu et al. [6] designed the system as a single hidden layer where such spline basis functions are used for the layer’s nodes and DEM data were learned as a neural network based on this. In this study, the B-spline NN was designed as 1971 hidden nodes using the spline basis function of order 3. This is compared with other fitting methods in the next section.
2.3. Support Vector Machine for Regression
Okwuashi and Ndehedehe [7] used SVMR to fit the DEM in the 225 m by 225 m field. The regression method using SVMR was compared with the artificial neural network and nearest neighbour method, and it showed the highest level of accuracy. SVMR uses the same principles as the SVM for classification with only a few minor differences. For regression, a margin of tolerance was set as the estimate of the SVMR, which would have already been determined from the problem [23]. SVMR should consider more design factors than SVM for classification. However, like SVM for classification, SVMR is learned to minimize the difference between the true value and the predicted value and maximize the margin between decision boundaries. For nonlinear SVMR, kernel functions map input data to higher dimensional feature space, simplifying it as a linear SVMR [23]. In our study, the Gaussian kernel function was used for the precision DTED regression model based on SVMR. In addition, the box constraint was used and set as 0.002178 to prevent overfitting. The optimization issue of the SVMR model is usually solved with the Lagrange dual formulation for computational convenience, and when the Karush-Kuhn-Tucker (KKT) complementarity conditions are satisfied, the optimal solution can be obtained. In this study, such KKT tolerance was set as and it was optimized using the iterative single data algorithm (ISDA), which updates a single multiplier per iteration.
3. DEM Fitting Method by ML-RBF-ELM
3.1. Conventional ELM
Since first being developed by Huang et al., the ELM has demonstrated good generalization ability and much faster learning than the backpropagation method of stochastic gradient descent and is therefore being studied in various fields. The ELM is a single hidden layer feedforward neural network. The input weights and hidden unit biases are chosen randomly and do not need to be adjusted [24]. As the output weights are analytically computed by the least squares method, this method provides a good generalization performance at an extremely fast learning speed. For distinct training sample , the output of ELM with hidden layer nodes is as in [8, 9]:
Here, is the activation function of the th hidden node, which is modelled as a RBF. is the th input data and is set to 2 in this study. and are the normalized longitude and latitude, respectively. is the th target (or true) data and means the normalized terrain height of the DEM in this study. and are the center and width, respectively, of the th RBF. By Huang et al. [24], the center and width are assigned randomly. means the Euclidean distance. is the number of hidden nodes and is set to 271 including a bias node in this study. The number of nodes was optimized to have the lowest training and test errors. The ELM can estimate any target function with zero error using equation (2). The relation is written as
Here, is the hidden layer output matrix and is the hidden layer set of the th sample data. The hidden layer set consists of activation functions. is the target matrix of the ELM, and is set to 1 because the normalized terrain elevation data is used as the true value of the training data. is the weight matrix of the output layer and the smallest norm least squares solution of equation (3) as follows:
Here, is the Moore-Penrose generalized inverse, and it is generally calculated as a singular value composition. In order to prevent the overfitting and improve the robustness, a ridge parameter is used [15, 16]. The optimization problem is described as
Here, is a regularization parameter and is set to 0.2. is the residual between the target value and estimated output value of the th sample. According to KKT conditions, the solution of equation (5) is calculated as follows [9]:
Here, is a -by- identity matrix.
4. ML-RBF-ELM
The conventional ELM is easy to configure and the learning speed is fast, but it cannot be applied to fields that process big data with very large nonlinearity or require high levels of accuracy. To solve these problems, many studies are being actively conducted on various advanced ELM algorithms. In the constraint in which real-time computation must be guaranteed, three hidden layers were used as shown in Figure 2 to maximize the training and test accuracies. The first hidden layer was designed as an ELM-based autoencoder (ELM-AE), and the remaining hidden layers were designed as conventional ELMs. As the main difference between the ELM-AE and conventional ELM, the ELM is a supervised NN and the outputs are the target values but the ELM-AE is an unsupervised NN and the outputs are equal to the inputs [15].
The output weight matrix of the 1st hidden layer output is calculated as follows [15]:
Here, and . is a -by- identity matrix and is the number of nodes in the 1st hidden layer. In this study, the ELM-AE was used to make the inputs more important in a rough terrain area than in a flat area. The 2nd and 3rd hidden layers shown in Figure 2 are similar to the conventional ELM. These layers were designed to represent the terrain height with refined accuracy. Unlike the conventional ELM, in the proposed ML-RBF-ELM method, the weights between the previous hidden layer and the current hidden layer were determined using K-means algorithms to enhance the training accuracy. The 1st and 2nd intermediate vectors are virtual vectors that represent the process in which the output weight matrices of the 1st and 2nd hidden layers explained in Algorithm 1 are learned. K-means clustering is used to partition sample data into clusters based on features in which each sample data belongs to the cluster with the nearest mean of a Gaussian function [25]. In this study, the mean and width values of the RBF nodes in the hidden layers are determined by the K-mean clustering based on the expectation-maximization algorithm [25]. The process of the proposed ML-RBF-ELM is given in Algorithm 1. The numbers of nodes in the three hidden layers, , , and , are set to 384, 2301, and 193, respectively.
| ||||||||||||||||||||||||||
4.1. Comparison of Various DEM Fitting Methods
Figure 3 illustrates the DEM fitting area with a total area of about km2 for DTED level 3 and km2 for DTED level 2. As mentioned above, we used 10 m resolution DTED level 3 and 30 m resolution DTED level 2 for the training and test data. The total obtainable number of sample data were 5391738 sets (161 MB) for level 3 and 5443200 sets (163 MB) for level 2. At a rate of 8 : 2, the sample data was used by dividing the training data (4313390 sets, 129 MB for level 3, and 4354560 sets, 130.4 MB for level 2) and test data (1078348 sets, 32 MB for level 3, and 1088640 sets, 32.6 MB for level 2). The sample data are acquired from the grid points such as , , , and in Figure 1. The learning was performed with an Intel® Core™ i7-6700K, CPU @4.00 GHz, RAM 16.0 GB computing environment. All simulations for the bilinear interpolation, B-spline NN, SVMR, conventional ELM, and the proposed ML-RBF-ELM were carried out in a MATLAB R2018a environment.
(a)
(b)
Table 1 compares the fitting accuracy and operation time according to the bilinear interpolation, B-spline NN, SVMR, conventional ELM, and the proposed ML-RBF-ELM for fitting DTED level 3. The fitting results of level 2 are also similar to those of Table 1. Each method was designed to maximize fitting accuracy. For the ML-RBF-ELM design, which was designed with three hidden layers, as shown in Table 1, the results confirmed that compared to the B-spline NN and conventional ELM designed as a single hidden layer, more nodes could be used in each hidden layer. Furthermore, we also confirmed that the ML-RBF-ELM, which was designed with three hidden layers, learned faster than the B-spline NN, which was designed as a single hidden layer. Likewise, it was evident that the ELM performed better in a shorter amount of time than the NN. The test time in Table 1 indicates the time needed for fitting one new sample to the terrain height. For the index to evaluate the fitting accuracy, the following mean value of the absolute percentage error was used. The value is calculated as follows:
Here, and are the estimated and target output value of the th sample data, respectively. is the total number of samples.
The training and test errors of the bilinear interpolation in Table 1 are approximate values. These values are about 25% of the absolute percentage errors when the terrain elevation value estimated at point was acquired from points , , , and in Figure 1. The conventional ELM had the smallest training and test time but the lowest fitting accuracy so it is not suitable for TRN. The SVMR did not satisfy the test time restrictions for real-time operation, so it is unsuitable. However, the proposed ML-RBF-ELM had the best fitting accuracy and its short test time allowed for real-time operation, so we judged that it would be suitable in TRN. Figures 4 and 5 are illustrations of the terrain height fitting results for the test data. We found that previous research results did not reach the accuracy of bilinear interpolation, as shown in Figure 4. Despite the advantages of not requiring a large memory space to store the DEM and the time needed for loading, if fitting accuracy results cannot be obtained at the same levels of bilinear interpolation, then the performance of TRN cannot be guaranteed. On the other hand, as evident in Figure 5, the proposed ML-RBF-ELM demonstrates a fitting accuracy similar to the bilinear interpolation of Figure 4. To store the weights of the proposed ML-RBF-ELM, we need a 105 KB size memory. This is only about 1/1533 of the 161 MB memory size for DTED level 3 and about 1/1552 of the 163 MB memory size for DTED level 2.
(a)
(b)
(c)
(d)
(a)
(b)
5. Performance of TRN Using ML-RBF-ELM
5.1. BKF-Based TRN and INS/TRN Integration
We verified the ML-RBF-ELM by conducting TRN with the ML-RBF-ELM and bilinear interpolation, which is superior to the fitting methods explained in the previous section. TRN simulations were performed using the interferometric radar altimeter (IRA) measurements and a bank of Kalman filters (BKF), which were used in Heli/SITAN [26]. Figure 6 presents a schematic diagram of our system with an inertial navigation system (INS), IRA, barometer, DEM, TRN S/W, and INS/TRN S/W. The TRN progresses in two steps, the time propagation and the measurement update. If the IRA measurements do not satisfy validity check conditions; the TRN performs only the time propagation. To acquire all the navigation information, including the position, velocity, and attitude, the INS/TRN-integrated navigation is needed. Therefore, the loosely coupled INS/TRN-integrated navigation was designed with an extended Kalman filter (EKF) and uses the position estimated by the BKF-based TRN for measurement. Details of the INS/TRN-integrated navigation are given in section 4.3 of reference [27]. If the estimates of TRN do not satisfy terrain validity check conditions, the INS/TRN performs only the time propagation. Table 2 shows the simulation conditions and BKF design parameters. On June 21, 2018, and January 21, 2019, two captive flight tests were carried out on a CESSNA 172 Skyhawk. The INS/GNSS-integrated navigation data and IRA outputs were used as the true trajectory and measurements in this simulation, respectively. We conducted Monte Carlo simulations 100 times using bilinear interpolation and the proposed ML-RBF-ELM. The proposed ML-RBF-ELM was implemented by using a S/W without the DEM, as shown in Figure 6.
The BKF is a linear filter that uses the multiple-model adaptive estimation (MMSE) technique to acquire navigation information sequentially, even in the case of a large initial position error. A state variable of one-dimensional Kalman filters comprising the BKF is the vertical bias error, which gradually changes. The BKF forms grids that are centered on the INS/TRN position estimate and assigns a Kalman filter to each grid. In this study, filters in the grid region corresponding to three times the initial position error were allocated. The time propagation for the th () filter state variable and its covariance estimate assigned to each grid is described as
Here, is the process noise covariance and is the sampling time. is the upward velocity of the aircraft in the current step. In this study, the IRA was used to measure the angle perpendicular to the direction of flight and the range from the aircraft to the nearest terrain. It then converted these measurements to three-dimension position information. As shown in Figure 7, the relative position vector of the nearest point from the aircraft is calculated as follows:
Here, , , and are the relative distances in the directions of longitude, latitude, and height, respectively. and are the range and look angle output, respectively, of the IRA. The virtual pitch angle and azimuth angle of the zero Doppler line are determined by the velocity of the aircraft [27].
The th filter gain and the measured values using the IRA measurement are as follows:
Here, and are the longitude and latitude, respectively, estimated by INS/TRN and is the measurement noise covariance. and are the longitude and latitude, respectively, at the nearest point as shown in Figure 7. and represent the radius of the curvature of the Earth ellipsoid in the north-south and east-west, respectively. is the terrain height of the th filter that is centered on the nearest point . The terrain height is calculated using the bilinear interpolation from the stored DEM or the regression model by the proposed ML-RBF-ELM. is the aircraft altitude measured by the barometer. The weighted residual squares (WRS) of the th filter at the th step are calculated as follows [26]:
The WRS indicates how well each filter in the grids matches the designed model. The smoothed WRS (SWRS) is expressed in the same way as the moving average of the WRS as follows [26]:
Here, and is the smoothing constant () and is determined by calculating how fast the INS estimate gets out of the filter exclusion region by the INS navigation error. In this study, was set as 0.018. The estimated covariance of each SWRS value in the filter exclusion region centered on the grid with the minimum SWRS value among total grids is as follows:
Using the estimated covariance, the longitude and latitude and estimated by the BKF are calculated as follows:
Here, and are the initial longitude and latitude, respectively. and represent the east and north axis position estimates of the th filter in the exclusion region. The terrain validity check conditions in Figure 6 are as follows:
Here, is the minimum SWRS value outside the beta area and represents how many times the filter retains the minimum SWRS value. and are their respective thresholds. Since the smallest value in the grid is , the terrain validity check conditions are always positive and the larger the difference, the rougher the terrain. The more often the above equation is satisfied, the more likely it is that the terrain is unique. The thresholds were each set as and 2 for the rough area and and 29 for the flat area through numerical simulations.
6. Simulation Results
Figure 8 compares the case of using the bilinear interpolation in the same simulation conditions (dotted lines) and the case of using the regression model fitted according to the proposed ML-RBF-ELM (solid lines). The position root mean square (RMS) errors in the north-south and east-west axes were smaller for the BKF-based TRN using the ML-RBF-ELM. Figure 8(a) shows such characteristics more clearly between 400 seconds and 500 seconds. Figure 8(b) also shows improved performance from 500 seconds onward. The 100 Monte Carlo simulation results of INS/TRN-integrated navigation, which uses the latitude and longitude information estimated by TRN as measurements, are shown in Figure 9. The performance difference between the bilinear interpolation and the ML-RBF-ELM in INS/TRN is smaller compared to TRN. However, it is evident that compared to using bilinear interpolation in Figures 9(a) and 9(c), using ML-RBF-ELM carries out TRN and INS/TRN more stably as shown in Figures 9(b) and 9(d). If the terrain height changes smoothly, the error of bilinear interpolation is small. However, the bilinear interpolation error may be bigger in terrains with higher nonlinearity. Meanwhile, the proposed ML-RBF-ELM is less affected by such terrain roughness or uniqueness and so the change in the overall RMS error may be smaller than bilinear interpolation. Table 3 shows the values of converting the horizontal position RMS errors of TRN and INS/TRN to circular error probability (CEP). Comparing TRN with the bilinear interpolation, TRN with the ML-RBF-ELM improved the result by about 2.335 m CEP in DTED level 3 and about 1.402 m CEP in DTED level 2. In other words, using ML-RBF-ELM, it is possible to fit a highly accurate regression model without the DEM. Moreover, TRN and INS/TRN using the proposed method achieved stable performance.
(a)
(b)
(a)
(b)
(c)
(d)
7. Conclusions
This study introduced DEM fitting methods, a key technology of TRN, that can be used as an alternative in environments where GNSS is unavailable. The DEM is a database that stores the terrain’s elevation data at a constant resolution. Therefore, large amounts of memory are needed to use the high-resolution DEM. Furthermore, the various interpolation methods for obtaining continuous elevation information from the discontinuous DEM includes fitting errors. Accordingly, there are currently many ongoing studies on DEM learning through various machine learning techniques. However, these efforts have only produced results on DEM learning in low-resolution or narrow fields and they have yet to achieve stable performance for TRN using these fitting methods. In this study, we used the ML-RBF-ELM to fit a 10 m resolution DEM in an area of km2 and a 30 m resolution DEM in an area of km2. The proposed regression model by the ML-RBF-ELM operates more quickly than the previous machine learning regression methods for the DEM and demonstrates similar fitting results with bilinear interpolation. Furthermore, we used BKF-based TRN, which showed improved TRN performance results over using the bilinear interpolation. Thus, our study verified that the proposed ML-RBF-ELM can make real-time operation of navigation and is suitable to TRN. It is clear that the proposed technology can be used by low-priced small unmanned aerial vehicles to execute TRN in environments where GNSS is not possible. Moreover, the proposed ML-RBF-ELM can be used in the field of simultaneous localization and map building (SLAM) because of the small training time.
Data Availability
The DTED level 3 full sample data and terrain reference navigation S/W used to support the findings of this study are not freely available, because of our institutional security policies. But the edited sample data, MRBF-ELM fitting S/W, and the simulation result data used to support the findings of this study are included within the supplementary information files (available here). The edited sample data size is smaller than original sample data size.
Conflicts of Interest
The authors declare that there is no conflict of interest regarding the publication of this paper.
Supplementary Materials
The list of the supplementary information files is as follows. The edited sample data: “TestDBData3.dat,” “TrainDBData3.dat” in “TrainDB_ML_RBF_ELM_Using AE2_180919” Folder MRBF-ELM fitting S/W: “TrainDB_ML_RBF_ELM_Using AE2_180919” Folder the simulation result data: “TRNSimulationOutput_180912” folder, and “CompareMatchingResults_180909” folder. (Supplementary Materials)