Abstract

With a series of complex lunar exploration missions carried out, reliable autonomous navigation approaches for lunar probes are in demand. In this paper, aiming to obtain the position of the probes with high accuracy and computational efficiency, a novel passive positioning scheme for lunar probes is presented. We first designed a reference system including four beacons with known positions. Then, we proposed arctangent temporal interferometry (ATI) method for ambiguity resolution in displacement measurement. Moreover, the quad-beacon displacement search (QBDS) algorithm is proposed for 2-dimension (2D) position search, which enables accurate and computational-efficient 2D positioning for lunar probes. The simulation validations show the effectiveness of the proposed positioning scheme and algorithms.

1. Introduction

In recent years, with the development of aerospace science and technology, there are increasing space probes launched for lunar, planetary, and interplanetary research, including Chang’e-5, InSight, and Tianwen-1 [14]. The moon is the closest celestial body to the earth, and it is also the first stop for humans to carry out deep-space exploration activities. In the past 60 years, more than a hundred spacecrafts have been launched for lunar exploration [5]. The probe navigation system has always been an important part of a lunar exploration mission [6, 7]. At present, lunar probes mainly rely on ground stations to provide telemetry and navigation services. However, as the moon exploration missions expand from the near to the far side of the moon, the ground telemetry is no longer available. Under such circumstance, reliable autonomous navigation approaches must be applied to ensure the probe has the ability to survive and explore during the mission cycle [8].

A number of autonomous navigation schemes have been proposed for lunar exploration missions. The inertial navigation is a commonly used method of autonomous navigation [9, 10]. It does not require any external information. However, its cumulative error leads to poor long-term accuracy. The stellar navigation method is another autonomous navigation method which can be used for lunar probes [1113]. It obtains the absolute position of the probe by measuring the angular differences between the earth and the navigation stars from the probe. Whereas the positioning accuracy using this method is limited by the accuracy of the sensors, which is not sufficiently high. A third method is to use relay satellites for radio navigation, which measures the distance and Doppler frequency to obtain the relative position and velocity information of the probe [14]. But unlike global navigation satellite system (GNSS) which is available all day, the visible time of lunar orbiting satellites is usually short as the number of satellites is very limited. Therefore, stable and continuous positioning services for lunar probes using relay satellites are infeasible at present. GNSS is also applied for lunar probe positioning, but the geometry of the visible satellites deteriorates greatly as the receiver leaves the earth [15]. There are attempts to build a GNSS-like system orbiting moon; however, the orbit accuracy and cost limits the application of such systems [16].

As a commonly used ranging method, carrier phase measurement has been utilized in various fields. For GNSS applications, the carrier phase can be measured to obtain the pseudo-range from the satellite to the receiver [17, 18]. The pseudorange has a very high precision but with an ambiguity due to the whole number of cycles between satellite and receiver is unknown. By resolving the integer ambiguity, the pseudorange with high precision can be utilized in high-precision positioning. Similar method is utilized for micromotion measurements using Doppler radar [1921]. By measuring the phase of the signal reflected and modulated by the moving target, the displacement of the target referring to the radar sensor is obtained. As the range difference exceeds the wavelength of the carrier, the phase difference is greater than , resulting the phase ambiguity. To solve this problem, a number of phase unwrapping approaches have been proposed, including conventional arctangent demodulation (CAD), extended differentiate and cross-multiplying (DACM), and arcsine demodulation [2023]. However, the disadvantages of these approaches include (1) high computational complexity and (2) significant accumulation error for fast-moving targets after approximation. These disadvantages limit the application of real-time high-precision ranging using carrier phase.

To solve the aforementioned problems, in this paper, a novel passive positioning scheme for lunar probes is proposed and investigated. A reference system including four beacons with known positions is designed which transmits continuous waves in a series of subbands from each beacon. Then, an accurate and computational-efficient algorithm for ambiguity resolution, which is called arctangent temporal interferometry (ATI), is presented for displacement measurement. After that, a target function for 2-dimension (2D) positioning based on the displacements referring to the beacons is constructed. A parallel 2D position search algorithm for lunar probes called quad-beacon displacement search (QBDS) was proposed, which enables accurate and computational-efficient 2D positioning results for lunar probes. The simulation results validate the effectiveness of the proposed scheme.

The remainder of this paper is organized as follows. The principle of displacement measurements using carrier phase and ATI algorithm is introduced in Section 2. The system configuration and proposed QBDS positioning algorithm of lunar probes is explained in Section 3. Section 4 gives the simulation results and analyzes the positioning performance. The conclusion is drawn in Section 5.

2. Displacement Measurement Using Carrier Phase

2.1. Signal Model

Consider a beacon transmitting a continuous wave (CW), neglecting the amplitude vibration, we can normalize the amplitude of the transmit signal to 1. Hence, the beacon transmits the signal with the following form, where and denote the carrier frequency of the signal and time, respectively. And denotes the arbitrary initial phase of the transmitter.

After transmission in the free space, the signal is received by the receiver on the lunar probe. The transmission induces a delay with respect to the transmit signal, that is, where is the amplitude of the received signal, is the noise term, and delay is determined by where denotes the instantaneous range from the beacon to the probe, and is the speed of light.

Then, the received signal is quadrature demodulated to baseband. Combining Equation (2) and Equation (3), the in-phase () and quadrature () outputs are represented by where is the arbitrary phase of the receiver, is the carrier wavelength, and denote the noise terms in / channels. We can observe from Equation (4) that the range can be resolved by extracting the phase of the received signal , that is,

Theoretically, we can obtain the range from the beacon to the probe with high accuracy via carrier phase measurement. And with these range information, the precise position of the receiver can be resolved as long as there are sufficient beacons. In this case, only two beacons are required at least for 2D positioning, and three beacons are required at least for 3D positioning (regardless of finite multiple solutions). Therefore, by placing several reference beacons with known exact coordinates on the moon face, the position of the lunar probes can be resolved by measuring the ranges from the beacons to the probe.

However, there are two practical problems for direct range measurement. First, the arbitrary initial phase of the transmitter and receiver is difficult to obtain. Meanwhile, the moving area of the probe is usually far larger than the scale of the microwave wavelength, which leads to phase ambiguity. In this scenario, the unambiguous range cannot be obtained directly.

In order to solve this problem, we apply range difference measurement instead of direct range measurement. Although we do not have the actual range, we can still resolve the variation of the range by continuously calculating the instantaneous phase of the received signal. The wrapped phase can be calculated by taking the arctangent of the received signal, that is,

To obtain unambiguous phase without jumping, unwrapping procedure is required, which is illustrated in Figure 1. A conventional unwrapping method is to perform direct phase shift, which is called conventional arctangent demodulation. Whenever the difference between consecutive phases is greater than or equal to , the following phases are shifted by adding multiples of until the difference is smaller than . Taking the first sampled phase as a reference, we have the unwrapped phase

Therefore, we can calculate the pseudorange from the beacon to the probe.

It is called pseudorange as the measurement has an unknown offset. However, this offset is not important as it is subtracted when calculating the range difference. Therefore, the range difference between the beacon and the probe from to can be determined by the following equation

The position of the probe at time and can then be resolved using the range differences from the probe to different beacons.

2.2. ATI Demodulation

As is shown above, the key to resolve range difference is to perform phase unwrapping. However, the aforementioned conventional arctangent demodulation algorithm has a series of drawbacks. First, the computational complexity is high as the entire data block from the current time to the end needs to be frequently shifted with period , which might be computational costly for the embedded system on the probe, as the computing resources are very limited. Also, this algorithm leads to serious accumulative error for a single wrong unwrapping point, where the latter data will have an entire phase offset.

In order to reduce the impact of wrong unwrapping, extended DACM algorithm is proposed. This algorithm calculates the phase increment between two adjacent samples and then accumulates the phase increments. Extended DACM algorithm greatly improves the performance facing the wrong unwrapping problem. The basic procedure of the extended DACM algorithm is as follows.

To obtain the phase increment, the continuous form is considered first. The angular frequency is the time derivative of the instantaneous phase, that is, where and denote the time derivative of and . The angular frequency can then be approximated in a discrete form, which is, where is the sampling interval. Here, we approximate the time derivative by the slope of adjacent samples. Then, we can derive the unambiguous phase increment , which has the expression

Taking , the unambiguous phase is represented as follows

where

In the extended DACM algorithm, the time derivative of and is approximated by the difference between two samples. Then, the phase increment is calculated by multiplying the instantaneous angular frequency and the time interval . These approximations may introduce additional error to the measurement, especially in low signal-to-noise ratio (SNR) scenarios.

In order to solve this problem, herein we propose an algorithm which calculates the phase increment between two adjacent samples regardless of the aforementioned assumptions, namely, ATI algorithm.

To simplify the phase calculation, we rewrite the baseband signal in complex form,

where is the imaginary unit. In discrete domain, it has the form

The phase increment can be resolved by taking the difference of and . Assuming that the sampling rate is high enough that the absolute value of phase increment is smaller than , the phase increment can be derived as the following equation:

Comparing Equation (16) with Equation (11), we can observe that the extended DACM does not consider the variation of signal amplitude over samples. If we assume , and the phase increment satisfies small angle assumption, i.e., when is sufficiently small, Equation (16) and Equation (11) can have a uniform form.

Then, similar to Equation (13), taking , the unambiguous phase can be represented by

where

As we already have the unambiguous phase, we can calculate the range difference according to Equation (9), that is,

3. 2D Positioning Scheme for Lunar Probes

3.1. System Configuration

As mentioned in Section II, the position of the probe at time and can be resolved using the range differences from the probe to different beacons as long as there are sufficient number of beacons.

Suppose the probe moves on a 2D plane, the trajectory starts from and ends at . As we can see, different from direct range positioning where only one coordinate needs to be resolved, there are two unknown coordinates, comprising four parameters . Therefore, if we can obtain at least four range differences, it is possible to determine the position of the probe.

In order to obtain a fine geometric dilution of precision (GDOP) and simplify the configuration, we set four beacons on the four vertices of a square, and the moving range of the probe is limited inside that square area. The positioning scene for lunar probe positioning is illustrated in Figure 2.

3.2. The Description of the Positioning Problem

Suppose the four beacons are distributed on the four vertices of a square with side-length , we can establish a Cartesian coordinate system referring to the beacons. The locations of the beacons are as follows: . The probe moves from the initial position to the terminal position , and the range difference from the initial position to the terminal position with respect to the four beacons is . We can locate the probe by solving the following equations

Equation (19) can be further written in a scalar form.

As we can see, by solving the four equations in Equation (20), the initial position and the terminal position can be estimated. However, Equation (20) is a set of nonlinear equations, which are difficult to solve the initial and terminal positions of the probe directly. In order to solve the initial and terminal positions of the probe, the following algorithm can be adapted.

Consider the probe is originally located at the initial position . For the four beacons , we can obtain the initial range from the probe to the beacons.

Then, we can then obtain the terminal range from the probe to the beacons, which can be calculated by where the range difference can be measured by applying ATI algorithm. Given two fixed reference points on the plane and the range from the third point to the two points, we can easily solve for the coordinates of the third point using law of cosines and we can divide the beacons into two separate groups to solve for the terminal position , respectively. For instance, Group A consists of beacons 1 and 2, and Group B consists of beacons 3 and 4.

Here, we derive the expression of terminal position calculated by beacons Group A. The geometry of the beacons and the probe is illustrated in Figure 3. The terminal angle of the probe with respect to Beacon 1 can be calculated by

Then, we derive the coordinates of terminal position . Theoretically, there are two solutions of the terminal position. The solution inside the feasible region is as follows.

Similarly, the coordinates of terminal position solved using beacons Group B are

If the initial position is correct, the terminal position calculated by two groups should be identical. Therefore, the problem is transformed into finding an initial position in the feasible region, so that the terminal positions and overlap. In other words, the distance between and can be used to measure the accuracy of the given initial position . Here, we define a target function , which is defined by the following equation

Combining Equations (21), (22), (24), (25), and (26), the cost function can be represented by a function of . As the specific form of the equation is very complicated, here, the cost function is simply denoted by

And the problem becomes an optimization problem, that is, where is the feasible region, which is the area enclosed by four beacons.

3.3. Newton-Raphson Method Based Algorithm

We can apply Newton-Raphson method to solve this optimization problem. Finding the minimum of the function requires finding the critical point which satisfies . After setting an initial solution , we can use the following formula for iteration. Each step generates a new iterate from the previous iterate according to where is the Hessian matrix defined by and is the gradient of defined by,

After each iteration, judge whether the difference between the two iterates is smaller than the threshold, i.e.,

If Equation (32) stands, output as the final solution of the initial position. Otherwise, continue the iteration. The whole process of initial position solution process is illustrated in Table 1.

3.4. QBDS Algorithm

Despite the initial position of the probe can be solved quickly by applying Newton-Raphson method, the calculation of the Hessian function is quite difficult. Also, Newton-Raphson method may not be able to obtain the global optimal solution of . Within the allowable range of accuracy, we can use a method based on grid search to solve the optimization problem in Equation (28), namely, QBDS.

As a simple but effective method, we can perform grid search to the initial position , and calculate the target function corresponding to each on-grid . By finding the minimum , we can solve the true initial position . Figure 4 gives an intuitive explanation of the QBDS process.

For convenience in presenting the figures, here, we define another target function . The problem becomes

Then, we can locate the probe in 2D plane when reaches its maximum.

We can observe from the QBDS algorithm that after calculating the range differences , the searching procedure can be run in a parallel form. This feature enables fast position resolution for real-time positioning, which helps the probe to survive in complex lunar environment.

4. Simulation Results and Analysis

In order to verify the feasibility and effectiveness of the theoretical part in this paper, a series of simulation is presented in this section.

4.1. Phase Unwrapping

To evaluate the phase unwrapping performance of ATI method, a scenario with linear range increment is considered. Suppose the lunar probe is originally located at the position of a beacon and moves away from the beacon with constant velocity  m/s. The carrier frequency is set  MHz. Three phase unwrapping algorithms including extended DACM, arcsine demodulation, and ATI are compared in the simulation. The curves of the mean square error (MSE) of the final position versus SNR are plotted in Figure 5, where the SNR is changing from 10 dB to 40 dB and each point is obtained utilizing 104 Monte-Carlo simulations.

As shown in Figure 5, the performance of ATI deteriorates for SNR lower than 12 dB. For SNR higher than 14 dB, ATI performs best among the three algorithms. All three algorithms approach similar SNR in phase unwrapping as the SNR approaches 40 dB. In this scenario, the algorithm for phase unwrapping has little effect on the result. The simulation validates the performance of the presented ATI algorithm for SNR higher than 14 dB.

4.2. Positioning

Suppose the lunar probe is restricted in a square area with length  m. We set four beacons on the four vertices of that area. Therefore, the locations of the beacons are, (m). In order to prevent signals from different beacons from interfering with each other, we adopt a frequency division multiplexing mode to manage the signals. We set the carrier frequency near 1000 MHz and a frequency interval of 1 MHz, so we can separate signals from four beacons with different band-pass filters. The carrier frequency of the beacons is set as follows:  MHz,  MHz,  MHz, and  MHz. As the Doppler frequency of the lunar probe is far smaller than the frequency interval, we can recover the received signal from each beacon without overlapping. The simulation parameters are shown in Table 2.

Three typical trajectories of the probe are investigated, those are (a) straight line; (b) semicircle; and (c) polyline, which are illustrated in Figure 6. The simulation time is 10 s, and the detailed parameters of the three trajectories are listed in Table 3.

To illustrate the influence of the SNR on range difference estimation and positioning, the SNR is set to 10 dB, 15 dB, and 20 dB, respectively. The range difference resolved by phase unwrapping and the distribution of the target function are presented in Figures 79. The positioning result of the terminal position indicated by the maximum point of the distribution of is shown in Table 4.

It can be observed that QBDS method works for all three trajectories. For , the large noise induces wrong unwrapping results, hence reduce positioning accuracy significantly. For and , no wrong unwrapping results, and the error of range difference is only induced by the phase noise within the period. As the SNR increases, the distribution map of becomes more concentrated, which may provide higher positioning accuracy for smaller search grids.

5. Conclusion

In conclusion, we proposed a positioning method and system for lunar probes. A novel phase unwrapping algorithm named ATI algorithm was presented, which calculates the unambiguous phase with higher precision and lower computational complexity. Then, a 2D positioning scheme for lunar probes is proposed, where we configure a set of square distributed beacons and apply QBDS algorithm for lunar probe positioning. Finally, the positioning method was evaluated by simulation. The feasibility and effectiveness of the method was validated by the simulation results.

Data Availability

All data, models, and code generated or used during the study appear in the submitted article.

Conflicts of Interest

The authors declare that they have no conflicts of interest.