Abstract
This research investigates the mutual interference between a GNSS signal refracted through the atmosphere and the same signal reflected from the surface below. For low-altitude occultations, the reflected signal can be present as multi-path when tracking the direct signal. This poses a problem when making remote sensing measurements using GNSS-RO or GNSS-R techniques. Here, this issue is addressed in the context of GNSS data collected by a mountaintop receiver on the Hawaiian island of Maui. First, an elevation-angle threshold, below which the reflected signal can be present as multipath on the direct signal, is calculated for the mountaintop geometry. Next, the first steps are taken to mitigate the mutual interference between the direct and reflected signals. A frequency-domain filtering approach and a time-domain smoothing approach are applied to the signal intensity in an effort to separate the multipath contribution and the direct signal contribution.
1 INTRODUCTION
As a signal propagates from the GNSS satellite to the receiver, it interacts with its environment—for example, the signal is refracted through the atmosphere and scattered from Earth’s surface below. Though these interactions are a nuisance when using GNSS for its originally intended purpose of navigation, they can be leveraged for remote sensing (Jin & Komjathy, 2010). GNSS radio occultation (GNSS-RO) obtains atmospheric parameters by measuring the amount of bending in the refracted signal (Kursinski, Hajj, Schofield, Linfield, & Hardy, 1997). GNSS reflectometry (GNSS-R) uses the reflected signal for purposes of scatterometry (e.g., to retrieve ocean surface winds or land vegetation characteristics) and altimetry (Zavorotny, Gleason, Cardellach, & Camps, 2014). One scenario of particular interest to the remote sensing of Earth’s atmosphere is the low-altitude occultation, enabling measurements of the moist lower troposphere, which has a high impact on weather patterns. In this case, the reflected signal will have a low-grazing angle (the angle between the local surface and the incident ray of the reflection). Low-grazing angle reflections are well-suited for carrier phase altimetry (Semmling et al., 2016). Furthermore, reflections at low-grazing angles occur significantly more frequently than those at high-grazing angles (Southwell & Dempster, 2018). Clearly, the ability to make GNSS-RO and GNSS-R measurements from low-altitude occultations is of importance.
However, these low-altitude occultations pose a technical challenge. In this scenario, it is possible for the path length difference between the direct and reflected signals to be small enough so that the reflected signal is present as multipath when tracking the direct signal. Thus, the possibility of making accurate remote sensing measurements from low-altitude occultations hinges on the ability to predict and mitigate the mutual interference between the direct and reflected signals.
The research described here addresses these two needs through the analysis of data recorded by an antenna dish placed on the summit of a mountain and directed towards GPS satellites rising over the horizon. The first goal is to calculate an elevation-angle threshold elth below which the reflected signal is present as multipath. To calculate elth, the code structure of the signals is considered, and the path length difference between the direct and reflected signals is calculated considering the altitude of the mountain-top receiver. The second goal is to devise a technique for separating the direct and reflected signals when operating below elth. Here, a frequency-domain filtering approach and a time-domain smoothing approach are applied to the tracked signal intensity (SI) in an effort to separate the effects of the direct and reflected signals. The smoothing approach, which outperforms the filtering approach at lower elevation angles, is further analyzed through a simulation study.
2 DATA COLLECTION
The experiment, illustrated in Figure 1, was performed in May 2017 near the summit of Haleakala on the Hawaiian island of Maui (latitude: 20°42’09’’ N, longitude: 156°15’24’’ W, altitude: 3060 m) (Morton et al., 2017). A 1.9-meter antenna dish (27 dB gain at GPS L1) was steered to track GNSS satellites near the horizon. In addition to signals from other GNSS constellations, three frequency bands of GPS signals were collected: L1 and L2 sampled at 5 MHz, and L5 sampled at 25 MHz. For low-elevation satellites, both the direct and reflected signals are within the 7-degree half power beam width of the antenna.
The signals are reflected from the dish and received by a horn antenna. The horizontal (H) and vertical (V) polarization outputs of the horn antenna are converted to RHCP and LHCP using 90-degree hybrid combiners, and raw data are recorded using six separate USRPs
USRP N200 software-defined radios are used to record raw data for the three frequency bands on RHCP and LHCP polarization channels. Although GPS signals typically flip polarization from RHCP to LHCP upon reflection from Earth’s surface, reflections for extremely low-grazing angles (smaller than 3 degrees) retain a significant amount of power within RHCP. Thus, the geometry of this experiment allows the direct and reflected signals to be recorded on the same polarization channel. Due to the reflection from the dish into a horn antenna, both the direct and reflected signals flip polarization from RHCP to LHCP just before reception. The data used in this paper come from tracking the data recorded on the LHCP channel (USRP 4, 5, and 6). The signals are tracked using a 0.5 chip correlator spacing and a 100 Hz time series of the SI for each frequency is computed in post processing from the I and Q correlator outputs. The data are divided into two datasets: dataset 1 contains GPS L1, L2, and L5 SI from a single rising occultation event (PRN 3 from approximately 8:00 to 8:50 UTC on May 26th) and is used to illustrate the methodology and qualitatively compare the filtering and smoothing approaches. Dataset 2 contains GPS L1-only SI for ten additional rising occultation events (recorded May 4th through 6th) and is used to facilitate the simulation study to quantify the performance of the smoothing approach.
3 METHODOLOGY
For a 0.5 chip spacing between the early, prompt, and late correlators in the tracking loop of the direct signal, multi-path effects can persist until the path length difference is 1.5 chips (Misra & Enge, 2010). So, elth is the elevation angle of the transmitter relative to the receiver for which the path length difference is 1.5 chips. When the transmitter is below this elevation, the reflected signal can be present as multipath in the tracking results of the direct signal. This threshold depends primarily on the height of the receiver above the surface and the chipping rate of the signal. As a first approximation, elth is determined by approximating the Earth as a sphere and assuming that the incident ray of the reflection is parallel to the direct signal path (see Figure 2).
Geometry used to determine the path length difference as a function of satellite elevation angle
Here, re is the radius of the Earth, el is the elevation angle of the transmitter relative to the receiver, and h is the receiver’s height above the surface. The path length difference is the difference between the length of the line segments a and b:
First, the angle θ is determined by numerically solving the equation
The length of the line segments are
The path length difference is therefore
The theoretical path length difference as a function of the satellite elevation angle is plotted in Figure 3.
Path length difference as a function of satellite elevation angle, calculated assuming spherical Earth and parallel incidence
The thresholds are set at 1.5 times the chip width: L1 and L2 at 439.57 m and L5 at 43.96 m. Based on these approximate results, elth is 4.03˚ for L1 and L2 and –0.41˚ for L5.
With these theoretical results to serve as a point of comparison, the relevant geometrical parameters are now determined for dataset 1 (see Figure 4). The path length difference is calculated assuming specular reflection. The specular point is the point on the Earth’s surface that minimizes the total reflected path length (Jales, 2012). This point is determined by approximating the Earth’s surface using the WGS84 ellipsoid and solving the minimization problem. Once the location of the specular point is known, the path length difference can be easily determined. Figure 4 also includes the time rate of change of the path length difference, which is later used to model the interferometric frequency introduced by the reflected signal acting as multipath.
Relevant geometrical parameters for the time interval of interest in dataset 1. The times are referenced to the beginning dataset
The value of elth is 3.94˚ for L1 and L2 and –0.43˚ for L5, which agrees well with the earlier results found using the spherical Earth approximation. According to these results, multipath should no longer be present in the tracking results of the direct signal after the 30-minute mark for L5 and after the 45-minute mark for L1 and L2.
Spectrograms of the tracked SI are used to investigate the mutual interference between the direct and reflected signals. To understand the spectral content of the SI, it is necessary to develop a model of the power of the received signal. Here, the composite received signal s(t) is modeled as a superposition of the direct and reflected signals (e.g., in Byun, Hajj, & Young, 2002):
where Ad(t) and φd(t) are the time-dependent amplitude and phase of the direct signal, and Ar(t) and φr(t) are the time-dependent amplitude and phase of the reflected signal. The phase terms are expressed in units of cycles. This model assumes that the signal is only reflected from the specular point. The power of s(t) is
The argument of the cosine term corresponds to the interferometric frequency that appears in the spectrograms. This frequency follows the time rate of change of the phase difference (equivalently, the time rate of change of the path length difference divided by the signal wavelength):
The time rate of change of the path length difference depicted in Figure 4 is normalized by the carrier wavelength of the corresponding signal so that it can be compared directly to the spectrogram. Figure 8 shows the spectrograms for dataset 1 with the theoretical interferometric frequency overlaid.
The other primary goal of this research is to devise a technique to separate the contributions to the SI of the direct and reflected signals. A simple approach is to apply a filter to the SI with a 0.5 Hz cutoff (determined visually from the spectrogram). Two different filters are considered: a third-order Butterworth filter and a fifth-order Chebyshev type I filter. Figure 5 illustrates the L1 SI spectrogram and the magnitude response of the filters used to remove or isolate the multipath contribution. To remove the multipath component and retain tropospheric features, a low-pass filter is used. To retain the multipath component and remove tropospheric features, a high-pass filter is used. In all cases, the filter is applied forwards and backwards to yield a linear phase response. The filter magnitude response depicted in Figure 5 is of the forwards-backwards filters (i.e., the square of the magnitude response of the filter applied once). A 1 dB ripple in the passband of the Chebyshev filter is arbitrarily chosen—this value could be tweaked if more specific requirements for the accuracy of the filtered result were known. In the end, the Chebyshev filter is chosen because of its sharper cutoff.
The L1 SI spectrogram from dataset 1 and the magnitude response of the forwards-backwards filters (Butterworth and Chebyshev type I)
The other technique devised for removal or isolation of the multipath contribution is a time-domain smoothing approach. Here, the peaks and troughs of the SI are located using a wavelet-based peak-finding algorithm (Du, Kibbe, & Lin, 2006). Splines are then used to interpolate the points of the maxima and minima to yield the upper and lower envelopes of the SI. As illustrated in Figure 6, which uses an example of the L1 SI from dataset 1, the smoothed SI is the average of the upper and lower envelopes.
The time-domain smoothing technique involves finding the peaks and troughs of the SI, interpolating curves through these points, and taking the average of the upper and lower envelopes
The later results show that the filtering technique successfully removes the multipath contribution when the interferometric frequency is clearly distinguishable. The smoothing approach shows promise at lower elevations where the filtering technique begins to fail. A simulation study using dataset 2 is performed to quantify the performance of the smoothing approach. First, the SI for a direct signal that incorporates realistic tropospheric scintillation is simulated by selecting a portion of the low-pass filtered SI where the multipath is clearly removed. Then, the SI of the combined signal is simulated using the aforementioned model of the composite received signal and computing its power. In this case, Ad(t) is the square root of the simulated direct SI, Ar(t) is Ad scaled by a constant factor determined by the Fresnel reflection coefficient, and the phase difference between φd(t) and φr(t) is determined by ΔR for the given event. The performance of the smoothing technique is measured by how well the simulated direct signal can be extracted from the simulated combined signal. This is quantified by calculating the percent error over a range of elevation-angle bins. The simulation procedure is illustrated in Figure 7.
Steps of the simulation procedure. This is illustrated here using the L1 SI from dataset 1; the direct signal is simulated using the low-pass filtered SI between minutes 28 and 38
4 RESULTS AND ANALYSIS
First, the single event from dataset 1 containing SI for GPS L1, L2, and L5 is used to qualitatively compare the interference mitigation techniques. The spectrograms of the raw SI and the theoretically determined interferometric frequency curves are depicted in Figure 8. (Although only the spectrograms shown here are for dataset 1, the L1 SI spectrograms for the ten events in dataset 2 exhibit the same behavior.) Recall that the multipath should drop out by 45 minutes for L1 and L2 and by 30 minutes for L5. In all three cases, the multipath contribution fades out before this point. For L5, the multipath contribution is hardly distinguishable from the fluctuations caused by tropospheric scintillation (the faint trend in L5 after 30 minutes that matches the interferometric frequency curve is likely due to sidelobes in the autocorrelation function of the code). So, a threshold elevation angle has been successfully established for the geometry of the mountaintop receiver.
Spectrograms of the raw SI for GPS L1, L2, and L5 signals for dataset 1
There are several other interesting features of the spectrograms that should be noted. Between 20 and 25 minutes, the GNSS satellite is barely above the horizon. During this time interval, the direct and reflected signals undergo a significant amount of bending, resulting in a discrepancy between the true interferometric frequency and the theoretical curve. A model that incorporated bending effects would be able to more faithfully reproduce the frequency of the multipath contribution. Also note the difficulty of distinguishing between the frequency fluctuations due to tropospheric scintillation and the multipath contribution between 20 and 25 minutes. After 25 minutes, the GNSS satellite is higher above the horizon, and the frequency of the multipath contribution is well-predicted by the interferometric frequency curve. At this point, it is also easier to distinguish between the multipath contribution and the direct signal contribution.
The effectiveness of the separation techniques is analyzed first by looking at the resulting SI spectrograms. Figure 9 shows the SI spectrograms after the low-pass Chebyshev filter is applied. The obvious multipath contribution has been removed after minute 25. However, between 20 and 25 minutes, the frequency of the multi-path contribution is below 0.5 Hz. Here, the filter cannot remove the multipath, and thus, the SI still contains contributions from the direct and reflected signals.
Spectrograms of the SI from dataset 1 after applying the low-pass Chebyshev type I filter
Figure 10 shows the spectrograms for the smoothed SI. The smoothing technique is analogous to low-pass filtering because it attempts to remove the multipath contribution. Again, the multipath contribution is removed after 25 minutes. Between 20 and 25 minutes, it appears that this approach has had more success in removing the multipath than the filtering approach. The extent of this success is analyzed later via simulation.
Spectrograms of the smoothed SI for dataset 1
Next, the performance of the mitigation techniques for multipath isolation is considered. The spectrograms of the SI after applying the high-pass Chebyshev filter are depicted in Figure 11. Again, the filtering approach is successful after 25 minutes where the multipath contribution is more easily distinguishable from the direct signal contribution. Between 20 and 25 minutes, the filter fails to isolate the multipath contribution, removing it instead.
Spectrograms of the SI for dataset 1 after applying the high-pass Chebyshev type I filter
The multipath contribution can also be isolated by subtracting the smoothed SI from the raw SI. Figure 12 depicts the spectrograms for the smoothing approach for multipath isolation. In this case, it appears that the multipath contribution is retained even between 20 and 25 minutes.
Spectrograms of SI minus smoothed SI for dataset 1
The difference between the filtering and smoothing approaches can be further illustrated by analyzing the SI time series. The time series between 22 and 22.5 minutes is shown in Figure 13(a). Here, the low-pass filter fails to remove the fluctuations corresponding to the multipath because the frequency of these fluctuations is smaller than 0.5 Hz. The smoothing approach, on the other hand, succeeds in removing the fluctuations. Figure 13(b) shows the L1 SI time series between 28 and 28.5 minutes. At this point, the SI fluctuations occur at a high enough frequency to be removed with the filtering approach. In this case, the results for the filtering and smoothing approaches are similar.
The filtering and smoothing approaches are compared by looking at the L1 SI time series from dataset 1 between (a) 22 and 22.5 minutes, and (b) 28 and 28.5 minutes
Finally, the ten events from dataset 2 are used to perform the simulation study to further explore the performance of the smoothing technique. Each data point in Figure 14 represents the average percent error between the smoothed result and the original simulated direct signal across a 0.1-degree elevation-angle bin. The thicker black curve is the average across all events. In the −1- to 0-degree elevation-angle range, the smoothing technique can recover the direct signal SI with less than 5% average error—this is a clear improvement over the filtering technique, which cannot reliably perform below 0 degrees elevation.
The percent error quantifies how well the smoothing technique can retrieve the direct SI. Each of the lighter colored curves represents the results for an event from dataset 2. The thicker black curve is the average across all events
5 CONCLUSION
In this paper, an elevation-angle threshold was defined and calculated as a means of determining when the reflected signal will be present as multipath in the tracking results of the direct signal. Furthermore, the first steps were taken in devising techniques to mitigate the mutual interference between these signal contributions that occurs below this threshold. In general, the filtering approach is only successful when the interferometric frequency is well-separated from the frequency contributions caused by tropospheric scintillation. The smoothing approach shows promise when the frequencies are not easily separable, a result corroborated by the simulation study. Still, further analysis is required to more fully characterize its performance.
In the future, a more robust mitigation technique will be developed and applied to other signal parameters in addition to the SI. For example, performing carrier phase altimetry requires knowing the phase difference between the direct and reflected signals. A more advanced technique would ideally be able to simultaneously estimate the amplitudes and phases of the direct and reflected signals that make up the composite received signal. The effects of tropospheric scintillation on the accuracy of the amplitude and phase estimation will also be considered.
HOW TO CITE THIS ARTICLE
Collett I, Jade Morton YT, Wang Y, Breitsch B. Characterization and mitigation of interference between GNSS radio occultation and reflectometry signals for low-altitude occultations. NAVIGATION. 2020;67:537–546. https://doi.org/10.1002/navi.375
ACKNOWLEDGMENTS
The data collection experiment was funded by AFRL grant #FA8650-14-1735. The data was collected by Steve Taylor and Harrison Bourne from the Satellite Navigation and Remote Sensing Lab at the University of Colorado Boulder. Dr. Frank van Graas from Ohio University and Mr. Neeraj Pujara from AFRL also assisted in the data collection experiment.
- Received May 23, 2019.
- Revision received April 20, 2020.
- Accepted May 10, 2020.
- © 2020 Institute of Navigation
This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.
Jump to section
Related Articles
Cited By...
- No citing articles found.