Single-Pass Soil Moisture Retrieval Using GNSS-R at L1 and L5 Bands: Results from Airborne Experiment

Abstract

Global Navigation Satellite System—Reflectometry (GNSS-R) has already proven its potential for retrieving a number of geophysical parameters, including soil moisture. However, single-pass GNSS-R soil moisture retrieval is still a challenge. This study presents a comparison of two different data sets acquired with the Microwave Interferometer Reflectometer (MIR), an airborne-based dual-band (L1/E1 and L5/E5a), multiconstellation (GPS and Galileo) GNSS-R instrument with two 19-element antenna arrays with four electronically steered beams each. The instrument was flown twice over the OzNet soil moisture monitoring network in southern New South Wales (Australia): the first flight was performed after a long period without rain, and the second one just after a rain event. In this work, the impact of surface roughness and vegetation attenuation in the reflectivity of the GNSS-R signal is assessed at both L1 and L5 bands. The work analyzes the reflectivity at different integration times, and finally, an artificial neural network is used to retrieve soil moisture from the reflectivity values. The algorithm is trained and compared to a 20-m resolution downscaled soil moisture estimate derived from SMOS soil moisture, Sentinel-2 normalized difference vegetation index (NDVI) data, and ECMWF Land Surface Temperature.

1. Introduction

Soil is a natural reservoir of water, making it the main supply store for plants to live. The surface storage is mainly depleted by the natural process of evaporation, percolation to lower layers in the soil, water uptake by plants, etc. Conversely, moderate-to-high surface soil moisture (SM) values increase flood [1] risks and affect soil erosion [2] by wind and rain. Consequently, monitoring the soil moisture content of this near-surface layer of soil is crucial for sustainable irrigation of crop fields (smart irrigation), forest fire risk prediction, assessment of vegetation senescence, and to have a better knowledge of the water cycle, which plays a key role in the climate feedback loops [3].

Soil moisture can be measured using in-situ probes, or by means of remote sensing techniques, for which several approaches have been shown to have the ability to retrieve surface soil moisture states at different spatiotemporal scales. Using L-band microwave radiometry, the European Space Agency (ESA) Soil Moisture and Ocean Salinity (SMOS) mission [4] and the National Aeronautics and Space Administration (NASA) Soil Moisture Active Passive (SMAP) mission [5] are providing soil moisture maps at a native resolution of ∼55 km [6] and 36 km [7], respectively. Thermal Infrared spectrometers have also been used to estimate soil moisture, although with lower accuracy than microwave radiometers [8]. Moreover, radar scatterometers (e.g., [9]) and Synthetic Aperture Radars (SAR) (e.g., [10,11,12]) at L- and C-bands have been used to derive soil moisture indexes. More recently, Global Navigation Satellite System—Reflectometry (GNSS-R) has been proven to allow for the retrieval of soil moisture from ground [13,14], airborne [15,16,17], and spaceborne [18,19,20,21,22,23,24] configurations.

GNSS-R offers the promise of an enhanced spatial resolution when compared to microwave radiometry sensors. The spatial resolution of GNSS-R receivers is mostly linked to the size of the first Fresnel zone [25], where for an airborne receiver at an altitude of ∼1000 m the spatial resolution ranges from 7.2 m/cos(${\theta }_i$) to 16.6 m/cos(${\theta }_i$) [26] and for a spaceborne receiver it scales with the squared root of the height, i.e., at 500 km height it is 22.4 times larger. Such high resolution from space has been analyzed in [27], where the GNSS-R signal collected by CyGNSS allowed the detection of riverbeds of a width of 200–250 m. However, current methods to retrieve SM from GNSS-R space-borne data are not providing such resolution. As discussed in [28], few reflections (less than 16%) contain a noticeable coherent component, thus showing that incoherent scattering is dominant over land. Under these conditions, the spatial resolution is degraded, as many contributions coming from the entire glistening zone are collected by the receiving antenna (i.e., 25–37 km as shown in [29]). Because of that, many algorithms to retrieve SM using GNSS-R data require averaging the GNSS-R observable to the SMAP native resolution of 36 km [23,30] or by applying spatial and temporal averaging (i.e., 22 km resolution in [31]). Other works have shown an enhanced spatial resolution as compared to the previous ones. As CyGNSS data is now tagging reflections containing a large coherent component [32], new algorithms are being developed providing large spatial resolutions up to 2 km [33] or 3 km [34] just including coherent reflections retrieved by CyGNSS. It is worth mentioning that the use of coherent reflections and an Artificial Neural Network (ANN) based algorithm is key to provide a high-resolution SM product [34].

Machine learning algorithms, and in particular ANNs, are now the latest approach to retrieve soil moisture from a wide range of remote sensing techniques. For instance, the SMOS soil moisture product assimilated by the European Centre for Medium-range Forecast (ECMWF) is obtained using an ANN [35]. However, ANNs challenge is to correctly train the algorithm: in this process, the key is not only to select the correct target but also the amount of ancillary data inputs used for the algorithm. Known ANN implementations are summarized in Table 1 from [30], and they require a large spatial averaging and the use of ancillary data due to the incoherent scattering that is produced over land. In most of the cases discussed in [30], the ancillary data used is either the normalized difference vegetation index (NDVI), the SMAP vegetation optical depth (VOD), or a combination of the soil texture and topography data.

Current reflectivity models (Equation (1) [36]) compensate for the vegetation attenuation and the surface roughness according to:

$$\Gamma (\theta )=R{(\theta )}^2{\gamma }^{2}exp(-4k^2{\sigma }_h^{2}co{s}^2(\theta )),$$

(1)

where $\Gamma$ is the reflectivity, $\theta$ stands for the local incidence angle, R corresponds to the amplitude of the Fresnel reflection coefficient, $\gamma$ is the transmissivity—which accounts for the vegetation attenuation—and it is modeled by the VOD or the NDVI as a proxy, k is the wavenumber (i.e., $\frac{2\pi }{\lambda }$) and ${\sigma }_h$ is the surface root-mean square (RMS) height.

If any of the terms included in Equation (1) is not estimated accurately, the soil moisture cannot be properly retrieved. This is shown experimentally in [17,37]: at L-band, soil moisture produces a change on the retrieved reflectivity of up to 17 dB for a range between 0 and 0.45 m${}^3/$m${}^3$, while under coherent scattering conditions, the surface roughness effect by itself may reduce the reflectivity up to 18 dB for a local RMS surface height variation of 0–4 cm [37] or up to 11 dB when a Kirchhoff Approximation simulator is applied to the multielevation surface [38]. Furthermore, vegetation may reduce the reflectivity by up to 11 dB for a VOD variation of 0–0.6. Therefore, it is critical to adequately estimate and include both parameters into the retrieval algorithm. While the vegetation impact can be, in principle, more easily corrected by means of the NDVI or the VOD, this is not the case for the surface roughness, as the effective surface RMS height is extremely complicated to be accurately modeled or retrieved. Previous studies from an airborne platform [17] showed that, even when estimating the surface roughness using laser profilers on ground, the reflectivity could not be properly corrected for, as the reflection of the GNSS signal actually takes place in the near soil surface, whose depth varies depending on the soil moisture content of the reflection area, and the signal wavelength [39]. Therefore, the effective RMS height also depends on the soil moisture content [40].

This study analyzes the effects of the surface roughness variations using GNSS-R data from the Microwave Interferometer Reflectometer (MIR) instrument [41,42], acquired during two flights across an extensively used agricultural area New South Wales, Australia. MIR is a state-of-the-art GNSS-R instrument designed to work as an interferometric GNSS-Reflectometer [43], but raw data was also stored at 32 MS/s at 1 bit. These data have been postprocessed following the conventional GNSS-R (cGNSS-R) technique at different integration times, as presented in [44,45] for the MIR flights over the ocean. This work shows the very first results of the MIR instrument over land, and it is organized as follows: Section 2 describes the data set used and the ancillary data used to calibrate/validate the instrument performance. Section 3 analyzes the statistics of the reflectivity of the different MIR flights over Yanco, where the surface roughness effect is evaluated at different integration times; its effect on the soil moisture retrieval is also evaluated. Finally, Section 4 shows the results of an ANN algorithm applied to the MIR data to retrieve soil moisture.

2. Data Description

MIR flew over the Yanco-designated portion covered by OzNet [46], New South Wales, Australia, on 1 May 2018, and on 18 June 2018 (see Figure 1). The two flights covered the same area, but the first flight was conducted after a long period without rain, under very dry soil conditions, while the second flight was conducted the day after a rainy day, in which the water content of the soil was substantially higher. These two flights are called herein as “Dry” and “Wet” flights, respectively.

2.1. Ground-Truth and Ancillary Data

MIR covered Yanco areas A and B (see Figure 1) in the two flights. As shown in Figure 2, the OzNet soil moisture of both days is significantly different. In the area covered by the plane track, the ground stations showed an average SM ∼0.05 m${}^3/$m${}^3$, and ∼0.27 m${}^3/$m${}^3$, for the dry and wet flights, respectively. However, considering only a limited set of in situ measurements is not enough to compare it with the GNSS-R measurements. For that reason, other remote sensing products have been used as reference data.

Two ancillary SM products are key for the correct interpretation of the data in this work. The first ones are the SMOS and SMAP soil moisture products. However, as the resolution is too large for our data set (∼55 and 36 km), a pixel downscaling [47,48] of the coarse resolution SMOS L3 soil moisture product has been enhanced down to 20 m, using SMOS brightness temperatures, ECMWF land surface temperature, and visible and near infra-red (VNIR) data from Sentinel-2. This technique has been validated with different soil moisture networks [49]. Moreover, the downscaled data have been validated against the OzNet SM network at both “Dry” and “Wet” flights, showing a bias of 0.01 m${}^3/$m${}^3$, with a root-mean-square deviation (RMSD) of 0.032 m${}^3/$m${}^3$ for the “Dry” flight, and a bias of 0.02 m${}^3/$m${}^3$, with a RMSD of 0.038 m${}^3/$m${}^3$ for the “Wet” flight. The RMSD for both maps is lower than the SMOS accuracy at native resolution (0.04 m${}^3/$m${}^3$), and it is lower than the RMSD presented in other studies [50]. Therefore, the presented downscaled product is suitable as ground-truth data.

As it can be seen in Figure 3 and Figure 4, the NDVI from Sentinel-2 presents some differences between the Dry and Wet flights, and the downscaled SMOS soil moisture product also presents large differences between both flights. However, the important values for this study are the ones collocated with the GNSS specular reflection points. For that reason, both the NDVI and the SM maps are interpolated to each of the specular points of the MIR instrument, and the histograms of both measurements for the two flights are shown in Figure 5.

2.2. GNSS-R Data

The MIR uplooking and downlooking antenna arrays are composed of 19-element dual-band patches. Each array creates 2 beams at L1/E1 and 2 beams at L5/E5a using an automatic beam-steering algorithm that compensates for the aircraft attitudes. Therefore, MIR is able to track up to four different GNSS satellites simultaneously with the uplooking antenna and four specular points with the downlooking antenna (two for each band). The instrument automatically tracks different satellites when the plane has a stable attitude in the Local-Vertical Local-Horizontal frame. In the case that a given satellite goes out of the antenna pattern, the algorithm automatically changes to a new one. Moreover, the instrument notifies whenever the plane is not having a stable attitude (i.e., during turns or maneuvers), and therefore flagged data need to be filtered out. The plane flew at an altitude of $h\sim$ 500 m at an average speed of $v\simeq$ 75 m/s. From this altitude, the size of the Fresnel zone $l_{F_r}$ is (Equation (2) [43,53]):

$$\begin{array}{c}{l}_{Fr}={\displaystyle \frac{\sqrt{\lambda R_r}}{cos({\theta }_{inc})}},\hfill \\ R_r={\displaystyle \frac{h}{cos({\theta }_{inc})}},\hfill \end{array}$$

(2)

where $\lambda$ is the signal wavelength, ($\lambda =19$ cm at L1/E1 and $\lambda =25$ cm at L5/E5a), and ${\theta }_{inc}$ denotes the incidence angle of the GNSS signal. Thus, for an incidence angle ${\theta }_{inc}=0^{\circ }$ (i.e., smaller Fresnel zone), $l_{F{r}_{L1}}=9.75$ m at L1/E1, and $l_{F{r}_{L5}}=10.75$ m at L5/E5a.

Based on the flight characteristics, the maximum integration time to prevent blurring of the first Fresnel zone is bounded by $T_{int}<{\displaystyle \frac{l_{F{r}_{L1}}}{v}}\simeq 130$ ms. For this study, and in order to have enough oversampling, $T_{int}=20$ ms has been selected. In that case, and depending on the local incidence angle, the number of samples overlapped in the same Fresnel zone varies from 6 to 11 samples, which corresponds to ${\theta }_{inc}=0^{\circ }$ and to ${\theta }_{inc}={45}^{\circ }$, respectively.

Because of the flight height and speed, the delay and Doppler spreads are negligible, and therefore the selected GNSS-R observable is the reflectivity as shown in Equation (1). Applying the Pseudo-Random Noise injection technique conceived in [54], the MIR instrument is calibrated as described in [55], including the antenna pattern compensation as in [56]. Furthermore, the reflectivity is calculated assuming only the peak of the waveform minus the noise according to Equation (3):

$$\Gamma ={\displaystyle \frac{P_{ref}-P_{N_{ref}}}{P_{dir}-P_{N_{dir}}}}\times {\displaystyle \frac{G_{up}({\theta }_{dir},{\varphi }_{dir})}{G_{dn}({\theta }_{ref},{\varphi }_{ref})}},$$

(3)

where P is the amplitude of the waveform at its maximum, $P_N$ is the noise power of the waveform, computed as the average of the delay bins before the peak of the Delay-Doppler Map, G is the gain of either the uplooking or the downlooking antenna, $\theta$ and $\varphi$ are the look angles to either the transmitting satellite (“dir”), and to the specular reflection point (“ref”). In this case, the GNSS satellite antenna pattern is not taken into account, neither the additional path losses, as the increased distance is negligible (500 m) with respect to the already traveled distance (>20,000 km).

The reflectivity is absolutely-calibrated using different water bodies that the MIR instrument crossed during both flights. As shown in Figure 6, the reflectivity over the selected water bodies is ∼−2.1 dB at incidence angle ∼20${}^{\circ }$. Figure 7 shows the histogram, normalized as a Probability Density Function (PDF), of the reflectivity. As it can be seen, mean reflectivity values have a difference between Dry and Wet flights of ∼7–8 dB. Considering Figure 21 from [37], the difference in received power is expected to be ∼4 dB for a difference of ∼0.22 m${}^3$/m${}^3$. However, Ref. [37] does not take into account the effect of the penetration depth of the GNSS reflection signal in different soil moisture conditions. In this work, the surface roughness plays a significant role in this averaging operation, as the average of a magnitude affected by an attenuation is a biased estimator of the actual magnitude. Additionally, as previously discussed, this roughness also depends on the soil moisture content [39], and therefore an analytic solution to that problem is very complex.

3. Soil Moisture Retrieval Using GNSS-R

As shown in Figure 23 from [37], and discussed in Section 1, herein a difference of just 4 cm in the local RMS height may produce a 18 dB drop in the power of the reflected GNSS signal. Furthermore, as the first Fresnel zone of the reflection is very small, the local surface roughness also has a high variability between overlapped nearby reflections.

3.1. Reflectivity Statistics Using Different Integration Times

The discussion of the previous section has shown that performing the average of the reflectivities (in linear units, and then taking its logarithm) of a certain area can reduce the effect of terrain inhomogeneity, and reduces the Speckle noise. In this section, the effect of increasing the integration time is addressed from a statistical point of view. In this case, consecutive incoherently-integrated reflectivity samples at $T_{int}$ = 20 ms are averaged (linearly, using a moving average filter, and then taking the logarithm of the resulting magnitude) up to 5 s. Different histograms for Beams 1 (L1) and 3 (L5) for both Dry and Wet flights are shown in Figure 8a,b, respectively for L1 and L5 cases. As it can be seen, as the effective integration time increases, the average value remains constant, but the standard deviation of the reflectivity decreases.

In addition, Figure 9 shows a selected track over a certain region where the local surface roughness of the crop field changes. As it can be seen, for shorter averaging times, a larger local variability of the signal is encountered, capturing small variations in the terrain, as in the middle of the image (outlined with a black box), where the signal crosses a small irrigation channel. However, as the effective integration time increases, small variabilities of the terrain (i.e., surface roughness) are averaged, and therefore the attenuation variability caused by the surface roughness is diluted, providing a negative bias in the reflectivity, which varies as a function of the integration time.

3.2. Surface Roughness Effect in Soil Moisture Retrievals

Since at high SM the reflectivity tends to saturate (Figure 21 of [37]), its sensitivity to SM decreases. In this case, if the scattering occurs close to the surface, the roughness effects on reflectivity, and consequently SM, are more noticeable. This is observed in Figure 9, in which the large variability in the reflectivity is mostly produced by local surface roughness variations. When the incoherent integration time is very short (i.e., 20 ms), the terrain contribution is not smoothed, and therefore the power changes caused by local roughness variations are highlighted.

In order to analyze in detail the effect of the surface roughness, the Wet flight over Yanco site B has been selected. The NDVI of the area is very low and almost constant, and the attenuation due to vegetation opacity can be neglected. The influence of the soil surface roughness on SM can be studied by isolating its corresponding term in Equation (1) and assuming a saturated reflectivity due to very high soil moisture values. Therefore, the surface RMS roughness can be estimated from Equation (4), as follows:

$${\sigma }_h=\sqrt{\frac{(-ln({10}^{\frac{\Gamma \left[dB\right]+5}{10}})}{(4k^2\times cos{({\theta }_{inc})}^2))}}.$$

(4)

The Fresnel reflection coefficient of Equation (1) is set to a constant value of −5 dB for a saturated SM value. As it is seen in [37], the GNSS-R reflectivity gets saturated for SM values larger than ∼0.25–0.3 m${}^3/$m${}^3$. Moreover, as it can be seen in Figure 7c from [21], these SM values produce a range of reflectivities from −6 to −4 dB. Thus, for the sake of simplicity, a reflectivity value of −5 dB has been selected for this study. Note that, this range of reflectivities is assuming a VOD < 0.2. As it is shown in Figure 2a [58], an L-band VOD < 0.2 is provided for NDVI values < 0.4.

Figure 10a shows the results of applying Equation (4) to the Yanco site B during the Wet flight, after removing all reflectivities from water bodies or with an NDVI > 0.4, so as to assume negligible attenuation due to vegetation. To perform this study, the selected beams at L1 and L5 are from the same GNSS spacecraft so that the specular reflection point is the same. In this case, for moist soils the penetration depth at both bands is limited to ∼2 cm [39]. As it can be seen, the estimated roughness at both bands is very similar, with an average value of ∼1.23 cm and a standard deviation of ∼0.68 cm at L1/E1, and an average value of ∼1.37 cm, and a standard deviation of ∼0.59 cm at L5/E5a, which is consistent with a slightly larger penetration depth at L5/E5a, and a reflection accuracy over a slightly flatter interface. In this case, the impact of an average roughness of 1.3 cm produces a degradation in the reflectivity of ∼1.4 dB, while a surface roughness of 2 cm produces a degradation up to 4 dB.

Finally, the roughness at L1 depending on the integration time is displayed in Figure 10b. As the effective integration time increases, the estimated surface roughness decreases to values lower than 1 cm. Averaging up to 1000 ms produces an average roughness of 0.88 cm, with a standard deviation of 0.45 cm, and averaging up to 5000 ms produces an average roughness of 0.78 cm and a standard deviation of 0.35 cm. Increasing the effective integration time reduces the surface roughness variability, but this introduces a bias in the estimation of the reflectivity.

4. Soil Moisture Retrieval Algorithms

Surface roughness is the variable limiting the accuracy of SM retrievals using GNSS-R data [17]. As introduced in Section 1, algorithms based on ANNs have proven to be very powerful tools to detect and solve nonlinear problems by minimizing the error of the algorithm output with respect to a known target. Despite that, up to now, most of the algorithms have shown a dependence on the ancillary data (e.g., [30]). In our proposed approach, the use of ancillary data has been reduced by adding statistical metrics of the reflectivity itself to the algorithm input, such as the standard deviation of the reflectivity. However, those statistics are computed from data collected at different time instants, and there is not yet an approach for single-pass retrieval. Despite that, this points out an interesting question: is there a relationship between the standard deviation of the reflectivity and the surface roughness?

4.1. Surface Roughness and Reflectivity Standard Deviation

To address this question, some examples are provided in both Dry and Wet flights over Yanco site B. The geolocated averaged $\Gamma$ and the moving standard deviation (movstd) of $\Gamma$ are computed and presented in Figure 11 and Figure 12.

Both figures present the $\Gamma$, averaged to 5000 ms, the movstd($\Gamma$) for the same integration interval. The movstd is calculated in linear units over N 20 ms integrated reflectivity measurements (i.e., for a $T_{int}=$ 5000 ms, N = 250), and then converted into dB units. Finally, the Noise-to-Signal Ratio (NSR), computed as movstd($\Gamma$) minus $\Gamma$, is presented in the third column. The color axis is evenly defined for the two flights and each of the three parameters.

As it can be seen in Figure 11, the two highlighted areas present a large decrease of the average reflectivity down to ∼−20 dB. The moving standard deviation is also affected, and the computed NSR in both areas increases. In this case, the increase in this NSR term can be linked to a rougher area, as in the bottom one. In this case, the reflectivity drop is linked to the loss of coherency caused by a vegetated area. However, it is important to remark that the average NSR is large: −0.78 dB. As it can be seen, lower NSR values correspond to areas with a larger reflectivity, and a lower reflectivity variation, which can be linked to smoother surfaces.

Moving to the Wet flight (Figure 12), the two previously highlighted areas are also shown, plus a third one in the middle of the map. First of all, it is important to remark that the average reflectivity is larger due to the higher soil moisture content of the flight. Moreover, the NSR is slightly lower than the Dry flight: −1.13 dB, and there are converging to the Speckle Noise NSR limit: −5.6 dB (p. 608 from [59]). This noise is the effect of the environmental conditions of the reflection scenario (i.e., surface facets, geometry, etc.), and it is a multiplicative noise that, among others, can be reduced by low-pass filtering, applying averages, or using neural networks [60]. The two selected sites of the Dry and Wet flights present a similar NSR, ∼1 dB for the top left one, and ∼2–3 dB for the one at the bottom right. Therefore, even though the soil moisture content is different, and the reflectivity value of the Dry flight is ∼7 dB lower than the Wet flight, the NSR is in the same order; therefore, the terrain inhomogeneity is similar in both cases.

Finally, Figure 13 presents a scatter density plot of the reflectivity computed at 20 ms, compared to the NSR computed at 5000 ms, for Dry and Wet cases, and at L1 and L5 bands. It is important to remark in the regime where the Speckle Noise is dominant (i.e., NSR close to −5.6 dB), the reflectivity tends to its average value, whilst for larger NSR values, the reflectivity displays a much larger variability due to terrain inhomogeneity and surface roughness.

The differences between the Dry and Wet flights are noticeable. In the Dry flight, the contribution from the surface roughness is dominant, and the NSR does not reach the −5.6 dB Speckle Noise SNR limit. On the contrary for the Wet flight, and especially at L5, the NSR presents a larger number of points with NSR lower than −4 dB, with some of the realizations in the −5.6 dB of the Speckle Noise SNR limit. The differences in NSR in both flights are linked to the surface roughness variation depending on the moisture content, and therefore the penetration depth of the incidence wave. The higher NSR of the Dry flight indicates a rougher terrain than the Wet flight. Furthermore, the flights at L5 also present a lower NSR than the ones at L1, due to the change in the penetration depth of the L5 signal, and the difference in the autocorrelation function width.

4.2. Artificial Neural Network

Due to the nonlinear behavior of $\Gamma$ and ${\sigma }_h$ in relation to SM, the soil moisture retrieval process does not have an analytical closed-form solution. However, the use of machine learning algorithms, and neural networks in particular, is a growing technique broadly used to solve nonlinear problems. In this case, selecting the proper inputs is crucial to accurately retrieve soil moisture. It has been shown in the previous section that the NSR is related in some way to the surface roughness, but as this parameter is computed from the subtraction of movstd($\Gamma$) and $\Gamma$, both parameters have been used separately in an ANN algorithm, therefore letting the network use both parameters independently. To do that, the following ANNs are proposed using the following four cases:

(1)

$\Gamma$, movstd($\Gamma$) as a proxy for ${\sigma }_h$, NDVI, and ${\theta }_{inc}$,

(2)

$\Gamma$, NDVI, and ${\theta }_{inc}$,

(3)

$\Gamma$, movstd($\Gamma$) as a proxy for ${\sigma }_h$, and ${\theta }_{inc}$,

(4)

$\Gamma$, and ${\theta }_{inc}$.

Furthermore, different integration times ($T_{int}$ = 0.1, 1, 2, and 5 s) are used, and the target output for all cases is the collocated and downscaled to 20 m SM data from SMOS/Sentinel-2 described in Section 2. Finally, the input data is split in two: one network is deployed for L1 and another for L5. Both networks are based on a three hidden-layer feed-forward network with 6 neurons each, and each data set (L1 and L5) is randomly divided into two parts: training and test. The training set is 20% of all the data available for both Dry and Wet flights, for the two Yanco A and B locations. The test set is the resulting 80%. This training set is the one used to train the network, randomly divided again in 70%, 15%, 15%, for training, validation, and testing of the network, respectively. In order to avoid neural network overfitting, early stopping and pruning techniques are applied to the trained network [61]. Finally, the results are applied to the overall data set. Note that the entire data set (i.e., both Yanco A and B areas, and all the possible NDVI values) are used in the algorithm to cover a wider reflectivity and terrain variability.

4.3. Results

Due to the environmental constraints of the two flights, the SM values of the Dry flight and the SM values of the Wet flight are quite distant. It is important to remark that the Dry flight was conducted after ∼1 month without rain and most of the area is not irrigated. Moreover, during the Dry flight SM values from 0.05 m${}^3/$m${}^3$ up to 0.1 m${}^3/$m${}^3$ are found. In contrast, the Wet flight was conducted after a strong rain event, and therefore most of the soil is very moist, in this case, values from ∼0.26 m${}^3/$m${}^3$ up to 0.33 m${}^3/$m${}^3$ are found.

Figure 14 and Figure 15 show the scatter density plot of the GNSS-R L1 SM output from the ANN, with respect to the SMOS/Sentinel-2 downscaled SM. The color axis is the density of points in logarithmic units. Analyzing Figure 14, it can be seen that the information provided by the NDVI increases the R parameter in any case. Moreover, the standard deviation of the error decreases, even without using the movstd($\Gamma$), but the best case is when this parameter is used. For short effective integration times, the R and the standard deviation of the error are slightly worse than the case for longer integration times. As the averaging increases, the estimation of the attenuation due to both the surface roughness and the vegetation are averaged, and the use of the movstd($\Gamma$) increases the correlation between the target and the ANN output. Furthermore, by looking at the shape of the output of the ANN, it can be seen the misclassification produced by the algorithm, where points falling into the Wet area, are classified as low soil moisture. This is mostly due to the effect of the reflectivity reduction due to surface roughness. In this case, by increasing the effective integration time (Figure 14 right, $T_{int}$ = 5 s), the attenuation is smoothed and the algorithm shows a better behavior, showing the lowest error. Finally, it is important to remark that, in the case where the reflectivity is used alone, the ANN output is the least accurate one, clearly indicating that additional information is required to retrieve soil moisture.

Moving to the L5 case (Figure 15), it can be seen that both the R and the standard deviation of the error are clearly better than in the L1 case. In any of the selected cases, the standard deviation of the error is ∼2–3 times lower than in the L1 case. Furthermore, even with a small averaging (i.e., $T_{int}=$ 0.1 s), the dispersion is smaller than at L1. Note that, these results are consistent with Figure 13, where the NSR at L5 was very close to the Speckle Noise limit, indicating that the L5 signal is less affected by surface roughness variations, especially during the Wet flight.

As the effective integration time increases, the ANN output shows lower error and a higher correlation coefficient with respect to the SMOS/Sentinel-2 downscaled soil moisture. Note that errors are drastically reduced for the case with larger averaging, where the standard deviation of the error is almost negligible (i.e., 0.016 m${}^3/$m${}^3$), showing a very low dispersion in both Dry and Wet flights.

4.4. Discussion

The presented results show a significant difference between L1 and L5 bands, with a standard deviation of the error at L1 being three times larger than at L5, despite the lower antenna directivity ($D_{L_1}$ = 21 dB, $D_{L_5}$ = 18 dB). This is due to a longer wavelength and a larger penetration depth, and by design, a much narrower autocorrelation function (30 m in space) at L5, which translates into a higher resolution. In this case, the peak of the L5 waveform contains contributions from a smaller glistening zone, increasing the coherency of the received signal. On the contrary, the L1 signal has a much larger autocorrelation function (300 m in space), and therefore contributions from a larger glistening zone are added in the L1 waveform, producing larger fluctuations than at L5.

Aside from the difference in bands, there is also a large variability depending on the selected integration time. In order to illustrate it, the same neural network is now deployed for $T_{int}=$ 0.1, 0.5, 1, 2, 3, 4, and 5 s. Results are shown in Figure 16.

Comparing cases 2 and 3 for large integration times, case 3 (i.e., using movstd($\Gamma$)) introduces larger errors than case 2 (i.e., using NDVI), which is not happening for lower integration times (i.e., 2000 ms at L1 and 1000 ms at L5). However, the combined use of the movstd($\Gamma$) parameter together with the NDVI provides the lowest error for all integration times. In this case, as the vegetated areas are quite small (see Figure 3), a very large integration time produces errors, as a $T_{int}=$ 5000 ms is equal to a specular point movement of ∼375 m. However, the lower the integration time, the larger the surface roughness effect. On the contrary, for the smallest integration time, 100 ms, where the plane movement is equal to the size of the first Fresnel zone, case 3 provides the same error as case 1, and a much smaller error than cases 3 and 4. Thanks to the oversampling in the along-track direction, several measurements of speckle noise and surface roughness are included in the recovery algorithm by means of the movstd($\Gamma$) term. However, for longer integration periods, e.g., $T_{int}=$ 2000 ms, the signal covers up to ∼15 Fresnel zones (i.e., $v_{plane}$ = 75 m/s, and $l_{F_r}\sim$10 m). The terrain inhomogeneity while crossing these areas together with the presence of different vegetated areas induces errors in the recovery algorithm, which can be only corrected using NDVI.

It is clear that there is a clear trade-off between spatial resolution and radiometric resolution, where it is not possible to achieve a good radiometric resolution and low root-mean-square error of the retrieved parameter and good spatial resolution at the same time. As shown in Figure 13 from [30], the standard deviation of their spaceborne GNSS-R SM retrieval algorithm is decreased to ±0.1 m${}^3/$m${}^3$ when more than 25 averages are performed, and the smaller the number of averages, the larger the error, and hence the worse the radiometric resolution.

As shown in this study, when increasing the integration time, and therefore lowering the spatial resolution, the standard deviation of the error decreases, but effects due to terrain changes in the along-track direction of the GNSS-R measurement induce errors that need to be corrected for using NDVI measurements. On the contrary, reducing the averaging leads to a much higher resolution, and the NDVI term is not providing additional information to the ANN algorithm. However, the standard deviation of the error at such low integration times is larger.

Just to remark that the radiometric and the spatial resolutions cannot be maximized at the same time and there will always be a trade-off (if no ancillary data is used) between the required SM error and the spatial resolution of the GNSS-R-derived SM product.

5. Conclusions

This study has analyzed the data collected by the MIR instrument during two flights (“Dry” and “Wet”) over the OzNet Yanco sites in New South Wales, Australia, during May–June 2018. The effect of increasing the averaging and its impact on the surface roughness estimation are addressed, showing that the effective integration time has to be increased to 5 s to neglect surface roughness effects. A statistical parameter based on the moving standard deviation over N samples of the reflectivity (movstd($\Gamma$)) has been presented as a proxy of the surface roughness effects when the averaging (N) is large enough. Finally, an ANN-based algorithm has been presented for different combinations of auxiliary data and reflectivity averages, for L1 and L5 cases. In both cases, the use of the movstd($\Gamma$) parameter reduces the error of the retrieved SM to 0.047 m${}^3/$m${}^3$ and 0.016 m${}^3/$m${}^3$ at L1 and L5 respectively, for a $T_{int}$ = 5000 ms. Furthermore, the L5 signal shows a larger correlation coefficient with the expected SM output than the L1 signal because of the higher penetration depth and the narrower autocorrelation function.