In the simulation experiment, to make use of the measured channel characteristic data, additional channel characteristics can be disturbed according to the method which is shown in Figure 2 [38].

In Figure 2, S i n ( t ) is the received if signal. S i n ( f ) is the output transformed to the frequency domain, S o u t ( f ) is the output after an equivalent nonlinear phase low-pass filter, and S o u t ( t ) is the output transformed to the time domain by inverse Fourier transform. The output expression of the signal after passing through the low-pass filter is: S o u t ( f ) = S i n ( f ) H ( f ) (1)

In engineering applications, the amplitude-phase characteristics of antennas and RF channels can be measured by signal sources and spectroradiometers. Figure 3 shows the measured amplitude-phase characteristics of four single antenna channels measured by signal source and spectrum analyzer. It can be seen from Figure 3A that the magnitude response of the actual channel is non-flat, while it can be seen from Figure 3B that the phase response of the actual channel is nonlinear. The amplitude-phase characteristics of the channel in the simulation experiment will be simulated by using the actual measured data.

The time-domain adaptive anti-jamming technique based on interference estimation is to make use of the difference in predictability between narrowband interference and wideband signal, so as to form the estimated value of narrowband interference, and then cancel it from the received signal to eliminate the narrowband interference.

According to the literature [31], the correlation function after the signal passes through the non-ideal channel is as follows: R ( τ ) = ∫ − b b s i n c 2 ( f ) ⋅ A ( f ) ⋅ exp { j 2 π f [ τ − τ g ( f ) ] } d f (2) where τ g ( f ) is the group delay of the equivalent low-pass filter.

It can be seen from Eq. 2 that when the navigation signal passes through the non-ideal channel, if the ideal constant amplitude and linear phase filter cannot be guaranteed, the correlation peak will no longer be symmetrical, resulting in distortion and thus affecting the ranging performance.

Figure 4 shows the correlation function when navigation signals pass through ideal channel and non-ideal channel respectively. But from the Figure 4, we cannot intuitively distinguish the difference between the two correlation peaks. Therefore, in the satellite navigation signal channel indicator system, S-curve Bias (SCB), an important indicator, is usually used to measure the symmetry of correlation peaks and quantitatively analyze the ranging bias of receivers [39, 40]. The definition of SCB is expressed as follows, where τ 0   、 τ 1   a n d   τ 2 meet the following constraints: τ S C B = τ 0 − | τ 1 − τ 2 | 2 , { τ 0 = a r g m a x R ( τ ) R ( τ 1 ) = R ( τ 2 ) | τ 1 − τ 2 | = Δ (3)

In Eq. 3, R ( τ ) is the correlation function between the received signal and the local signal, and Δ is the interval of the correlator.

As shown in Figure 5, under ideal and non-ideal channel characteristics, the maximum pseudocode ranging deviation exceeds 0.5 ns? Once again it is proved that the error of ranging is caused by the navigation signal passing through the non-ideal channel.

In the single-antenna-based interference suppression algorithm, we adopt the time-domain interference suppression based on linear estimation with mature technical basic theory. In this paper, combined with the current situation of engineering implementation, the classical minimum mean square error (LMS) algorithm [41, 42] is selected to provide convenience for the following analysis and application.

In addition, the filter structure adopted in this paper is the common odd-order bilateral tapped transverse filter, namely the filter with interpolation structure. It uses both past and future data to estimate x ( n ) based on the data set { x ( 0 ) ⋯ x ( n − 1 ) x ( n + 1 ) ⋯ x ( N − 1 ) } .

The order of the filter is (2M + 1). In addition, to prevent the filter coefficient from converging to 0 and facilitate calculation, the intermediate weight of the filter is fixed as 1. In the subsequent studies of this paper, the odd-order bilateral tapped lateral filter is taken as an example for analysis. For an even-order unilateral tapped filter, the local pseudocode signal can be obtained by appropriate delay.

The flow chart of the robust ranging algorithm is shown in Figure 8.

First of all, in order to eliminate the superposition effect of related functions brought by the anti-jamming module in the iterative process, we will design the filter coefficient for preprocessing the local pseudocode according to the anti-jamming filter coefficient, and the specific calculation formula will be given later. The signal output y s through the preprocessing filter is expressed as follows: y s ( n ) = α H s = ∑ i = − L L α i s ( n + i ) ( L ≤ M ) (11)

Then, the navigation signal passing through the anti-jamming filter is correlated with the preprocessed local signal, and its correlation function is expressed as follows: R ′ ( τ ) = ∑ n = − M M y ( n ) y s * ( τ − n ) (12)

Substituting Eqs 7, 11 into Eq. 12 can be calculated as follows: R ′ ( τ ) = ∑ n = 1 N [ ∑ k = − M M w k x ( n + k ) ] [ ∑ i = − L L α i s ( τ − n − i ) ] * (13)

Further calculation of the above related functions can be obtained as follows: R ′ ( τ ) = ( ∑ i = 2 M − L M w i α i − M − L ) R x s ( τ + ( M + L ) τ 0 ) + ( ∑ i = 2 M − L − 1 M w i α i − M − L + 1 ) R x s ( τ + ( M + L − 1 ) τ 0 ) + ⋅ ⋅ ⋅ + ( ∑ i = M − 2 L + 1 M w i α i + 1 ) R x s ( τ + τ 0 ) + ( ∑ i = − L L w i α i ) R x s ( τ ) + ( ∑ i = − M 2 L − M − 1 w i α i − 1 ) R x s ( τ − τ 0 ) + ⋅ ⋅ ⋅ + ( ∑ i = − M L − 2 M + 1 w i α i + M + L − 1 ) R x s ( τ − ( M + L − 1 ) τ 0 ) + ( ∑ i = − M L − 2 M w i α i + M + L ) R x s ( τ − ( M + L ) τ 0 ) (14)

Further considering that when the order of the filter preprocessing the local pseudocode is equal to the order of the anti-jamming filter (that is, when 2L+1 = 2M + 1), the Eq. 14 can be further simplified to obtain: R ′ ( τ ) = ∑ j = M − 2 L + 1 , k = 1 k = M + L j = 2 M − L , ( ∑ i = j M w i α i − k ) R x s ( τ + k τ 0 ) + ( ∑ i = − L L w i α i ) R x s ( τ ) + ∑ j = L − 2 M , k = − M − L k = − 1 j = 2 L − M − 1 , ( ∑ i = − M j w i α i + k ) R x s ( τ − k τ 0 ) (15) When the order of the local pseudocode preprocessing filter is equal to the order of the anti-jamming filter (that is, when 2L+1<2M + 1), the Eq. 14 can be further simplified to obtain: R ′ ( τ ) = ∑ j = M − 2 L + 1 , k = M k = M + L , j = M ( ∑ i = j M w i α i − k ) R x s ( τ + k τ 0 ) + ∑ k = 1 M − 1 ( ∑ i = k − L k + L w i α i − k ) R x s ( τ + k τ 0 ) + ( ∑ i = − L L w i α i ) R x s ( τ ) + ∑ k = 1 − M − 1 ( ∑ i = − k − L L − k w i α i + k ) R x s ( τ − k   τ 0 ) + ∑ j = − M , k = − M − L k = − M , j = 2 L − M − 1 ( ∑ i = − M j w i α i + k ) R x s ( τ − k τ 0 ) (16)

In Section 3 of this paper, it has been deduced and proved that the aggravation of correlation peak distortion is caused by the superposition effect of correlation functions at different moments after the signal passes through the anti-jamming module. Therefore, in Eq. 15 and Eq. 16, we preprocessed the local signal through the filter coefficient α i , so that the coefficient before R x s ( τ ) of the principal component of the correlation function is 1, and the pre-coefficient of R x s ( τ ± k τ 0 ) at other different delay times is 0. The equations are as follows: { ∑ i = − L L w i α i = 1 ,   k = 0 ∑ i w i α i ± k = 0 ,   k ≠ 0 (17)

In the form of matrix, the solution to α i can be expressed as: [ w − 2 L w − 2 L + 1 ⋯ w − 1 w 0 w − 2 L + 1 w − 2 L + 2 ⋯ w 0 w 1 ⋮ ⋮ ⋰ ⋮ ⋮ w − 1 w 0 ⋯ w 2 L − 2 w 2 L − 1 w 0 w 1 ⋯ w 2 L − 1 w 2 L ] [ α − L α − L + 1 ⋮ α L − 1 α L ] = [ 0 ⋮ 1 ⋮ 0 ] (18)

According to Eq. 16, it not only ensures that the principal component R x s ( τ ) of the correlation function R ′ ( τ ) will not change, but also eliminates the superposition effect of the correlation function caused by iteration due to the anti-jamming process, and ensures robust ranging of correlation peak in the capture and tracking stage. At this time, the peak value of the correlation function closest to the peak value of R x s ( τ ) should be two code chips, and the preprocessing filter order is 2L+1. The formula of L is as follows, where ⌊ ⋅ ⌋ stands for rounding down. L = ⌊ 2 F s F c ⌋ (19)

In the navigation receiver, Early minus Later (E-L) code phase discriminator is used for ranging [28, 29]. In Eq. 19, if the delay component of the correlation function only affects the principal component of the correlation function beyond the E-L code interval, the ranging result will still not be affected. Therefore, when the preprocessed local signals are correlated with the signals passing through the anti-jamming filter, the robust ranging performance can be maintained by eliminating other correlation function delay components within the E-L code interval of the correlation function principal component. The code interval is 2∆, then the optimal order of the preprocessing filter through which the local signal passes is 2 L 0 + 1 , where L 0 is expressed as follows: L 0 = ⌊ ( 1 + Δ ) F s F c ⌋ (20)

In Figure 9, when the sampling rate is 20 MHz and the code rate is 10.23MHz, according to Eq. 20, it can be seen that the optimal preprocessing filter order is five orders. It can be seen from the figure that in scenarios with different narrowband interference and correlator intervals, when the preprocessing filter order reaches the theoretically optimal order, the ranging value reaches a stable state, and the optimal order is effective.

Carrier-to-noise ratio (CNR) is defined as the ratio of signal power to the noise power in a unit bandwidth, and is an important indicator for evaluating signal quality in navigation reception [43]. In order to further derive the loss limit of the CNR of the robust ranging algorithm proposed in Section 4.1, the signal energy loss due to the superposition effect of the correlation function is analyzed from the energy perspective of the correlation domain. In the correlation domain, in the absence of interference and noise, the area of the correlation function in [ − T c , T c ] is the signal energy. After passing through the anti-jamming filter, the total area of the correlation function at other time delays caused by the iterative effect is the signal loss energy. Regardless of interference and noise, the energy of the navigation signal correlation function is expressed as follows: E = ∫ τ = − T c T c | w 0 R s s ( τ ) | 2 d τ (21) where, T c represents the code width.

The total energy loss of the pure navigation signal after passing through the anti-jamming filter can be deduced as: E L o s s = ∑ k = − M M k ≠ 0 ∫ τ = k τ 0 − T c k τ 0 + T c | w k R s s ( τ − k τ 0 ) | 2 d τ (22)

Then the power of signal loss can be expressed as: P L o s s = E L o s s / 2 T c E / 2 T c = ∑ k = − M k = M , k ≠ 0 | w k | 2 | w 0 | 2 (23)

Because of the complexity of the expression of the filter coefficients in the time domain, the coefficients are derived from the ideal band-stop filter in the frequency domain. Assume that the symbols of the input signal, weight value, and output signal in frequency domain anti-jamming are X ( k ) , H ( k ) , and Y ( k ) , the corresponding symbols in the time domain are x ( n ) , h ( n ) and y ( n ) . Since the frequency domain and the time domain satisfy the Fourier transform relationship, filtering in the frequency domain can be equivalent to filtering in the time domain. The frequency domain and time domain anti-jamming weights can be expressed as: Y ( k ) = X ( k ) H ( k ) (24) y ( n ) = x ( n ) ⊗ h ( n ) (25)

Focus on analyzing the frequency domain filter H ( k ) , the relationship between H ( k ) and h ( n ) is as follows: P L o s s = E L o s s / 2 T c E / 2 T c = ∑ k = − M k = M , k ≠ 0 | w k | 2 | w 0 | 2 (26)

Assuming zero is set on the frequency interval of [ k 1   k 2 ) and the remaining interval remains 1, then the filter coefficient in Eq. 26 can be further expressed as: h ( n ) = i f f t [ H ( k ) ] = 1 N ∑ k = 0 k 1 − 1 e j 2 π N k n + 1 N ∑ k 2 + 1 N − 1 e j 2 π N k n = 1 N ∑ k = 0 N − 1 e j 2 π N k n − 1 N ∑ k 1 + 1 k 2 e j 2 π N k n = { N − ( k 2 − k 1 ) N n = 0 1 N e j 2 π N k 2 n − e j 2 π N k 1 n 1 − e j 2 π N n n ≠ 0 (27)

In order to prevent the filter coefficient from converging to 0 and to facilitate calculation, the middle tap is usually constrained to one when using a bilateral tapped transverse filter. Then in the time domain, as in Eq. 27, the modulus ratio of the tapped coefficient at n ≠ 0 to the middle tapped coefficient is as follows: w c o m _ a b s = | 1 N e j 2 π N k 2 n − e j 2 π N k 1 n 1 − e j 2 π N n | (28)

Using Euler’s formula to expand and further simplify: w c o m _ a b s = 1 N [ cos ( 2 k 2 π N ) − cos ( 2 k 1 π N ) ] 2 + [ sin ( 2 k 2 π N ) − sin ( 2 k 1 π N ) ] 2 [ 1 − cos ( 2 n π N ) ] 2 + sin 2 ( 2 n π N ) = 1 N sin 2 ( ( k 1 − k 2 ) π N ) sin 2 ( n π N ) ≤ 1 N 1 sin 2 ( n π N ) = 1 N 1 | sin ( n π N ) | (29) As shown in Eq. 29, w c o m _ a b s is regarded as a binary function of n and N. When N ∈ [ 2 , ∞ ) , n ∈ [ 1 , N − 1 ] , w c o m _ a b s gets the maximum value at N = 2 , n = 1 . By substituting the values of n and N into Eq. 31, the modulus ratio of tap coefficient at n ≠ 0 to the middle tap coefficient can satisfy the following relation: w c o m _ a b s = | w n ≠ 0 w 0 | ≤ 1 2 (30)

The signal energy loss is maximum when the delay component of the correlation function is outside the interval [ − T c , T c ] (equivalent to when the delay correlation function does not overlap with the principal component of the correlation function). Under the condition of channel characteristics, taking the third-order bilateral tapped transverse filter as an example, the maximum energy loss of correlation function of pure navigation signal before and after passing through the anti-jamming filter is shown in Figure 10:

In Figure 10, the red correlation function at the top is the correlation function before the pure navigation signal passed through the anti-jamming module under the channel characteristics, and the area of the correlation function within the interval [ − T c , T c ] is the signal energy. Therefore, under the N-order filter coefficient, the maximum signal power loss caused by the superposition effect can be obtained as follows.

Substitute Eq. 29 into Eq. 23: L o s s C F S _ MAX = E L o s s / 2 T c E / 2 T c = ∑ n = 1 N − 1 [ 1 N 1 | sin ( n π N ) | ] 2 (31)

By numerical simulation of Eq. 31, it can be seen from Figure 11 that when N is greater than order 10, the CNR loss caused by signal superposition tends to a fixed value of 1.8 dB (one decimal place is reserved). It is concluded that the maximum signal power loss caused by the adaptive correction algorithm is 1.8 dB.

The parameters for this simulation are listed in Table 1. Where ISR represents the interference-to-signal ratio, and CNR represents the carrier-to-signal ratio, and F s represents sampling frequency, F c represents code rate. In addition, in the simulation, the anti-jamming filter adopts the lateral filter with a double tap. The interference types are single frequency interference and narrowband interference with different bandwidths.

It should be noted that in the following experimental scenarios, the channel characteristics adopt the measured amplitude and phase characteristics shown in Figure 3. In addition, it should be noted that in the following simulation scenarios, the ranging deviation values of the vertical coordinates are based on the ranging values under the condition of no interference.

In Figure 12, according to formula (20), it can be known that the optimal preprocessing filter order is fifth, and from Figure 14, the distorted correlation peaks after passing through the traditional LMS algorithm are all repaired to different degrees after passing through preprocessing filters of different orders. When the preprocessing filter order is greater than or equal to the optimal order, the different delay components of the correlation function cumulative in the correlator interval can be eliminated, so that the distorted correlation peak can be adaptively corrected.

In order to better evaluate the ranging deviation at different correlator intervals, the SCB curves of the improved algorithm and the traditional LMS algorithm under different narrow-band interference are drawn. As shown in Figure 13, the blue curve is the deviation caused by the actual channel characteristics without interference. Compared with the traditional LMS algorithm, the proposed algorithm significantly reduces the ranging deviation, and the ranging value of the signal after passing through the non-ideal channel characteristics is also controlled within the index range.

In Figure 14, compared with the ranging value of the received signal after passing through the non-ideal channel characteristic, the ranging deviation of the traditional LMS algorithm under different narrowband interference increases with the increase of the interference bandwidth, and it can reach about 1.5ns under 2 MHz narrowband interference. The improved algorithm proposed in this paper corrects the deviation to within 0.5 ns, which conforms to the ranging performance index of the ranging receiver and achieves the effect of robust ranging.

In Figure 15, whether under the traditional LMS algorithm or the improved algorithm proposed in this paper, with the increase of interference bandwidth, the loss of the SNR also increases. The SNR loss of the improved algorithm is slightly higher than that of the traditional LMS algorithm, and its loss is less than 1.8 dB in narrowband interference scenarios with different bandwidths, which is consistent with the theoretical analysis in section 4.3.