Coherent Combining and Long Coherent Integration for BOC Signal Acquisition under Strong Interference

3.1 Enabling Technique 1: Long Coherent Integration in Stages

As shown in Figure 4, the sampling rate is f_s and each short-term correlation involves N samples over an interval of T₁ = N∆t_s with ∆t_s = 1/f_s. It can be shown that the m-th short-term correlation between t_m and t_m+1, denoted by c_m, can be written as:

1

where A is the signal amplitude and ϕ₀ is the phase difference between the incoming signal and carrier replica at t_m.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 4

Division of long integration (C_M) into multiple short-term correlations (c_m, m = 1, 2, …, M)

The short-term correlation c_m may be small if the signal is weak with a small A, or relatively small as compared to the interference signal not shown in (1). The correlation function R(∆τ_m) in Equation (1) is ideally of a triangular shape of unity height and suffers from an irreversible loss R(∆τ_m) = 1 − ∆τ_m/T_c < 1 due to the timing error |∆τ_m| ≠ 0. In acquisition, a tolerable error is typically chosen to be , where T_c is the code chip duration. This choice is determined by the code structure.

Similarly, the sinc-function in Equation (1) represents another irreversible loss due to frequency error ∆f. To limit the loss, the tolerable error is typically chosen such that . Clearly, it is determined by frequency search spacing and integration length, which our algorithm can control.

If the long integration is performed directly over T₂, the required frequency resolution is 1/(2T₂), inversely proportional to the integration interval, which significantly increases the number of search steps. In addition, the input data rate is very high at f_s. A more practical approach is to perform long coherent integration in two stages. In the first stage, short-term correlations would be performed for the input samples at high rate (f_s) over T₁ at coarse frequency search bins with spacing of 1/(2T₁). In the second stage, the short-term correlations, now at a lower data rate (1/T₁ instead of f_s), would be further added up to T₂ at refined frequency search bins with spacing of 1/(2T₂) as shown in Figure 5.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 5

Coarse frequency search (1st stage) and refined frequency search (2nd Stage)

3.2 Enabling Technique 2: Phase Rotation for Coherent Integration

The individual short-term correlations c_m in Equation (1) can be added up over M periods into the long coherent integration C_M so as to boost the weak amplitude A and make up the irreversible losses due to ∆τ and ∆f. C_M can be shown to have the following form:

2

Due to the frequency error ∆f, each short-term correlation is subject to an unknown phase rotation, which is the last term of Equation (2). We can obtain an estimate of this frequency error through frequency search. As shown in Figure 5, the refined frequency search bin around a coarse frequency search bin provides an estimate of the underlying frequency error as:

3

Then, the phase-rotated long coherent integration can be written as:

4

The phase rotation in Equation (4) can be implemented in hardware as a carrier numerically-controlled oscillator (NCO) for a given frequency . Alternatively, it can be implemented in software using the fast Fourier transform (FFT), which is equivalent to multiplying the short-term correlations with complex exponentials (phase rotation) for compensation at a refined frequency. The use of the FFT to perform the phase rotation in Equation (4), which is equivalent to Doppler frequency removal, is applied in double-FFT high-sensitivity acquisition (Jiménez et al., 2006), block-accumulating coherent integration over extended interval (BACIX; Yang & Han, 2006), and synthesized Doppler frequency hypothesis testing (SDHT; Kong, 2015). The FFT can be replaced with the chirp or fractional Fourier transform to deal with high dynamics (Wu et al., 2015).

3.3 Enabling Technique 3: Compensation for Code Migration

Driven by a non-zero code Doppler, the code phase of an incoming signal migrates from one search point to a neighboring one while the local code replica is fixed at the initial point. In the tracking mode, the code replica is adjusted constantly to catch up with the incoming signal. In the acquisition stage, however, a code replica is prepared per frequency bin with a potential maximum frequency error being a half of that bin. The code Doppler error is not an issue as long as the code migration stays within a grid point. For any integration interval over which code migrates across search points, less and less signal accumulation appears at the original search point but more and more on a neighboring one, which delimits an effective coherent integration for each given residual Doppler.

To exemplify the effect at L1 = 1,575.42 MHz, consider an example with T₁ = 1 ms for a carrier Doppler error of f_d = 250 Hz, the code Doppler error at a chipping rate of f_chip = 5.115 mcps is f_chip·f_d/L1 = f_d/308, which translates to a code migration rate of 0.8116 chips/sec. In other words, it takes about 0.2, 0.4, and 0.8 sec for the peak to migrate by ¼, ½, and one full chip. Residual Doppler causes the complete migration across a code phase search bin after 1.6 sec, and each code phase search point dominates about 0.8 sec as shown in Figure 6.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 6

Illustration of code migration and its effect on correlation peak values

In Figure 6(a), five consecutive code phases are color-coded and a signal starts at Code Phase 3 (red) and sweeps through Code Phase 4 (black) and then Code Phase 5 (pink). Code Phase 4 (the black curve) takes over Code Phase 3 (the red curve) at 0.65 sec, which is, in turn, taken over by Code Phase 5 (the pink curve) at 1.25 sec as shown in Figure 6(b). In this example, the long integration period is 2 sec and the refined Doppler search bin is 0.25 Hz. With a 5-ms increment being the data bit interval, the contribution to the cumulative peak of coherent integration per data bit is shown in Figure 6(c) for different code phases.

The general triangle shape in Figure 6(b) is the rising and falling of correlations when the code migrates through different phases. Note that on top of the triangle shape are stair-like variations with a stair width of 0.08 sec. There are eight stair steps over each ½ chips, which correspond to eight samples averaged prior to a simple 8:1 decimation in down-sampling (f_s = 16f_chip). A smooth response is expected if a low-pass filtering is used.

There are different methods to compensate for code Doppler as summarized below (Foucras et al., 2014):

• Method 1: Code Doppler compensation on spreading code sequences: One technique is to generate a local spreading code for each possible incoming Doppler frequency (Ziedan, 2006). The method is very expensive in terms of computation and storage. In another method, for every Doppler frequency estimate, an interpolation is made on both the received signal and the local replica to ensure that there is exactly one full spreading code period present in the coherent integration interval (Psiaki, 2001). For FFT-based correlation, the code delay induced by code Doppler can be compensated for in the frequency domain by multiplying the FFT of a local spreading code by a complex exponential corresponding to the delay (Akopian, 2001; Krasner, 1997).

• Method 2: Code Doppler compensation on correlator outputs: One technique is to choose the pickup points in coherent integration results to compensate for phase sliding based on determining the amount of code phase shift from the estimated code Doppler (Jiao et al., 2012). Another technique is the de-skewed short-term correlator (Betz, 1984), implemented in DirAc, which resamples/interpolates the short correlators for alignment. Migration of correlation peak, i.e., shifting the correlation function (Guo et al., 2017) can be compensated for in the frequency domain using the Keystone transformation (Li et al., 2006) in much the same way as range migration in delay-Doppler imaging of moving targets by radar.

In the present implementation, we use the approach of re-indexing short-term correlations based on the present code Doppler estimate. For a second stage frequency search bin around the coarse frequency search bin , the corresponding code Doppler (the code Doppler difference between the incoming signal f_d and the local replica driven by ) is given by:

5

where f_chip is the nominal chip rate.

Denote the code phase of the first chip in each short-term correlation along the long coherent integration interval by n_chip(m), m = 1, 2, …, M, which evolves from the initial code phase n_chip(1) driven by the code Doppler as:

6a

6b

In terms of code search bins, it can be written as:

6c

where b_chip is the search bin width in chips. For example, when the correlator spacing is a sampling interval, b_chip = f_chip/f_s = 5.115/56.32 ≈ 0.1 chips/bin where f_s = 56.32 MHz is the sampling rate. The short-term correlations at n_bin(m) for m = 1, 2, …, M are then selected for the search bin to perform the residual frequency compensation in coherent integration (4).

3.4 Enabling Technique 4: Compensation for Frequency Divergence between Sidebands

The upper and lower sideband signals, when treated separately as two BPSK signals, have a different center frequency. As such, they are subject to different Doppler frequency shift and frequency-dependent dispersive effects (Rushanan & Winters, 2010). In particular, when the interval for coherent integration is long enough to allow the code phase to migrate across a large portion of a code chip, the code migration introduces an amplitude modulation to the short-term correlations. The amplitude modulation can also shift the frequency of the upper and lower sidebands.

The spectra of the upper and lower sidebands are complex conjugates of each other. The above-mentioned frequency shifts are actually in opposite directions, leading to a frequency offset or frequency divergence between the upper and lower sidebands. The frequency offset becomes detrimental to signal combining when it is larger than the refined frequency resolution. Frequency alignment is therefore necessary to remove the frequency offset between the sidebands so as to reap more gain out of sideband combining.

As shown in Figure 7(a), the upper and lower sidebands of an IF signal are bandpass-filtered and down-converted to the baseband separately. Once in the baseband, they are processed independently and then combined. Figure 7(b) shows an example of combining upper and lower sidebands with frequency divergence and its compensation. In this example, T₁ = 1 ms (the short-term correlation interval), ∆f = 500 Hz (the frequency bin), ∆τ = 0.5 chips (the code phase spacing), f_d = 250 Hz (the maximum frequency error for T₁), T₂ = 5 ms (the bit interval), and T₃ = 1,000 ms (the long coherent integration interval with 0.5 Hz required for post-correlation frequency compensation).

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 7

Frequency divergence and compensation

The Doppler shifts experienced by the upper and lower sidebands, denoted by and , respectively, can be written as:

7a

7b

7c

where c is the speed of light, v is the range rate, f_L is the carrier frequency (L1 or L2), f_sc is the subcarrier frequency, and f_d is the Doppler shift experienced at the carrier (the band center).

From Equations (7a) and (7b), the difference in Doppler shifts (frequency divergence) experienced by the upper and lower sidebands is:

8

Note that the Doppler frequency shift experienced by each sideband is only a few Hz (e.g., 1.625 Hz for f_d = 250 Hz at L1), which is not visible when T_LCI < 308 ms because ∆f⁽²⁾ > 1.625 Hz. However, it is visible when the refined frequency search bin spacing is on the order of 1 Hz (T_LCI = 500 ms) and separable on the order of 0.5 Hz (T_LCI = 1,000 ms). Clearly, the frequency divergence in this example exceeds the frequency resolution and thus needs to be accounted for prior to combining.

Considering the frequency divergence at different sidebands, bringing Equations (7a) and (7b) into Equation (8) gives:

9a

9b

9c

An additional term appears in the compensation (4) for the upper sideband signal as:

10

Similar equations for the lower sideband signal can be written as:

11a

11b

In the above analysis, the group delay between the sidebands (due to signal propagation and/or receiver antenna/RF artifacts) is not considered. However, for a large sideband separation, the group delay can no longer be ignored and needs to be factored in for coherent combining. Hence, a search is added in the delay space similar to coherent combining of L1 and L2 as described in Section 3.6.

3.5 Enabling Technique 5: Estimation and Correction of Symbol Signs

Modern GNSS signals contain data and data-less (pilot) signals. In Galileo E5, for instance, data and pilot signals are separately modulated on I- and Q-components. In GPS, on the other hand, the time division data multiplexing (TDDM) scheme is used for M-code (Dafesh et al., 2002) in which the data bits and data-less bits are put on the odd and even chips separately.

It is straightforward to apply coherent integration to pilot signals, i.e., data-less even chips (Dafesh et al., 2002). A conservative approach can apply non-coherent integration to data-modulated signals (odd chips) since the signs of modulated data symbols are unknown. Yet, such an approach results in a significant squaring loss, which can be detrimental for weak signal and strong interference scenarios. To boost the signal strength, it is possible to apply sign-guessing/correction when going from bit interval T₂ to long coherent interval T₃ (Soloviev & Dickman, 2011, 2014; Soloviev et al., 2008, 2012; Yang & Han, 2006, 2007).

One way to guess and correct the bits is to find their signs such that the accumulated signal power is maximized:

12

where the superscripts o and e stand for the odd chips (data-modulated) and even chips (data-less) while the superscripts u and l stand for the upper and lower sidebands, respectively; ci_N is the value of chips integrated over N data bits, and is the guessed sign for the data bit b_i.

Simulation results indicate that the successful rate of sign guessing depends strongly on the signal energy per data bit. Maximization over multiple data bits may reduce the bit error rate. Yet, it cannot guarantee the uniqueness of each guessed bit sequence because of the effect of residual Doppler errors. At the low data rate of 50 bps, the approach seems useful but becomes less practical at the high data rate of 200 bps or for extremely weak signals.

For GNSS signals with high symbol rate data and secondary codes that can change their sign within each primary code period, advanced methods that are tailored to individual code inner structures can be employed (Borio, 2008; Leclère et al., 2017; Svaton & Vejražka, 2020).

3.6 Enabling Technique 6: Coherent Combining of Multi-Frequency Signals

A natural extension of the above approach is to combine signals at multiple frequencies to further improve the acquisition performance. Signals at individual frequencies are first accumulated independently using the combination of the coherent and non-coherent integration strategies previously described.

The accumulation results in the form of signal energy vs. code phase and Doppler shift are then combined. However, ionospheric delays can create differences between code phases at L1 and L2 in excess of the code search grid size. Hence, an adjustment is performed to avoid energy losses.

The delay search space is formed to include the worst-case ionospheric delay between L1 and L2 signals. For each search bin, L2 accumulation results are first adjusted via the re-indexing of code phase search bins. For example, to adjust for the negative delay of one code chip, the code bin index is incremented by two (for a half-chip spacing). The re-indexing of frequency bins is not required.

Then, L1 and adjusted L2 accumulated signal energies are added non-coherently. The adjustment bin that maximizes the overall signal energy is used to determine the acquired code phase and Doppler frequency. Coherent integration is possible but it requires an additional search to compensate for phase errors between L1 and L2, which is a trade-off between performance gains vs. extra computations. An example of coherent combining E5aI, E5aQ, E5bI, and E5bQ in real data can be found in Yang and Soloviev (2016).

4.1 Ability of Long Coherent Integration

The ability to perform long coherent integration over one second is demonstrated in the example shown in Figure 8. A BOC(10, 5) signal under an open-sky condition is generated by a Spirent GPS RF Simulator and down-converted to baseband and sampled by a transform-domain instrumentation GPS receiver (TRIGR) that uses a stress compensated (SC)-cut quartz OCXO at 10 MHz (Gunawardena, 2007). The true frequency of the signal is marked with a black arrow in Figure 8(a).

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 8

Example of long coherent integration over 1 second at three frequency search bins

A coarse frequency search grid is chosen to have three bins spaced by 250 Hz. Bin 1 (purple) is smaller while Bin 3 (red) is larger than Bin 2 (green) by 250 Hz. The short-term 1-ms correlations have an effective frequency coverage of ±500 Hz. As a result, all three frequency bins are responsive to the signal, but the strongest is from Bin 2, which is the closest to the true signal. The refined frequency searches from the three coarse bins all converge to the same estimate as shown in Figure 8(b).

However, code migration occurs for the three coarse bins as shown in Figure 8(c). The code search grid has a spacing about 1/10 of a chip. Since Bin 2 is closest to the true signal, the frequency error is small for code and carrier replicas generated at this frequency bin. There is practically no code migration and the correlation peak stays at the same code phase index (the green curve) when integrated longer and longer over time.

On the other hand, since the frequency of Bin 1 is smaller than the true signal, its frequency error is negative. As the signal is stretched, the code phase is delayed. The short-term correlation peak moves to a smaller and smaller index for subsequent correlation intervals (the purple curve). Similarly, the frequency of Bin 3 is larger than the true signal, its frequency error is positive. As the signal gets squeezed, the code phase is advanced. The short-term correlation peak moves to a larger and larger index for subsequent correlation intervals (the red curve).

The correlation peak-to-noise floor average power ratio (PAR) as a function of pre-detection integration (PDI) = 1 ms to 1,000 ms is shown in Figure 8(d) at the three coarse frequency bins. The results demonstrate long coherent integration over 1 sec. Note that the step jumps at PDI = 64, 256, and 512 in Figure 8(d) are caused by frequency resolution improvement. In the implementation, the phase rotation specified in Equation (4) is implemented with the FFT, zero-padded to the next power of two of its length. A higher resolution means a smaller loss due to reduced frequency errors.

Another example explaining code migration due to frequency error is shown in Figure 9. As shown, the three frequency bins all capture part of the signal energy. For the 1-ms correlations, their frequency response has the first nulls at ±1 kHz. Any two adjacent bins, which can be as large as 500 Hz apart, have a non-trivial gain to the incoming signal.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 9

An illustrating example of code migration

However, the correlation peak values for different frequency bins do not always coincide in code phase, but rather are separated more and more as the integration interval increases as shown by partial sums of correlation plotted against code phases at PDI = 100, 500, and 1,000 ms in Figure 9(a), 9(b), and 9(c), respectively.

Figure 9(d) illustrates the code migration process. The code Doppler rate is 0.811 chips/sec for a frequency error of 500 Hz. At those frequency bins, the code replicas have to start from larger offset code phases and then are shifted by their respective code Doppler so as to align with the incoming code phases.

Ideally, the signal strength increases linearly with the integration time. The linearity, however, may not hold up in reality due to implementation impairments. Figures 10(a) and 10(b) show the partial sums of coherently integrated odd and even chips for upper and lower sidebands on L1 and L2 for PRN7, respectively.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 10

Impairments affecting effective gain of long coherent integration

First note that the partial sums of the lower even (LE) are higher than the upper even (UE) in both L1 and L2, which are attributed to channel asymmetry. The RF front-end (RF FE) used in our data collection was not designed for the bandwidth of BOC(10, 5); furthermore, the bandwidth of an upper sideband extraction filter is narrower than that for a lower sideband. For L1 in Figure 10(a), the lower even (LE, indicated in orange) and lower odd (LO, indicated in yellow) have about the same magnitude whereas the upper odd (UO, indicated in purple) is slightly higher than the upper even (UE, indicated in blue).

It is rather clear in Figure 10(b) that the partial sums start to flatten out, which is due to uncompensated frequency error and clock instability that start to manifest in long integration. In addition to frequency drift and code migration, this saturation effect is another limiting factor as to how long the coherent integration can last. Due to saturation, the small increase in signal strength may not justify the additional amount of resources needed to extend the coherent integration further—a diminishing return aspect to be taken into account in the design trade-offs.

4.2 Analysis of Correlation Functions over Search Space

Figure 11 shows the correlation functions integrated over one second over a time-frequency search grid consisting of three coarse frequency bins (500-Hz spacing) expanded to 3,000 fine frequency bins (0.5-Hz spacing), and 169 code phases (1/10-chip spacing). It was under an open-sky condition for PRN = 7 on L1, which had the true Doppler f_d ≈ −477 Hz. Three coarse frequency search bins were at . The residual frequency error was 77 Hz from the second bin.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 11

1-s correlations over time-frequency search grid

The even chips on lower and upper sidebands were coherently integrated over 1 s as shown in Figure 11(a) and 11(b), respectively. For the odd chips on lower and upper sidebands, they were coherently integrated to 5 ms (a data bit duration) and then non-coherently to 1 s as shown in Figure 11(c) and 11(d), respectively. Additional examples for PRNs on L1 and L2 from different initial conditions can be found in Yang et al. (2020b).

Figure 12 compares the correlations of even chips, odd chips, and their combinations. In the code-phase domain, Figure 12(a) shows the correlation functions of even chips that were coherently integrated to 1 second for all frequency search bins. Figure 12(b) shows the correlation functions of odd chips bitwise coherently and then non-coherently integrated to 1 second.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 12

Coherent vs. coherent/non-coherent integrations of even and odd chips

As shown, both the even and odd chip correlations had about the same peak value, size, and shape. However, the noise floor of the coherently integrated even chips in Figure 12(a) proved to be much lower than that of odd chips in Figure 12(b). In contrast, Figure 12(c) shows the non-coherent integration of the even and odd chips from both L1 and L2, which had a peak four times higher.

In the frequency domain, Figure 12(d) shows the correlation functions of even chips coherently integrated to 1 second for all searched code phases. Figure 12(e) shows the correlation functions of odd chips bitwise coherently and then non-coherently integrated to 1 second. As shown, both the even and odd chip correlations had about the same peak value. However, the noise floor of the coherently integrated even chips in Figure 12(d) was much lower than that of odd in Figure 12(e).

Besides, the coherently integrated even chips in Figure 12(d) had a very narrow peak with the first nulls at ±0.5 Hz. In contrast, the non-coherently integrated odds had an elevated noise floor and its mainlobe was wide at ±200 Hz corresponding to the bit duration of 5 ms. It indirectly points to the benefits of coherent over non-coherent integration. Figure 12(f) shows the non-coherent integration of even and odd chips from both L1 and L2, which had a peak four times higher. It also shows secondary peaks in frequency, which is explained next.

4.3 Analysis of Secondary Peaks

Secondary peaks develop in the special case in which the true signal appears in the middle of two adjacent frequency bins. Consider PRN = 7 on L1 whose true Doppler is about f_d = −477 Hz and four coarse frequency search bins are chosen at . The true signal would be 250 Hz away from the middle two bins.

Figure 13(a) shows the non-coherent integration of the bitwise correlations of odd chips combined from upper and lower sidebands where a strong secondary peak appeared along the code phase. The color-coded correlation surface is shown in Figure 13(b), and the cuts along the frequency axis for two code phases and the cuts along the code phase axis for two frequency bins are shown in Figure 13(c) and 13(d), respectively. Note that the two USB curves are lower than the LSB ones due to channel asymmetry as explained in Figure 10.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 13

Bitwise correlations of odd chips on upper and lower sidebands: Secondary peaks when the signal straddles two frequency bins

The peak along the frequency in Figure 13(c) was made of two halves: The left half peaked at an index of 2,000 (−476.5 Hz) and had its left null at 1,600 (200-Hz offset) while the right half peaked at 2,001 (−477 Hz) and had its right null at 2,400 (also 200-Hz offset). The two halves together made a peak width of ±200 Hz, corresponding to the data bit of 5 ms, as expected. Note that only a search at a negative frequency (refined frequency bin index smaller than 2,000) could produce the correct alignment for code phase 135 and similarly only a positive frequency (larger than index 2,001) could produce the correct alignment for code phase 118, which explains the sharp drop around 2,000.

However, two peaks of comparable amplitude appeared along the code phase axis as shown in Figure 13(d). Namely, the 1-s non-coherent integration of bitwise coherent correlations of USB and LSB odd chips at Frequency Bin 2,001 provided the correct peak at code phase 118, while the non-coherent integration at Frequency Bin 2,000 produced a false peak at 135, which was off by 18 samples.

There are several reasons for this phenomenon. First, the frequency response of the bitwise coherent integration was so wide that it allowed many frequency bins to “compete,” so to speak. Second, when the true signal falls exactly in the middle of two adjacent coarse frequency search bins, it can be reached from the opposite direction, giving rise to spurious peaks. Third, a frequency search bin on the wrong side accumulates from a code phase also on the wrong side. For the example in Figure 13, the separation between the true and false peaks was twice the code phase drift over 1 s for a frequency error of 250 Hz, that is, 2f_chip/L1 × f_d in chips or 2f_s/L1 × f_d = 2 × 56.32/1575.42 = 17.87 samples. Further explanations and additional examples about spurious peaks can be found in Yang et al. (2020).

4.4 Embedded Acquisition Performance under Interference

The generation of test scenarios subject to interference is shown in Figure 14. Interference signals were injected digitally to the BOC(10, 5) data (an M-prime code) generated by a Spirent GPS RF Simulator. The RF signal was recorded by a special GPS receiver, TRIGR (Gunawardena, 2007), from which the noise power was estimated as was the signal power, determined with the estimated C/N₀ and verified against the value in the configuration file of Spirent GPS RF Simulator.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 14

Block diagram of interference signal injection

As shown, two types of interference signals were injected. One was a band-limited broadband (BLBB) noise jammer within ±10 MHz, which wipes out the band center where C/A-code and P(Y) code mainly reside so as to emulate a blue force electronic attack (Dominguez et al., 2016; Fishman & Betz, 2000; Holmes et al., 1999; Leclère et al., 2016; Ward, 2003b). A random noise was generated and passed through a low-pass filter and the samples were scaled so that the output noise power matched the desired jamming-to-signal ratio (J/S) of 10 dB. The other interference signal was a matched-spectrum jammer (MSJ) as another BOC(10, 5) with its strength adjustable up to J/S = 51 dB.

An automatic gain control was applied to the combined signal to avoid signal clipping and maintain it within the level for type conversion from double (64 bits) to integer (16 bits), which was then adaptively quantized to three bits for hardware correlation on a Xilinx Zynq (FPGA + ARM) card.

As shown in Figure 14, the 1-ms short-term correlations generated by FPGA-implemented parallel correlators as part of the first stage coarse frequency search are sent to the embedded ARM processor for coherent integration with the second stage refined frequency search. At the same time, the 1-ms short-term correlations are off-loaded for post-processing and verification.

Due to hardware limitations (e.g., the asymmetry of upper and lower sideband channels shown in Figure 10), only the even chips on the lower sideband were processed to illustrate the embedded SDR ability for long coherent integration. The performance of direct acquisition under interference is summarized in Figure 15. The correlation peak-to-noise floor peak power ratio (PPR) in dB (the y-axis) is shown as a function of PDI = 1 ms to 1,000 ms (the x-axis) for eight interference conditions with J/S = 10, 20, 30, 35, 40, 45, 48, and 51 dB, respectively, plus a noise-only case.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 15

Peak-to-peak ratio of CCLCI vs. PDI for different J/S

During each test run, the short-term correlations were dumped into a data file for offline verification. These short-term correlations were coherently integrated for different PDI, shown as the solid lines and color-coded for different J/S values in Figure 15. The embedded 1-s long coherent integration results output by the hardware are marked with crosses (×) in Figure 15. At PDI = 1,000 ms, the offline processing results matched well with those of embedded processing, demonstrating the CCLCI performance on hardware.

The time and frequency parameters obtained from the initial 1-s correlation, which was based on repeated replay of the signal samples stored in memory on the Xilinx Zynq card, were projected for live signals in subsequent 1-s correlations. The correlation peak was again clearly above the noise floor, indicating that the signal was indeed detected and the signal’s time and frequency estimates were valid. However, the second 1-s correlation peak was not as sharp as the first one. It is understandable because the long coherent integration is very sensitive to frequency errors (need to be better than 0.5 Hz for 1-s correlation) and over the integration intervals, the signal parameters had changed due to signal dynamics and/or combined satellite and receiver clock drifts.

4.5 Combined CI and NCI Implementation Trade-offs

The ability of CCLCI for 1-s coherent integration was demonstrated in embedded tests. It is of further interest to compare CCLCI with standard acquisition schemes based on short coherent integration (CI) cascaded with non-coherent integrations (NCI) as suggested by Pany (2020) to better understand the performance gain of CCLCI especially during radio frequency interference (RFI) and to provide insights for implementation trade-offs. Among the possible post-correlation or pre-detection integration schemes (Gómez-Casco et al., 2020), we consider here three different CI and NCI combinations with CCLCI as the baseline as:

• Scheme 0: CI = 1,000 ms and NCI = 0

• Scheme 1: CI = 5 ms (corresponding to high data rate of 200 Hz) and NCI = 200

• Scheme 2: CI = 20 ms (corresponding to low data rate of 50 Hz) and NCI = 50

• Scheme 3: CI = 100 ms and NCI = 10

The same embedded 1-ms correlations as Scheme 0 in Figure 15 were processed for Schemes 1, 2, and 3 and the results are shown in Figure 16(a), 16(b), and 16(c), respectively.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 16

PPR vs. PDI of combined coherent and non-coherent integration schemes under different J/S

The case without signal (noise alone) was used to establish the baseline with the PPR, which was used as a metric rather than the usual PAR for performance evaluation. It is a choice of convenience (easy implementation and fast calculation). Typically, the PPR is lower than the PAR, thus this choice was on the conservative side. As a measured quantity, PAR is an estimate of the underlying SNR and may be used interchangeably with SNR in analysis. A simplistic analysis of the relative performance is presented below.

First, look at the noise-only case where no signal is present. For Scheme 0, the output PPR of CCLCI stabilized at 10 dB. For the combined CI/NCI schemes, it went from 4 ~ 5 dB down to 1 ~ 2 dB and became smaller for shorter CI and, asymptotically, it went to 0 dB for larger NCIs. This can be understood as NCIs provide an estimate of the noise power (variance) for all search points and it thus equalizes the peaks and valleys in the long run.

For Scheme 0 with coherent integration, the PPR proved to grow up even though at a gradually slower rates and after 500 ms, the gain may not have warranted the expensive computations in view of the diminishing return effect. For the combined CI/NCI schemes, their growth rates were smaller and even more so for shorter CIs. It practically stopped at 200 ms for CI = 5 ms in Scheme 1, 400 ms for CI = 20 ms in Scheme 2, and 600 ms for CI = 100 ms in Scheme 3, respectively. This indicates a practical limit on how long NCIs can be useful, especially for weak signals.

If a 5-dB or more margin in a PPR was required in test threshold for reliable detection, Scheme 0 could stop at 500 ms for J/S = 51 dB. Scheme 1 with CI = 5 ms could not detect a signal with J/S > 51 dB, but those with < 45 dB after 250 ms. Similarly, Scheme 2 with CI = 20 ms could not detect a signal with J/S > 51 dB, but those with < 40 dB after 200 ms. Finally, Scheme 3 with CI = 100 ms could not detect J/S > 51 dB (barely after 800 ms), but < 48 dB after 500 ms. This clearly shows the superior detection performance of Scheme 0 over the CI/NCI schemes.

Consider now the output of the PPR of the three combined CI/NCI schemes under J/S = 10 dB. The end PPR at 1 second was 21 dB for Scheme 1, 25 dB for Scheme 2, and 31 dB for Scheme 3, respectively. From CI = 5 ms, 20 ms, and 100 ms to 1,000 ms, the time factor was 200, 50, and 10 and the corresponding coherent processing gain was 23, 17, and 10 dB, respectively. Their difference relative to Scheme 0 with PPR = 42 dB was about 21, 17, and 11 dB, respectively, which are well matched within a few dB.

As shown in Figure 16, the combined CI/NCI schemes developed a dip in their initial PPR values for the cases with shorter CIs and larger J/S. This dip also appeared in the noise-only case. The PPR of stronger signal cases picked up again, but the noise-only case did not regrow as some weaker signal cases didn’t. These dips are attributed to the effect of squaring loss (SL) during the initial short coherent integration in which noise and jamming dominated until the accumulated signal started to become meaningful (Betz & Cerruti, 2020; Esteves et al., 2016; Lowe, 1999; Rodríguez et al., 2004; Strassle et al., 2007; Woo, 2000; Yang et al., 2020a).

However, the PPR of some cases even after a larger NCI did not exceed their very initial level. In these cases, it was not the non-coherent integration of the signal but the averaging-out of noise and jamming that determined the detection performance. As such, a variable detection threshold could be set, which would decrease as a combined CI/NCI scheme progressed in order to locate the weak signals, quite counter-intuitively due to the nonlinear effect.

Roughly speaking, J/S = 51 dB is equivalent to C/N₀ = 20 dB-Hz within 10.23 MHz. Since the coherent integration only uses even chips, the effective C/N₀ was reduced by 3 dB to 17 dB-Hz. When the CI = 1,000 ms, the equivalent SNR in the 1-Hz bandwidth was 17 dB, which is the PPR of Scheme 0 at the end showed in Figure 15. For Scheme 1, the equivalent SNR in the 200-Hz bandwidth for CI = 5 ms was −6 dB, which was boosted to 5.5 dB after the NCI = 200. Due to the SL of 3 dB, the equivalent SNR was 2.5 dB while the experimental value was about 1 dB as shown in Figure 16(a).

For Scheme 2, the equivalent SNR in the 50-Hz bandwidth for CI = 20 ms was 0 dB, which was boosted to 6.5 dB after NCI = 20. After the SL of 3 dB, the equivalent SNR was about 3.5 dB, which is near the experimental value of the 3 dB shown in Figure 16(b). Finally, for Scheme 3, the equivalent SNR in the 10-Hz bandwidth for CI = 100 ms was 7 dB, which was boosted to 12 dB after the NCI = 10. Due to the squaring loss of 3 dB, the equivalent SNR was about 7 dB vs. the experimental 6.5 dB shown in Figure 16(c).

On the other hand, for low J/S cases in which the interference was comparable to the noise level, the output SNR in ratio can be written as (s/n)_out = P_s/(P_n + P_j) = (P_s/P_n)/(1 + P_j/P_n) ≈ (P_s/P_n)/2 = (s/n)_in/2 when P_j/P_n ≈ 1. In dB, (S/N)_out ≈ (S/N)_in − 3. For Scheme 0, (S/N)_in = 45 dB in the 1-Hz bandwidth after the 1-s coherent integration under an open-sky condition of C/N₀ = 45 dB-Hz before jamming insertion. As a result, (S/N)_out = 45 − 3 = 42 dB, which matches the top first curve well at J/S = 10 dB shown in Figure 15.

Similarly for Scheme 2, the SNR in the 50-Hz bandwidth for CI = 20 ms was 45 − 17 = 28 dB, which was boosted to 6.5 dB after NCI = 20 (no SL at high SNRs), leading to (S/N)_in = 34.5 dB, and (S/N)_out = (S/N)_in − 3 = 34.5 − 3 = 31.5 dB, which is close to the top first curve at J/S = 10 dB as shown in Figure 16(b). The analysis also holds for Schemes 1 and 3.

When the input J/S went from 10 dB to 51 dB, an increase of 41 dB (the first and eighth curves from the top in Figures 15 and 16), the output PPR for Scheme 0 went down from 42 dB to 17 dB, a reduction of 25 dB as shown in Figure 15. Similarly, a reduction of 9, 22, and 25 dB for Schemes 1, 2, and 3 was observed, respectively. There was no linear relationship between the change of the input J/S and that of the output PPR, but they can be accounted for by such factors as S/(N+J), CI, , and SL, as analyzed above.

In the initial demonstration, the test data with a smaller number of frequency bins and hundred code phases were considered for a ballpark analysis. Test runs may one day be extended to cover wider time and frequency search intervals, further taking into consideration such implementation losses as quantization errors on FPGA and arithmetic rounding-off errors on ARM, among others. Theoretical analysis is of interest using the system detection and false alarm probabilities (Borio et al., 2006a) as well as extreme value theory (Gumbel, 2004; Turunen, 2007) and the advanced concepts of effective C/N₀ (Borio et al., 2006b) and spectral separation coefficients (SSC; Betz, 2000, 2001).

4.6 Receiver Operating Characteristics (ROC)

The receiver operating characteristics (ROC) curve is defined as the probability of detection (P_d) as a function of the probability of false alarm (P_fa). In our case of signal acquisition, P_d is therefore the probability of detecting a desired signal under the H₁ hypothesis that the signal is present. P_fa is the probability of declaring the presence of a signal under the hypothesis H₀ in which there are only noise and interference without the signal.

A total of 400 runs, 200 for H₁ and 200 for H₀, were made in the Monte Carlo simulation to estimate the ROC curves for the CCLCI method. A single sideband signal of BOC(10, 5) was subject to a white noise with SNR = −10 dB and an MSJ of J/S = 53 dB. The signal was sampled at f_s = 4f_c, f_c = 20.46 Msps.

The initial timing uncertainty was 10 chips, which was covered by 40 correlators and spaced by one sample for the code phase search. The initial frequency uncertainty was ±1 kHz, which was covered by five frequency search bins with the coarse frequency search spacing of 500 Hz for the short correlation interval of T₁ = 1 ms. The actual frequency was at f_d = 249.75 Hz, which is close to the largest frequency error possible of 250 Hz between two frequency search bins.

In the simulation, PDI varied from 600 ms to 1,000 ms with an increment of 100 ms as shown in Figure 17. The corresponding refined frequency search spacing was chosen to be 1/(2 × PDI). At PDI = 1,000 ms, the spacing was 0.5 Hz, and the residual frequency error for f_d = 249.75 Hz was −0. 25 Hz relative to the coarse bin at 500 Hz, which represents a worst-case scenario with the largest possible frequency error.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 17

Estimated ROC curves for different PDIs under SNR = −10 dB and J/S = 53 dB from Monte Carlo simulations

Two test statistics were generated at the search grid points for detection: one was the correlation peak-to-noise floor average power ratio (PAR) and the other was the peak-to-noise floor peak power ratio (PPR), both excluding the mainlobe (±2 and ±4 in frequency and time directions, respectively) around the point under test. The comparison of the two plots in Figure 17, though simplistic, indicates that PAR outperforms PPR in terms of ROC curves, which is consistent with the observations made in Geiger et al. (2012).

In our embedded implementation shown in Figure 15, keeping the two largest peaks when stepping through frequency bins and code phases proved to be faster and easier than obtaining the average across the entire time and frequency search zone, thus is the logic behind the choice of PPR over PAR, a trade-off of performance for complexity.

4.7 Clock Stability for Handover from Acquisition to Tracking

Local oscillator stability (or lack of it) remains an important performance-limiting factor in long coherent integration. The phase noise (short-term stability) effectively has two effects on correlation, namely, the expected correlation peak decreases and the peak variance increases with noise (Serna et al., 2010). Due to the de-coherence of the carrier, the signal energy accumulation no longer increases linearly with time, limiting the ultimate sensitivity level one can expect from long coherent integration (Watson et al., 2006). As an example, the saturation of partial sums as shown in Figure 10(b) is attributed to the combined effect of uncompensated frequency error and clock instability that manifests in long integration.

In Dominguez et al. (2016), three types of clocks, namely, TCXO (h₀ = 9.43 × 10^–20, h_–2 = 3.8 × 10^–21), CSAC (h₀ = 7.2 × 10^–21, h_–2 = 2.7 × 10^–27), and OCXO (h₀ = 3.4 × 10^–22, h_–2 = 1.3 × 10^–24), were analyzed for their effects on coherent integration. TCXO could be used to integrate up to 100 ms coherently but not satisfactorily when the coherent integration time is 500 ms and beyond. CSAC could be used for coherent integration up to 500 ms, but shows considerable degradation at 1,000 ms. In contrast, OCXO could integrate coherently up to 1,000 ms with inconsequential effects from the phase noise.

For the embedded implementation described in the previous section, the samples of the incoming signal over 1 second were stored in a double data rate (DDR) memory and repeatedly processed faster than real-time for searching through the time-frequency uncertainty zone, which would have otherwise taken hundreds of seconds to find the first satellite. Once detected, an estimate of the elapsed time from the start of data collection was then propagated forward to the live signal. Clearly, a concern arises about the long-term stability on the order of several hundred seconds; that is, the oscillator drift over the search period may cause time deviation (TDEV) of an appreciable size with respect to the correlator acquisition window.

If the return from playback to live samples misses the signal (i.e., cannot find it where it was predicted), one has to search around that spot anew, albeit with a much-reduced uncertainty. TDEV can thus be expressed in terms of the modified Allan deviation as TDEV ≈ τσ_y(τ) where τ is the observation time. As an example, a Bliley LP62 OXCO has a TDEV on the order of 10^–9 s for τ ~ 100 s. With 1,024 correlators spaced at 0.454 chips at 5.115 MHz, the acquisition window amounts to 9.1 × 10^–5 s. For this case, the 3-sigma timing propagation error remains well within the acquisition window. In general, the goal of the system design (the tolerance that can be accommodated) was to ensure that the three-sigma predicted interval would be well within the parallel correlators so that the code and carrier tracking loops could be closed once a signal is acquired.