An advanced error state Kalman filter (ESKF)-based terrain contour matching (TERCOM) method for tracking an aerial vehicle using a low-cost digital elevation map

State of the art review

The terrain contour matching (TERCOM) is a class of TRN algorithm which uses pre-stored terrain elevation data to determine the current location of the aircraft. The history of the TERCOM system dates back fundamentally to the design of a cruise missile system in the late 1950s. The cruise missile system is a pilotless, continuously powered and a dispensable vehicle, which is specifically used for delivering the nuclear devices. One of the major advantages of the TERCOM technique is that it can achieve faster localization with accuracy in cases of large initial positioning error (Peng et al., 2018). The disadvantage is that the technique has poor real-time performance. Another TRN algorithm is the Sandia Inertial Terrain Aided Navigation (SITAN), which uses the Extended Kalman Filtering (EKF) method (Jayaramu et al., 2021). SITAN provides the capability of real-time positioning, however, the major disadvantage is that it produces high false-fixes in cases of large initial positioning errors (Dai & Kang, 2014).

In the recent years, hybrid TAN strategies have also been reported combining TERCOM and Sandia Inertial Terrain Aided Navigation (SITAN) algorithms to exploit the advantages of both the technologies in order to meet the military requirements. The Terrain Profile Matching (TERPROM) is a hybrid two-phase TRN method which uses TERCOM strategy in the acquisition mode and SITAN algorithm in the tracking mode (Eroglu & Yilmaz, 2014). Apart from position estimation, the system offers a number of features including Predictive Ground Collision Avoidance (PGCA), Obstruction Warnings (OW), Terrain Following (TF), and Passive Ranging (PR) (Lee & Bang, 2018a). The disadvantage of the scheme is that it takes longer computational time due to batch processing in the TERCOM method.

The accuracy and reliability are one of the key factors characterizing a navigation system. In Yoo et al. (2012), an improved TERCOM algorithm using the velocity correction method was proposed. The observability analysis was performed and a selective error method was employed to correct the velocity errors in the INS in addition to the position errors. In Zhao (2012), a TERCOM strategy based on three stage logic was proposed in order to mitigate the false-fix issues in the TERCOM algorithm. In Wu, Fei & Li (2012), a multi-path Terrain Contour Matching algorithm using the laser scanning altimetry replacing the traditional radar altimeter sensor. The results of the scheme verified that 2D altimetry can effectively meet the localization requirement of the UAV. The disadvantage is that the scheme is applicable for low altitude UAV flight missions only. In Yan et al. (2018), a TRN system based on particle filter (PF) and adaptive scaling method for Digital Elevation Maps (DEM) for improving the navigation precision of the UAV was proposed. The system is capable of providing good navigation accuracy in cases of high UAV altitudes. The disadvantage is that the performance of the system degrades in case of flat and repetitive terrain. Two major issues with the PF-based approach are particle degradation and particle impoverishment (Yousuf & Kadri, 2024). In order to mitigate this issue, a TRN method based on intelligent particle filter (IPF) was proposed in Chai, Li & Qiao (2022). In this method, genetic algorithm (GA) was incorporated in the PF resampling stage to improve the robustness of the IPF. In Veselý et al. (2022), a system based on point mass filter (PMF) was proposed to deal with the non-linear filtering problem arising in a TRN based navigation. The PMF-based TRN showed high estimation performance, however the major disadvantage was the computational efficiency of the method. To overcome the drawbacks of the conventional PMF, the authors in Duník et al. (2019) proposed a TRN system based on Rao-Blackwellised point-mass filter. The system preserves the advantages of point mass filter providing high estimation accuracy and predictable computational complexity. Similarly, in Park & Park (2020), a Two Stage Point Mass Filtering (TSPF) state augmentation approach for TRN application was presented. In this method, the non-linear states, i.e., the UAV latitude and longitude were estimated using the PMF whereas the linear state variable, i.e., the altitude was estimated using a Linear Kalman Filter (KF). The results demonstrated that the TSPMF technique was more computationally efficient compared to the RBPMF method while the estimation performance was almost similar to RBPMF (Veselý et al., 2022). In Liu, Wang & Yao (2014), a TRN method based on Extended Kalman Filter and B-Spline Neural Network was presented for UAV navigation. The combined strategy provided better estimation results compared to EKF based TRN. The disadvantage of the BSNN network is that it is trained by gradient-based methods which may fall into local minima during the learning process. In Eroglu & Yilmaz (2014), a TERCOM method utilizing the Long Short Term Memory (LSTM) deep learning scheme was introduced. The LSTM network was trained to learn all possible flight data in order to estimate the current location of the UAV. A map reader approach which uses a filtered map was employed in training the LSTM. However, the authors report issues in the presented scheme from the aspect of real-world implementations. Similarly, in Lee & Bang (2018b), a TRN scheme combining Rao-Blackwellized Particle Filter (RBPF) and LSTM network was presented to minimize the navigation accuracy of a UAV. The scheme provided precise navigation compared to the conventional RBPF based TRN approach. The disadvantage was that study related to the real-time implementation of the method was not presented. More recently, in order to improve the position accuracy of a TRN system, the concept of fuzzy logic has been proposed. In Liu, Zhang & Huang (2021), a TRN was proposed in which fuzzy logic was applied to estimate the distribution of variance of particles in a Particle Filter to improve the navigation accuracy of an Autonomous Underwater Vehicle (AUV). It was reported that the proposed method in this work achieved good tracking accuracy. In Zhao, Xu & Weiqiang (2021), a multi-index grey fuzzy decision making method was proposed for terrain navigability analysis for an AUV. In this method, a point mass filter (PMF) was used to simulate the navigability of the matching areas having different terrain characteristics. The results of this method indicate feasibility of the proposed fuzzy evaluation scheme in selecting the underwater matching areas for better AUV navigation.

Furthermore, A robust Kalman filter (RKF) that handles measurement outliers through kernel density estimation (KDE) is presented in this article in Gao et al. (2025). A new RKF is created by using Bayesian estimation and a logarithmic Gaussian kernel to model sudden noise changes. Experiments and simulations demonstrate how well it works to increase vehicle navigation accuracy in outlier scenarios. A set-membership hybrid Kalman filter (SM-HKF) is proposed in this article in reference (Zhao et al., 2020) for nonlinear state estimation under both unknown but bounded (UBB) and stochastic errors. It obtains an optimal Kalman gain that takes into consideration all uncertainties by linearizing the system and combining errors using the Minkowski sum. Simulations demonstrate that when it comes to managing both systematic and stochastic errors, SM-HKF performs better than the extended Kalman filter (EKF). In Gao et al. (2021), the authors suggest a distributed optimal fusion approach based on cubature rules for UAV navigation with integrated miniature inertial measurement unit (MIMU), global navigation satellite system (GNSS) and celestial navigation system (CNS). By detecting and forecasting kinematic model errors using Mahalanobis distance, it tackles issues caused by nonlinearity. Subsystems employ modified cubature Kalman filters, and the outputs are combined to estimate the globally optimal state. Experiments and simulations verify increased navigation accuracy. In order to address measurement uncertainty in nonlinear estimation, the article in reference (Hu et al., 2024) suggests an indirect fuzzy robust Cubature Kalman Filter (CKF). It modifies CKF via a scaling matrix using a fuzzy inference system (FIS) with normalized inputs and triangular membership functions. This improves convergence and robustness. Simulations demonstrate better results than conventional fuzzy-based filters. An advanced cubature information filter for tracking multiple wideband sources indoors is presented in the article in Zhang et al. (2019). To increase tracking accuracy in difficult environments, it integrates a distributed noise statistics estimator. The technique provides a reliable indoor source tracking solution and improves performance in noisy environments. A novel Gaussian filtering technique is presented in the article cited (Gao et al., 2020), which overcomes the drawbacks of conventional methods that depend on known system noise characteristics. Utilizing maximum, a-posterior theory-based random weighting estimations, the approach adaptively modifies weights to reduce the effects of noise. Experiments and simulations demonstrate that this strategy increases estimation accuracy over conventional techniques. To overcome the drawbacks of conventional Kalman filters, which rely on known system models and noise statistics, the authors in Xia et al. (2020) suggest the fitting H-infinity Kalman filter (FHKF). FHKF increases the robustness and stability of nonlinear uncertain systems by estimating coefficients via least weighted squares using a fitting transformation. Its exceptional performance in preserving accuracy and stability is demonstrated by simulations and real-world INS/GPS experiments.

Also, an innovative orthogonality-based robust unscented Kalman filter (IO-RUKF) for hypersonic vehicle navigation is proposed in this article in Hu et al. (2019). By employing hypothesis testing to identify anomalies and a robust factor to modify the Kalman gain, it mitigates measurement errors such as outliers and non-Gaussian noise. Simulations confirm that the method increases unscented Kalman filter (UKF) robustness without compromising accuracy. Using the maximum likelihood principle, this article in Hu et al. (2020) suggests an adaptive UKF for vehicular INS/GPS integration that estimates process noise covariance. It improves robustness against process noise uncertainty by updating noise estimates online through the introduction of a fixed-length memory window. Experiments and simulations verify that it performs better than the typical UKF. In order to manage system nonlinearity and uncertain noise statistics, this article in Meng et al. (2016) suggests an adaptive UKF for direct INS/GNSS integration. It estimates and updates process and measurement noise covariances online using covariance matching. Compared to standard and adaptive-robust UKF methods, simulations and experiments demonstrate increased accuracy and robustness. In Hu, Gao & Zhong (2015), a refined strong tracking unscented Kalman filter (RSTUKF) is used to propose a direct INS/GNSS integration method. Through assumption testing, it detects kinematic model errors and modifies the predicted covariance using a suboptimal fading factor. In the absence of model error, this maintains UKF accuracy while improving robustness. Experiments and simulations verify its efficacy.

Contribution of the work

This article presents a novel integration of fuzzy logic with the error state Kalman filter (ESKF) to significantly enhance the positioning accuracy and reduce the computational cost of the TERCOM navigation algorithm. Unlike traditional approaches that rely on full DEM scans for terrain correlation, we introduce a fuzzy-based heuristic mean absolute deviation (FH-MAD) method that intelligently restricts the search space using real-time aircraft dynamics, leading to faster and more efficient localization. Specifically, the key contributions of this work are as follows:

• Fuzzy-Enhanced Terrain Matching: we develop a fuzzy logic-based extension of the MAD algorithm (FH-MAD), which adaptively prioritizes terrain correlation regions based on flight parameters, resulting in improved positioning decisions.

• Region-of-Interest (ROI) Selection using Flight Attitude: the proposed method dynamically selects a rectangular sub-region of the DEM using aircraft yaw and roll angle conditions, significantly reducing the computational burden compared to exhaustive DEM searches.

• Improved Localization Accuracy: the FH-MAD integrated with ESKF demonstrates superior accuracy in UAV localization, validated through rigorous performance metrics under varying flight scenarios.

• Significant Computational Gains: our method achieves a notable reduction in execution time when compared with the conventional MAD-based TERCOM algorithm, making it viable for real-time onboard processing.

• Robust ESKF Design for Maneuvering UAVs: we provide a carefully designed and tuned ESKF framework tailored to handle dynamic flight conditions, ensuring consistent performance during complex UAV maneuvers.

• Real-Time Application Readiness: the proposed approach is well-suited for real-time UAV navigation in GPS-denied environments, offering a practical and scalable solution for autonomous aerial systems.

The rest of the article is structures as follows: in “System Overview”, a brief overview of the proposed TRN method is presented. In “Proposed Method”, a detailed discussion on the proposed hybrid strategy based on fuzzy-based heuristic and ESKF is provided. The simulation tests and analysis is provided in “Simulation Results”. Finally, the conclusions are given in “Discussion”.

Data collection

For the proposed Fuzzy Logic Correlation Scheme, the X-Plane-11 software was used for the aircraft flight simulation. X-Plane is one of the most advanced flight simulators used in industry as well as for research purpose. X-Plane differentiates itself from other simulators by implementing an aerodynamic model called blade element theory (Yu et al., 2020; Gusti Agung Agastya & Kusumoputro, 2019). Traditionally, flight simulators emulate the real-world performance of an aircraft by using empirical data in predefined lookup tables to determine aerodynamic forces such as lift or drag, which vary with differing flight conditions. The Skywalker X8 model was designed in the X-Plane and was used to generate the flight test data. The sensor data obtained from X-Plane was exported in the MATLAB environment for further processing. The Mission planner software was used to specify the waypoints and fly the mission in the X-Plane as shown in Fig. 2. The Digital Elevation Map (DEM) of the same area was obtained from NASA Earthdata (n.d). An example of the listing of all the eighteen variables of interest generated during the flight simulation is shown in Table 1 below. These variables are the ‘inputs’ to the proposed navigation scheme.

INS navigation algorithm and transformations

The Inertial Navigation System is the combination of three major modules, i.e., (1) the inertial measurement unit (IMU) (2) the co-ordinate transformation module and (3) the two-step integrator module. An IMU measures the three Degree of Freedom (DoF) angular rate and magnitudes of specific forces (Shaukat, Moinuddin & Otero, 2021). These measurements are initially obtained in the body frame (i.e., the sensor frame) axes. The INS model selected is based on strap-down INS (SINS) technology. In this scheme, the integration is performed using the trapezoidal integration (TI) method to obtain aircraft position. In this work, the east-north-up (ENU) co-ordinate system is selected as the world co-ordinate system. The conversion from ENU frame to the latitude-longitude-altitude (LLA) frame is a two-step process (1) conversion from ENU to Earth-Centered, Earth-Fixed (ECEF) and (2) conversion from the ECEF to LLA frame. The transformation from ECEF to LLA frame can be performed by the non-iterative transformation method defined in You (2000) or the close form method called the Jijie Zhu’s algorithm (Osen, 2017). Moreover, it should be noted that the IMU measures the relative acceleration, also known as the specific force. It does not measure the true kinematic acceleration directly due to presence of the earth’s gravitational acceleration effect modelled as:

(1) $$f_i{b}^b={\dot{R}}_i^b{\ddot{p}}_i-g^b.$$

In the above equation, ${\ddot{p}}_i^b$ represents acceleration in the inertial frame and $g^b$ represents the gravity in the body frame. ${\dot{R}}_i^b$ is the rate of change in the rotation matrix represented in the body frame (b). And $f_{ib}^b$ represents the specific force.

DEM processing and conventional correlation scheme

In the proposed framework, the Mean Absolute Deviation (MAD) algorithm is used for the purpose of matching, where the matching area selection is based on Fuzzy Logic approach detailed in the next section. The MAD algorithm correlates the terrain heights estimated from the system with the DEM based topographical values. This correlation is described by the following equation:

(34) $$MA{D}_{n,m}=-\frac{1}{N}\sum _{i=1}^N|h_{meas}^{(i)}-h_{DB}(n+i,m)|.$$

In the above equation, $h_{m}eas$ and $h_{DB}$ are representing the measured terrain heights and the stored terrain database respectively. The variable ‘N’ represents the number of samples. Also, n and m represent the rows and columns of the terrain matrix stored in the digital terrain database. The position fix is obtained using the minimum MAD value given by:

(35) $$P_{MAD}(n,m)=P(\arg \underset{n,m}{min}(MA{D}_{n,m})).$$

Applying the correlation process, the terrain height with the minimum values are considered as appropriate candidates. The steps comprise of comparing and matching the difference between the barometer and the radar altimeter measurements (which is essentially the terrain height) with a pre-stored DEM. The DEM models topography as well as the topographical change in a geographical region (Wilson & Ramirez-Serrano, 2014). The DEM is structured as a two-dimensional matrices of elevation values assigned as the DEM matrix Z. In a DEM cell, the geographical location is computed according to a specific co-ordinate system assigning a projected location in the two-dimensional XY plane. The terrain roughness and uniqueness factor is also crucial for the successful application of the TERCOM system. TERCOM results are more appreciable for unique terrains (Chen, Xu & Ding, 2022). The terrain roughness characteristics of the region for which the simulations were conducted are defined by the two parameters, i.e., the Sigma-Z and Sigma-T values. These parameters are employed as the terrain roughness indicators defined in Lu, Jian & Xiaowen (2019), Raković, Simonović & Grbović (2020). The terrain roughness described by Sigma-T is represented by the standard deviation of terrain elevation values. Areas with very flat terrains and lakes have smaller Sigma-T values. Similar to Sigma-T values, Sigma-Z represents the roughness in the terrain. Another parameter is the terrain correlation length. The terrain correlation length describes the separation between two rows or columns of the terrain matrix which are required to reduce their normalized auto-correlation function to e−1 = 0.3697. In this article, only the first two parameters, i.e., Sigma-T and Sigma-Z are adopted for the analysis of terrain roughness characteristics.

Fuzzy logic based ROI selection for TERCOM

The conventional mean absolute devotion (MAD) algorithm is found to be computationally expensive in the traditional TERCOM-based strategy. In the proposed method, a heuristic based fuzzy logic system (FLS) has been designed to select the matching area in order to improve the shortcomings of the conventional MAD algorithm. The matching area selection is performed using the aircraft heading and roll angle parameters as inputs to the FLS. Since in practical situations, not only the yaw angle can change during the actual flight, the roll angle may also vary which affects the matching region due to aircraft maneuver.

The FLS module calculates the co-ordinates of the rectangular region (defining the matching area) based on the instantaneous heading angle (Y) and roll angle (R) information from the aircraft orientation sensor. The proposed fuzzy logic system takes the yaw (Y) and roll angle (R) as inputs and provides the co-ordinates (E1, E2, N1, N2) that define the rectangular region of the matching area based on the established fuzzy rules. The inputs and output membership functions are selected to be triangular membership functions. Because of their simplicity, effectiveness, and speed of computation, triangular membership functions are suited for real-time systems with constrained resources (Varshney & Goyal, 2023). They are simple to interpret and tune due to their intuitive design and low memory requirements. Control and decision-making tasks can benefit from these functions’ ability to accurately approximate gradual transitions between fuzzy sets and their compatibility with fuzzy rule-based systems due to their simplicity (Khairuddin et al., 2021; Samonto et al., 2021). The FLS system is designed as a 2-input and 4-output multi-input-multi-output (MIMO) system. Figure 3 shows the relationship between the inputs variables (Y and R) and output variables (E1, E2, N1 and N2). For the proposed system, the domain of the yaw angle parameter is considered from [0, 360] and the aircraft roll angle is considered within the domain [−45, 45]. In Fig. 3, the four regions for aircraft maneuver are based on the aircraft heading angle parameter. For each region, a separate FLS module is designed and stored in the system. In Fig. 3 the Region-I of the aircraft motion is defined from 0 to 90 degrees with respect to the yaw angle. The input membership functions for this case are shown in Fig. 4. In Fig. 3, the parameter ‘L’ represents a parameter related to the terrain characteristics over which the aircraft is flying. For simplification, this parameter is set to a constant parameter in this article. Based on parameter L and yaw and $\mathbf{\text{θ}}$, the computation of the output variables E1, E2, N1 and N2 is possible, however, the uncertainty arises when there are changes in the roll angle parameter as well as sensor noise parameters. Therefore, while constructing the output membership functions (MFs), the domain of outputs is also defined considering the errors in the north and east directions, i.e., ${\mathbf{\text{σ}}}_{\mathbf{E}}$ and ${\mathbf{\text{σ}}}_{\mathbf{N}}$ respectively. The output membership functions are shown in Figs. 5 and 6 respectively. The output membership functions are designed incorporating the possible changes in position in the East and North frame due to the effect of roll angle parameter. With reference to Fig. 3, the description of the output fuzzy variables is provided in Table 3. The corresponding set of fuzzy rules are provided in Table 4. The fuzzy rules for the remaining three regions are constructed in a similar manner. According to Fig. 3, the rules are constructed according to the following logic:

R1: IF the yaw angle is small and roll angle is zero THEN set E1 and E2 high and N1 and N2 low

With reference to the above logic, we define the input subset of three-degree division, ${\mathrm{S}}_{\mathrm{i}\mathrm{n}1}=\{\mathrm{L}\mathrm{o}\mathrm{w},\mathrm{M}\mathrm{e}\mathrm{d}\mathrm{i}\mathrm{u}\mathrm{m},\mathrm{H}\mathrm{i}\mathrm{g}\mathrm{h}\}$ for input-1 (yaw angle) and ${\mathrm{S}}_{\mathrm{i}\mathrm{n}2}=\{\mathrm{L}\mathrm{o}\mathrm{w},\mathrm{M}\mathrm{e}\mathrm{d}\mathrm{i}\mathrm{u}\mathrm{m},\mathrm{H}\mathrm{i}\mathrm{g}\mathrm{h}\}$ for input-2 (i.e., the roll angle corresponding to a change from negative to positive). The output subset is defined with a four-degree division given by ${\mathrm{S}}_{\mathrm{o}\mathrm{u}\mathrm{t}}=\{\mathrm{Z}\mathrm{e}\mathrm{r}\mathrm{o},\mathrm{L}\mathrm{o}\mathrm{w},\mathrm{M}\mathrm{e}\mathrm{d}\mathrm{i}\mathrm{u}\mathrm{m},\mathrm{H}\mathrm{i}\mathrm{g}\mathrm{h}\}$, corresponding to the minimum and maximum values for the rectangular boundaries, i.e., $\mathrm{E}1,\mathrm{E}2,\mathrm{N}1,\mathrm{N}2$.

Test Case-I

For the first test case the DEM having dimensions of 1,500 × 635 is shown in Fig. 7 where the bright colors in the contour plot represent higher altitudes. The corresponding three dimensional plot has also been shown in Fig. 7, where the brightest region (i.e., the yellow color) represents the highest elevation values. The region of complete aircraft maneuver is shown in Fig. 8 in red color in the region defined by latitude between 24.7 to 24.9 degrees and the longitudes between 66.9 to 67.4 degrees. In Fig. 8, the aircraft starts its flight initially from the $\mathrm{p}\mathrm{o}\mathrm{i}\mathrm{n}\mathrm{t}=[{24.88}^{\circ },{66.93}^{\circ }]$, continues its flight to the $\mathrm{p}\mathrm{o}\mathrm{i}\mathrm{n}\mathrm{t}=[{24.96}^{\circ },{66.97}^{\circ }]$ and returns back towards its initial point. The flight path between these position values are selected as the region of interest, i.e., the ROI which is highlighted by the ellipse in Fig. 8 and the magnified region is shown in Fig. 8 where the green color represents the region where the testing using the modified MAD and ESKF algorithm was performed. In Fig. 9, it was observed that in this region, the latitude and longitude of the aircraft were both increasing. The GPS measurements are used as reference for comparing the position obtained from the proposed method. In Fig. 9, the latitude and longitude obtained using the INS are shown, in which the drifts can be observed in the INS position estimates, i.e., the negative drift error in the longitude position generated by the INS. This affects the position information of the overall trajectory followed by the aircraft. In Fig. 10, the ESKF prediction in the test region are shown, the reference positions (i.e., the latitude and longitude) are represented by blue color and the estimated positions are shown in red. In Fig. 11 the variation in the aircraft velocity and 3D acceleration components are shown respectively. Figure 11 shows that the aircraft velocity increases from the time t = 600–800 s with the maximum velocity being 6 $\mathrm{m}/\mathrm{s}$. The aircraft yaw angle and the roll angle variations for Test-Case-I are shown in Fig. 11 respectively. The yaw angle follows a decreasing trend from 0.62 radians to 0.52 radians, whereas the roll angle can be considered approximately 0. In the proposed FL based MAD scheme, this corresponds to the case where only one angle, i.e., the yaw angle of the aircraft is changing and the effect of roll angle is not being considered. The ESKF prediction error is shown in Fig. 12, where the error range in the final prediction are 0.02–0.096 degrees respectively showing the effectiveness of the system matrix design of the proposed ESKF method and the optimal tuning of the filter parameters.

Test Case II

In the second test case the DEM having dimensions of 1,501 × 1,801 is shown in Fig. 13, where the bright colors in the contour plot represent higher altitudes. These values are considered as the terrain topographical database values. The corresponding three dimensional plot is shown in Fig. 13, where the brightest region (i.e., the yellow color) represents the highest elevation. In contrast to the Test Region-I, the region for Test-II is composed of more terrain feature. The region of complete aircraft maneuver is shown in Fig. 14 respectively. The flight path is shown in red color defined by latitude between 33.573 degrees to 33.581 degrees and the longitudes between 72.785–72.835 degrees. The aircraft starts its flight initially from the point = [33.57, 72.79] and continues its flight to the point [33.58, 72.83]. The flight path between these position values are selected as the ROI as shown in Fig. 14 where the region defined in the green color represents the region where the testing using the modified MAD and ESKF algorithm was performed. In Fig. 15, the GPS-based positions for this test case are shown. It can be observed that in this region the latitude of the aircraft is decreasing while the longitude of the aircraft is increasing. For this case, the test data was selected from time t = 12,000–12,200 s. Similar to Test Case-I, the aircraft sensors (i.e., the radar altimeter and the barometer) were used to compare the terrain height feature from the second pre-stored DEM. In Fig. 15, the latitude and longitude obtained using the INS is shown, in which the drifts can be observed in the INS position estimates, i.e., compared with the GPS reference position shown in Fig. 15, the negative drift error in the longitude position is generated by the INS which affects the position information of the aircraft. In Fig. 16, the estimates of the position outputs, i.e., the latitude and longitude are shown obtained from the proposed scheme. It can be observed that the drifts are significantly reduced when compared with the INS-based estimations. The velocity and acceleration of the aircraft for this test case are shown in Fig. 17 respectively. The velocity range is from 14.0–14.5 m/s. The velocity remains almost constant for t = 12,000–12,140 s and increases after t = 12,140 s. The velocity of the aircraft remains constant i.e., V = 14.44 m/s seconds for the rest of the aircraft motion from the time t = 12,160–12,200 s. The range of the aircraft heading i.e., yaw during this time is 4.766–4.768 radians and the roll angle is within the range of 0.19–0.22 radians as shown in Fig. 17 respectively. In this test case, the yaw angle variation is such that it remains constant till t = 12,140 s, increases rapidly after time t = 12,140 s and then decreases. While the roll angle starts decreasing slowly and a rapid change can be observed after time t = 12,140. This corresponds to the case where the effects of variation in both the yaw and the roll angle of the aircraft are considered in the proposed FL based MAD scheme. The predicted errors from the ESKF strategy are shown in Fig. 18. The latitude errors are obtained in the range 0.055–0.5 degrees. While the longitude errors are between 0.098 and 0.1047 degrees depicting minimum error in the acceptable ranges. This highlights the effectiveness of the system matrix design of the proposed ESKF method and the optimal tuning of the filter parameters and the appropriate design of the fuzzy logic system.

The values of the terrain roughness parameters for the two test regions are shown in Table 8. Where it can be observed that Test Region-II has more features compared to Test Region-I. The values of the metrics MAE, MSE and RMSE computed for two test cases are shown in Table 9 (with conventional MAD) and Table 10 (with fuzzy heuristic MAD). From the Table 9, it can be observed that the error values with the conventional MAD are lower for Test-Case-I compared to the position errors for Test-Case-II despite higher Sigma-Z and Sigma-T values. The reason is as follows: although, in general, the more unique terrain features, the higher is the positioning accuracy. However, the accuracy also depends on other factors such DEM-size and profile length selection (Yousuf & Kadri, 2024). In this article, for MAD correlation, the profile length ‘l’ is selected to be l = 10 for both test cases. The greater the profile length, the greater is the accuracy. Also, we used different DEM sizes for the two test cases. Based on the DEM size, no exact idea can be made as the accuracy may increase or decrease with minimizing or maximizing the DEM size. In Table 10, the error metrics for the proposed hybrid ESKF and FH-MAD are given. It can be observed that proposed method performs well in terms of position accuracy for the two DEMs having different terrain features with the errors being bounded within an acceptable range. However, the accuracy for Test Region-II is less compared to Test Region-I for the reason described above. Another reason could be the selection of the terrain related parameter ‘L’. This parameter is depending on terrain characteristics. However, it is assumed to be a constant parameter in this article and has been selected empirically. The computational time with the proposed MAD algorithm and the conventional MAD algorithm for the two test cases for Region-I and Region-II is compared in Fig. 19.

Comparison between conventional MAD and the proposed FH-MAD algorithm

The mean absolute difference (MAD) algorithm is frequently used for terrain matching in the standard terrain contour matching (TERCOM) system. The traditional MAD algorithm calculates the absolute difference between the corresponding elevation values at each position by sliding a reference terrain map over the elevation data obtained by the on-board sensor (usually from an altimeter or radar). The location that it produces the lowest MAD value is regarded as the most likely current position of the vehicle. The mean of these absolute differences is computed for every possible position. Because it must assess every potential match within a particular search area, regardless of terrain features or matching suitability, this thorough comparison process is computationally demanding. Particularly in real-time or resource-constrained applications, the deterministic approach may result in longer processing times and lower efficiency due to its inability to adjust to dynamic or ambiguous terrain.

The DEM size used in Test Case-II is greater than the DEM size in Test-Case-I. The conventional MAD therefore requires more time searching the DEM for Test-Case II compared to Test-Case-I. Similar trend can be observed with the proposed FH-MAD algorithm for the two regions. However, it can be clearly observed from the graph that the computational time is significantly reduced with the proposed FH-MAD correlation algorithm satisfying one of the major requirements for real-time performance. Furthermore, it must be highlighted that the quality of the predicted estimates from a TAN algorithm mainly depends on four factors: (1) the INS (2) the accuracy of the DEM (3) the terrain features and (4) the correctness of the altimeter measurements. The resolution of the DEM used in this work is 1 arc s. The predicted error for the latitude is up to five decimal places (i.e., 1 m). Also, the predicted error for the longitude is also up to five decimal places (i.e., 1 m). The major reason for this error, despite the convergence trend of the proposed ESKF algorithm, is the non-availability of a high-quality DEM. Therefore, the future work is focused towards the task of incorporating a high-quality DEM with the proposed hybrid technique in order to test the levels of accuracies of the predicted target region.