University of Southern California Dissertations and Theses

Electromagnetic scattering models for satellite remote sensing of soil moisture using reflectometry from microwave signals of opportunity

(USC Thesis Other)

Electromagnetic scattering models for satellite remote sensing of soil moisture using reflectometry from microwave signals of opportunity

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content ELECTROMAGNETICSCATTERINGMODELS FORSATELLITEREMOTESENSINGOFSOIL MOISTUREUSINGREFLECTOMETRYFROM MICROWAVESIGNALSOFOPPORTUNITY JAMES D. CAMPBELL A DISSERTATION PRESENTED TO THE FACULTY OF THE USC GRADUATE SCHOOL IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY (ELECTRICAL ENGINEERING) MING HSIEH DEPARTMENT OF ELECTRICAL AND COMPUTER ENGINEERING UNIVERSITY OF SOUTHERN CALIFORNIA LOS ANGELES, CA, USA DECEMBER 2019 Contents Abstract iv Acknowledgements v Notation vi 1 Introduction 1 2 Review of Related Work 4 2.1 Remote Sensing Principles . . . . . . . . . . . . . . . . . . . . . . 4 2.2 Missions with Operational Soil Moisture Products . . . . . . . . 4 2.3 Soil Moisture from Other Non-Re ectometry Missions . . . . . . 6 2.4 GNSS Re ectometry . . . . . . . . . . . . . . . . . . . . . . . . . 7 2.5 Re ectometry with Communication Satellite Signals . . . . . . . 10 2.6 Theoretical Models for Rough Surface Scattering . . . . . . . . . 11 3 Sensitivity Study 14 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 3.2 Statistical Modeling . . . . . . . . . . . . . . . . . . . . . . . . . 15 3.2.1 Observation Model . . . . . . . . . . . . . . . . . . . . . . 15 3.2.2 Cost Function and Error Covariance Matrix . . . . . . . . 16 3.2.3 Partial Derivatives . . . . . . . . . . . . . . . . . . . . . . 17 3.2.4 Inner Product Representation . . . . . . . . . . . . . . . . 18 3.2.5 Soil Moisture Retrieval Performance . . . . . . . . . . . . 18 3.3 Numerical Results . . . . . . . . . . . . . . . . . . . . . . . . . . 20 3.3.1 Parameters . . . . . . . . . . . . . . . . . . . . . . . . . . 20 3.3.2 Sensitivities . . . . . . . . . . . . . . . . . . . . . . . . . . 21 3.3.3 Calibration Accuracy Requirement . . . . . . . . . . . . . 22 3.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 3.5 Future Directions . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 4 Models for Planar Rough Surfaces 26 4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 4.2 Theoretical Modeling . . . . . . . . . . . . . . . . . . . . . . . . . 26 ii 4.2.1 Total Power . . . . . . . . . . . . . . . . . . . . . . . . . . 26 4.2.2 Bistatic Radar Cross Section . . . . . . . . . . . . . . . . 27 4.2.3 Woodward Ambiguity Function . . . . . . . . . . . . . . . 27 4.2.4 Coherent Component . . . . . . . . . . . . . . . . . . . . 28 4.2.5 Noncoherent Component . . . . . . . . . . . . . . . . . . 28 4.2.6 Change of Polarimetric Basis from Linear to Circular for NBRCS . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31 4.2.7 Fresnel Footprint . . . . . . . . . . . . . . . . . . . . . . . 32 4.2.7.1 Elliptical Approximation of First Fresnel Zone . 32 4.2.7.2 Fresnel Footprint of Noncoherent Integration . . 34 4.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 4.3.1 Kati Thanda-Lake Eyre . . . . . . . . . . . . . . . . . . . 36 4.3.1.1 Dataset Selection . . . . . . . . . . . . . . . . . 36 4.3.1.2 Coherent Modeling . . . . . . . . . . . . . . . . 37 4.3.1.3 Noncoherent Modeling . . . . . . . . . . . . . . . 40 4.3.2 SoilSCAPE Sites . . . . . . . . . . . . . . . . . . . . . . . 44 4.3.2.1 Background . . . . . . . . . . . . . . . . . . . . . 44 4.3.2.2 Fresnel Footprint . . . . . . . . . . . . . . . . . . 46 4.3.2.3 Tonzi Ranch Dataset Identication . . . . . . . . 47 4.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53 5 Models for Rough Surfaces with Topography 55 5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55 5.2 Theoretical Model . . . . . . . . . . . . . . . . . . . . . . . . . . 56 5.2.1 Deterministic Modeling . . . . . . . . . . . . . . . . . . . 56 5.2.2 Stochastic Modeling . . . . . . . . . . . . . . . . . . . . . 61 5.2.3 CYGNSS-Specic Modeling . . . . . . . . . . . . . . . . . 64 5.2.3.1 Surface Bistatic Radar Cross Section . . . . . . 64 5.2.3.2 Modeling the Woodward Ambiguity Function . . 64 5.3 Validation Case Study . . . . . . . . . . . . . . . . . . . . . . . . 65 5.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70 5.5 Future Directions . . . . . . . . . . . . . . . . . . . . . . . . . . . 71 5.5.1 Land Cover . . . . . . . . . . . . . . . . . . . . . . . . . . 71 5.5.2 Calibration and Site Selection . . . . . . . . . . . . . . . . 71 5.5.3 Geometrical Optics Limit . . . . . . . . . . . . . . . . . . 71 5.5.4 Fluctuations . . . . . . . . . . . . . . . . . . . . . . . . . 72 5.5.5 Heterogeneities of Non-DEM Parameters . . . . . . . . . 72 5.5.6 Gradient and Gradient Error Estimation . . . . . . . . . . 72 5.5.7 Computational Eciency . . . . . . . . . . . . . . . . . . 73 A Partial Derivatives of Fresnel Re ectivity 74 Acronyms 77 References 82 iii Abstract Global land surface observations of soil moisture from space are key to under- standing water, energy, and carbon cycles. As such, soil moisture observations support many applications of human interest, including weather and climate forecasting, ood and landslide prediction, aquifer modeling, drought analysis, crop productivity evaluation, wildre danger assessment, and human health. Cur- rent soil moisture observation missions rely on large, expensive Earth-observing instruments, such as radiometers and scatterometers, which tend to lack spatial and temporal resolution. An alternative approach for Earth observation is the use of re ectometry from signals of opportunity (SoOp-R). In particular, SoOps from global navigation satellite systems and digital communication satellites are continually re ected from Earth's surface, picking up information about soil moisture and other biogeophysical parameters in the process. Due to their low weight and power, SoOp-R receivers can y on small satellites with the potential for improved resolution and revisit rates at a lower cost. The aim of the present research is to develop physics-based models for satellite remote sensing of land surfaces using SoOp-R. We begin by investigating the theoretical accuracy of a time series method for soil moisture retrieval without relying on any non-SoOp-R soil moisture products to calibrate each resolution cell. The method exploits the incidence angle diversity that is inherent in SoOp-R data. The method assumes coherent scattering of electromagnetic waves in the presence of instrument calibration errors, vegetation cover uncertainty, and surface roughness uncertainty. We conclude that the SoOp-R-only retrieval problem for current missions, which have only a single frequency and a single polarization, is poorly conditioned in general but may be feasible in the case of a bare rough surface. Next, we add a noncoherent scattering model based on the small perturbation method (SPM) and we attempt to validate the models with satellite data. Initial results from a large dry lakebed using the coherent model are promising. However, we nd that contribution of the noncoherent SPM is always negligible relative to that of the coherent model within their common regime of validity. Furthermore, for validation sites having in-situ sensors, we nd that both models fail due to topography within the regions of active scattering. Finally, we develop a new model for the scattering of electromagnetic waves from land surfaces with topography. The model makes use of a digital elevation model and represents an extension of a well-known existing model for scattering from a rough sea surface. Initial validation results are encouraging, and areas for future work are identied. iv Acknowledgements I am grateful to my faculty advisor M. Moghaddam for providing a stimulating research environment, introducing me to a global research community, and guiding me throughout my Ph.D. journey. I thank M. Moghaddam and my other dissertation committee members, A. Tabatabaeenejad and K. T. Sanders, for their thoughtful feedback and suggestions. My colleague A. Melebari contributed to the implementation of the topography- based model, including optimization by parallel computing and the coding of data quality checks and gradient calculations. R. H. Chen provided soil parame- ters from Soil SURvey GeOgraphic (SSURGO), and J. B. Whitcomb furnished gradient subroutines. I thank V. U. Zavorotny and J. T. Johnson for their technical comments on the electromagnetic scattering models that I presented at the Cyclone Global Navigation Satellite System (CYGNSS) science team meeting in Ann Arbor, Michigan, on June 5, 2019. I am indebted to A. J. O'Brien and S. Gleason for their guidance on how to register the delay and Doppler of CYGNSS data. My Ph.D. work was funded through a fellowship from the Engineering Advanced Study Program (ASP) of my employer, Raytheon Space and Airborne Systems. My former supervisor V. C. Kirk introduced me to the ASP and motivated me to apply. My former ASP advisor, the late W. P. Ballance, provided invaluable advice, strategy, and encouragement through my rst few years. Finally, I thank my ASP advisor R. D. Schaefer for facilitating the renewal of my fellowship for the nal phase of my research. This research was performed at the University of Southern California, sup- ported in part by the National Science Foundation under grant number 1643004 and in part by the National Aeronautics and Space Administration under grant number 80NSSC18K0704. Although I received no direct support from these contracts, they helped to cover my travel expenses for attending engineering conferences and science team meetings where I presented my progress. The contracts also provided partial support for my faculty advisor, who oversaw the research. v Notation Cartesian vectors are denoted by a right over arrow, and unit vectors by a hat. Non-Cartesian vectors, such as state vectors from control theory, are written in boldface. Estimated non-Cartesian vectors are written in boldface with a hat, in accordance with control systems convention. Thus, the hat can have two dierent meanings, depending on context. Matrices are written in uppercase boldface. Italic symbols denote variables or quantities, while upright symbols indicate non-variables or labels. The asterisk denotes complex conjugation. A superscript T indicates matrix transposition. Angle brackets denote either inner product or ensemble mean, depending on context. Over bars suggest averaging or approximation by a constant over some region. While the notation for much of the original literature on wave scattering from rough surfaces comes from the eld of acoustics, we prefer to use notation from the eld of electromagnetics for consistency with other work in bistatic radar. We note that, in the electronic version of this manuscript, many items are hyperlinked for the reader's convenience, including acronyms, table-of-contents entries, digital object identiers (DOIs), web addresses, and references to equa- tions, citations, gures, tables, section numbers, etc. vi Chapter 1 Introduction Soil moisture is a fundamental quantity in the study of our planet, providing observability into processes of evapotranspiration and groundwater recharge. These processes, in turn, in uence cycles of water, energy and carbon on re- gional and global scales. Applications supported by soil moisture observation include weather and climate forecasting, ood and landslide prediction, aquifer modeling, drought analysis, crop productivity evaluation, and human health [1]. Additionally, soil moisture shows promise for wildre danger assessment [2]. The advancement of dependable soil moisture retrieval methods is therefore of considerable human interest. Soil moisture observation requirements are not adequately satised by existing or planned satellite missions. The two most capable resources for accurate global soil moisture determination are currently the European Space Agency (ESA) Soil Moisture and Ocean Salinity (SMOS) mission and the National Aeronautics and Space Administration (NASA) Soil Moisture Active Passive (SMAP) mission. Both SMOS and SMAP carry microwave L-band radiometers capable of measuring near-surface soil moisture through moderate vegetation cover to an accuracy of 0:04 m 3 m 3 with spatial resolutions of about 40 km and return intervals of 2 to 3 days. However, spatial resolution needs to be increased to sub-kilometer scales and revisit rates need to be increased to at least once per day in order to achieve current Earth science goals [1]. Additionally, increased sensing depth is desired since surface soil moisture measurements from L-band microwave radiometry currently need to be extrap- olated down to the root zone to observe plant transpiration. Furthermore, neither SMOS nor SMAP can penetrate dense vegetation. Finally, no follow-on or replacement missions have been planned to extend the soil moisture data record once the existing missions reach the end of their useful lives. The value of soil moisture data grows with the length of the data record since meaningful monthly and seasonal statistics cannot be achieved without an extensive satellite 1 record. In short, soil moisture observation faces challenges of spatial and tem- poral resolution, vegetation cover penetration, sensing depth, and data record continuity [1]. An alternative approach, which has been receiving increasing attention in recent years and which has potential to address these challenges, is the use of re ectometry from signals of opportunity (SoOp-R). Sources of SoOps include global navigation satellite systems (GNSSs) and communications satellites. Re- ceivers can detect the re ection of SoOps from Earth's surface, thereby operating as bistatic radars. In particular, SoOp-R with GNSS sources and receivers is known as global navigation satellite system re ectometry (GNSS-R) [3]. For smooth surfaces where the scattering of SoOps is dominated by specular re ection, subkilometer resolution from low Earth orbit (LEO) can potentially be achieved [4]. Although such high resolution measurements are spatially sparse, they could potentially support spatiotemporal upscaling of SMOS and SMAP data to address the aforementioned resolution and revisit shortfalls [5{7]. Additionally, re ectometry is less sensitive to vegetation cover and more sensitive to soil moisture than conventional backscattering measurements from monostatic radar [8]. Furthermore, the associated signal-to-noise ratio (SNR) improvement opens the possibility of increasing the depth of sensing and building a soil moisture prole with depth. Finally, the use of SoOps can help to address data record continuity since many SoOp sources, such as GNSSs, are well funded and highly persistent. Various methods and results for retrieving soil moisture content from GNSS-R data have already been reported. One method uses a strictly noncoherent surface scattering model based on Hagfors' law [9]. Others assume that coherent scattering from a planar rough surface dominates [6, 10{12]. Here, the term noncoherent refers to scattering of the signal in a variety of directions with random uctuations, while the term coherent indicates re ection in the specular direction without uctuations. Thus, consensus on an appropriate model for the scattering of GNSS signals from land surfaces has not yet been reached. Furthermore, all GNSS-R-based retrieval methods thus far have had to rely heavily on soil moisture products from non-GNSS-R missions to calibrate each GNSS-R retrieval cell. The aim of the present research is to develop improved physics-based models for satellite remote sensing of land surfaces using GNSS-R. These models could be used to develop inversion algorithms to retrieve biogeophysical parameters, such as soil moisture, and reduce reliance on soil moisture products from non- GNSS-R missions. The models could also be used for sensitivity analyses and diagnostic purposes. We focus on GNSS signals since re ectometry missions using communication satellite signals are not yet operational as of 2019. Nevertheless, we anticipate that many of our results could be extended to SoOp-R using communication satellite signals. The remainder of this document is organized as follows. Chapter 2 provides a review of related work. In Chapter 3, we investigate the theoretical accuracy of a time series method for soil moisture retrieval that exploits the incidence angle diversity inherent in GNSS-R data. The method assumes coherent scattering of 2 electromagnetic waves in the presence of instrument calibration errors, vegetation cover uncertainty, and surface roughness uncertainty. In Chapter 4, we add a noncoherent scattering model based on the small perturbation method (SPM) and we attempt to validate the models with satellite data. In Chapter 5, motivated by deciencies in the previous chapter's models, we develop a new model for GNSS-R measurements over land surfaces with topography. The model makes use of a digital elevation model (DEM) and represents an extension of a well-known existing model for the scattering of GNSS signals from a rough sea surface. 3 Chapter 2 Review of Related Work Sections 2.1, 2.2, and 2.3 provide general background on soil moisture observation. They cover remote sensing principles, operational missions, and related non- re ectometry missions, respectively. Section 2.4 describes re ectometry from GNSSs, and Section 2.5 introduces re ectometry from communications satellites. Section 2.6 gives a background on GNSS measurement models. 2.1 Remote Sensing Principles The dielectric constant of soil is sensitive to moisture content in the microwave regime. Both passive and active microwave sensors can observe this sensitivity. Radiometers are passive instruments that measure electromagnetic radiation power. They can detect thermal emissions from the soil and radiation from the atmosphere scattered from the soil. Both the thermal emission and the scattering processes are sensitive to soil moisture through the soil dielectric constant. Radars, including scatterometers, imagers, and re ectometers, are active sensors that consist of a transmitter and a receiver. When microwaves emitted by the transmitter impinge on a soil, some of the power is transmitted into the soil and absorbed and some is scattered away from the soil and detected by the receiver. This scattering process also is sensitive to soil moisture through the soil dielectric constant [13]. 2.2 Missions with Operational Soil Moisture Products The 2017-2027 Decadal Survey for Earth Science and Applications from Space (ESAS 2017) identies soil moisture as a targeted observable supporting multiple science objectives [1]. Additionally, the implementation plan of the Global Climate Observing System (GCOS) by the World Meteorological Organization 4 (WMO) of the United Nations (UN) includes soil moisture among the essential climate variables (ECVs) [14]. Both ESAS 2017 and the GCOS plan identify the following four missions and sensors as supporting soil moisture observation: 1. SMAP 2. SMOS 3. Advanced SCATterometer (ASCAT) 4. Advanced Microwave Scanning Radiometer 2 (AMSR2) Of these, ESAS 2017 identies the NASA SMAP mission [15] and the ESA SMOS mission [16,17] as the most capable. Each of SMAP and SMOS ies a passive microwave L-band radiometer designed for soil moisture observation. Radiometers have high sensitivity but low spatial resolution. L band has the advantage of being able to operate in darkness and penetrate cloud cover and up to moderate levels of vegetation cover. L band is also legally protected from radio frequency interference (RFI) for the purposes of radio astronomy and passive Earth sensing [15,17]. The current resolution of the SMAP and SMOS radiometers is about 40 km with revisit times of 2 to 3 days. The SMAP mission also includes an active radar intended to provide a high resolution capability to complement the radiometry. However, the SMAP radar failed in July 2015 due to an anomaly with the power supply for the radar's high power amplier, leaving the science mission to continue with the radiometer alone. Consequently, the SMAP mission has been unable to achieve its anticipated combined active-passive resolution of 1 to 10 km [1]. The SMOS mission, which launched November 2, 2009, is already well beyond its planned duration of 3 years, while the SMAP mission, which launched January 31, 2015, is expected to continue until at least 2022. No enhanced replacement missions have been planned yet [1]. In addition to SMOS and SMAP, two other missions currently provide oper- ational soil moisture products. The ASCAT on the Meteorological Operational (MetOp) satellites of the European Organization for the Exploitation of Meteoro- logical Satellites (EUMETSAT) operates in C band and in vertical polarization with nominal resolution of 50 km and with research resolution of up to 25 km. With three MetOps currently operational, launched in 2006, 2012, and 2018, respectively, revisit time is sub-daily. Although ASCAT was not designed with soil moisture in mind and its frequency in C band was initially deemed less optimal for soil moisture than L band, various soil moisture applications have been demonstrated successfully [18{21]. The AMSR2 on the Japanese Aerospace Exploration Agency (JAXA) Global Change Observation Mission-Water 1 (GCOM-W1) satellite launched in May 2012. AMSR2 is a passive radiometer with dual polarization and multiple frequency bands. The C- and X-bands are used for soil moisture retrieval with resolution of 50 km and revisit times of 1 to 2 days [22{24]. The predecessor of AMSR2 was the Advanced Microwave Scanning Radiometer for EOS (AMSR-E) 5 on the Aqua satellite of the NASA Earth Observing System (EOS) program. AMSR-E is noted here only because it is referenced in the GCOS plan instead of its successor. 2.3 Soil Moisture from Other Non-Re ectometry Missions In addition to the operational soil moisture products provided by the four missions identied in the previous section, which include three radiometers and one scatterometer, soil moisture can also be retrieved from other types of sensors or from combinations of sensors. In the microwave regime, other sensors include synthetic aperture radar (SAR) and re ectometers (to be discussed in the next section). The advantage of SAR is high resolution, but the disadvantage is that revisit times tend to be too infrequent for most soil moisture applications. The ESA Sentinel-1 constellation [25], however, is the rst SAR mission with a 6-day repeat cycle, providing suciently frequent revisits for eective data assimilation. Sentinel-1 currently consists of two satellites carrying C-band SAR. In a recent study, joint assimilation of Sentinel-1 and SMAP observations produced better surface soil moisture estimates than either of the two sensors separately [26]. The ESA Climate Change Initiative (CCI) has assembled nearly four decades of global soil moisture records. The CCI project blends data from radiometers and scatterometers to generate a combined soil moisture product. The project seeks to incorporate SAR and thermal infrared data in the future for improved spatial resolution [27]. Although not a satellite mission, the NASA Earth Venture 1 (EV-1) Airborne Microwave Observatory of Subcanopy and Subsurface (AirMOSS) investiga- tion [28] is worth highlighting since, unlike the other sensors covered thus far which observe only the near-surface layer of soil, the frequency of the AirMOSS sensor in P band is suciently low that it penetrates down into the root zone. This penetration, combined with the multiple polarizations of the AirMOSS SAR data, facilitates the retrieval of soil moisture not only at the surface but also as a depth prole into the root zone. The root zone soil moisture (RZSM) is important for prediction of carbon uxes due to respiration and photosynthesis [29]. The AirMOSS radar has also been used to support the NASA Arctic Boreal Vulnerability Experiment (ABoVE). In arctic soils containing permafrost, the active layer is the top layer that thaws during the summer and freezes again during the autumn. The active layer thickness (ALT) is of particular interest since, if it increases from year to year, large amounts of greenhouse gases from decaying organic material could be released into the atmosphere. A recent study suggests that a multifrequency approach consisting of P-band and L-band data together for ALT retrieval outperforms approaches with only one of the two frequencies. The study illustrates the value of incorporating multiple frequencies, multiple polarizations, and multiple seasons to enhance observability of multiple unknown geophysical parameters and to improve retrieval accuracy [30]. 6 2.4 GNSS Re ectometry The concept of GNSS-R was rst proposed by Hall and Cordey for multistatic scatterometry from the sea surface in 1988 [31]. Although they concluded that the concept was not feasible with Global Positioning System (GPS) sources, their conclusion would soon be challenged. In July of 1991, a group of engineers with Dassault Electronique detected anomalous behavior with a GPS receiver that they were testing with an Alpha Jet military aircraft over the Atlantic Ocean at the French Ranges near Bordeaux. Instead of the usual hemispheric GPS antenna, the engineers had elected to install a wrap-around antenna on the nose of the Alpha Jet, which is capable of maneuvers of up to 6 Gs, to allow the GPS receiver to continue tracking during sharp turns with large roll angles. The group expected no multipath problems over the ocean because their calculations based on the Rayleigh criterion for smooth surfaces showed that, due to the roughness of the sea surface, no specular re ection should be possible [32]. However, upon investigating the anomalous behavior, the engineers were surprised to nd evidence that the GPS receiver had indeed locked on to a sea multipath signal over the Atlantic Ocean when reacquiring a satellite track lost during a maneuver. They suspected that the receiver had continued to track the sea multipath as the aircraft made its way back to the coast. They also discovered that the navigation solution generated by the GPS receiver during the period of suspected multipath was well outside the uncertainty bounds of the solution provided by independent radar tracking at the same time [32]. To investigate this phenomenon, the engineers carried out experiments to collect raw intermediate frequency (IF) data over the sea. They also developed a basic computer simulation of the multipath scattering from the rough sea surface. In spite of certain discrepancies between their simulation and experimental results, analysis of the experimental raw IF data conrmed the engineers' suspicions. In particular, they found that the delay-locked loop (DLL) of the GPS receiver was easily able to track GPS signals re ected from the rough sea surface. Because the total path length of the sea multipath signal is greater than the length of the direct path, the navigation solution of a GPS receiver becomes inaccurate when the receiver locks on to multipath, thus explaining the anomalous behavior of the GPS receiver observed during the July 1991 test ight [32]. We note that the Dassault simulation results imply a dependency between the power level received in the sea multipath channel and the roughness of the sea surface caused by wind. As it turns out, this dependency would be the motivation for future GNSS-R satellite missions to observe wind speed over the sea as illustrated in Fig. 2.1. Additionally, the Dassault engineers perceived that the time delay of the multipath signal relative to the direct path signal was dependent on the height of the receiver above the sea surface. This principle would prompt future eorts to observe changes in sea surface heights using GNSS-R satellite altimetry, as proposed by Mart n-Neira in 1993 [33]. 7 11 | Pa g e Estimation of Ocean Surface Wind Speed Using the onboard DDMI, each CYGNSS microsatellite observatory will receive direct signals from GPS satellites, as well as signals reflected off the ocean surface. As depicted in the image below, the direct signal is transmitted from the orbiting GPS satellite and received by the single zenith-pointing antenna (i.e., top side), while the scattered GPS signal scattered off the ocean surface is received by the two nadir-pointing antennas (i.e., bottom side). The direct signal is used to pinpoint the location of the observatory, while the reflected signals respond to ocean surface roughness, from which the ocean surface wind speed is derived. The specular point is the location on the ocean surface where all of the scattering originates if the surface is perfectly smooth. With a roughened surface, the scattering originates from a diffuse region called the glistening zone that is centered on the specular point. Image credit: University of Michigan The power in the signals scattered by the ocean surface and received by the DDMI is used to create a series of Delay Doppler Maps, examples of which are shown below. The y-axis represents the time delay between the direct (from the GPS) and scattered (from the ocean) received signals, while the x-axis represents the shift in frequency between the direct and scattered received signals. The two axes are normalized with respect to the delay and Doppler shift at the specular point, the spot on the ocean surface where the scattered signal strength is largest. The wind speed is estimated from the DDM by relating the region of strongest scattering (the dark red region) to the ocean surface roughness. A smooth ocean surface will reflect a GPS signal directly up toward the CYGNSS observatory, producing a strong received signal. A roughened ocean will result in more diffuse scattering of the signal in all directions, resulting in a weaker received signal. Therefore, strong received signals represent a smooth ocean surface and calm wind conditions, while weak received signals Figure 2.1: Cyclone Global Navigation Satellite System (CYGNSS) observation of direct and scattered signals. Image credit: University of Michigan. In 1997, NASA conducted an experiment consisting of a series of aircraft ights using a specialized GPS receiver to verify the cross-correlation function over the Atlantic Ocean. The experiment demonstrated that GNSS signal re ections can sense ocean surface roughness and related wind speed [34,35]. In 2002, Jet Propulsion Laboratory (JPL) researchers published results in which they succeeded in detecting and characterizing a GPS signal in calibra- tion data collected over the Pacic Ocean near the Galapagos Islands by the Spaceborne Imaging Radar C-Band (SIR-C) on NASA Space Shuttle Endeavour in October 1994. This nd represents the rst ever spaceborne observation of GNSS signal re ections from Earth at steep incidence [36]. Also in 2002, NASA ew a specialized GPS receiver/re ectometer on the National Science Foundation (NSF)/National Center for Atmospheric Research (NCAR) C-130 Hercules aircraft [37] over cropland near Ames, Iowa, as part of the Soil Moisture Experiment 2002 (SMEX02). Results from this rst-ever controlled GNSS-R experiment for soil moisture remote sensing indicated that GNSS-R is indeed sensitive to variations in soil moisture [38,39]. In October 2003, the United Kingdom Disaster Monitoring Constellation (UK-DMC)-1 satellite carried the rst GNSS-R into space as an experimental payload developed through a partnership between Surrey Satellite Technology Limited (SSTL) and National Oceanography Centre (NOC) [40]. The UK-DMC data rst demonstrated the feasibility of spaceborne GNSS-R for global sensing of ocean, land, and ice [41{44]. For the ocean wind application in particular, the UK-DMC data provided convincing validation of the GNSS-R techniques [45]. 8 In July 2014, the United Kingdom TechDemoSat-1 (TDS-1) satellite carried experimental payloads into orbit, including the Space GNSS Receiver Remote Sensing Instrument (SGR-ReSI) to demonstrate GNSS-R technology. Funded by the United Kingdom Space Agency (UKSA) and InnovateUK and built by SSTL, the TDS-1 SGR-ReSI supports the GPS L1 band with potential to be reprogrammed for L2C, Galileo, and GLObal NAvigation Satellite System (GLONASS) bands. In 2015, TDS-1 SGR-ReSI acquired the rst-ever spaceborne GNSS-R observations of hurricane wind conditions. The observations went up to and through the eyewall without loss of data, providing wind proles consistent with the expected structures around the eye of tropical cyclones [46{ 48]. Additionally, TDS-1 observations over the Punjab region of India and Pakistan showed a re ected signal sensitivity of 7 dB to temporal changes in soil moisture [49]. Recently, Galileo E1 and GPS L2C signals were recovered, making TDS-1 the rst multi-constellation GNSS-R [50]. In December 2018, routine data collection from the SGR-ReSI on TDS-1 ceased in preparation for deorbit. In December 2016, NASA launched the eight microsatellites comprising the Earth Venture Mission 1 (EVM-1) Cyclone Global Navigation Satellite System (CYGNSS), as depicted in Fig. 2.2. With development led by the University of Michigan and the Southwest Research Institute (SwRI), the CYGNSS mission measures wind speeds within tropical cyclones using sea-re ected GPS signals collected by an updated version of the SSTL SGR-ReSI, which is also referred to as the delay-Doppler mapping instrument (DDMI), carried as payload on each of the eight satellites. The DDMI operates at only one frequency, GPS L1, and one re ectometry polarization, left-hand circular polarization (LHCP). Unlike all other missions described above, the CYGNSS orbit is not polar but rather inclined at an angle of 35° relative to the equatorial plane. Thus, CYGNSS coverage is limited to, and optimized for, tropical regions [51,52]. Observations within the CYGNSS coverage zone have the potential to ll the current spatial and temporal gap in existing soil moisture products [6,7]. The rst version of a CYGNSS soil moisture product was released in May 2019 [10,11]. A number of additional GNSS-R-equipped satellite missions are being planned, including a GNSS-R payload for Spire radio occultation (RO) nanosatellites [53], the Flexible Microwave Payloads (FMPLs) of the Universidad Polit ecnica Catalunya (UPC) NanoSat Lab [54,55], the ESA GNSS rE ectometry, Radio Occultation, and Scatterometry on board the International Space Station (GEROS-ISS) [56], and the ESA Passive RE ecTometry and DosimeTrY (PRETTY) mission [50], and a possible GNSS-R payload for the LEO and medium Earth orbit (MEO) segments of the Kepler system proposed by the Deutsches Zentrum f ur Luft- und Raumfahrt (DLR) [57]. Furthermore, GNSSs themselves are proliferating, providing additional op- tions for GNSS-R sources. Besides the U.S. GPS, the ESA Galileo, and the Russian GLONASS GNSSs referenced above, other current and future GNSSs include the Indian GPS Aided Geo Augmented Navigation (GAGAN) system, the Indian Regional Navigation Satellite System (IRNSS), the Japanese Quasi- 9 7 | Pa g e First Satellite Signal Acquisition • First contact: Approximately three hours after microsatellite separation • Data downlink services will be provided by the SSC Space USA, Inc. Spacecraft Orbit • The CYGNSS Microsatellite will employ a 510km circular low-Earth orbit (LEO) at an inclination of 35˚ with the equator. Launch Site • Cape Canaveral Air Force Station, Cape Canaveral, FL, USA Launch Operations • NASA Launch Services Program, Kennedy Space Center, Merritt Island, FL, USA Mission Management • The Principal Investigator is based at the University of Michigan and has teamed up with the Southwest Research Institute for project management and implementation. Mission Operations Center • Southwest Research Institute, Boulder, CO, USA Science Operations Center • University of Michigan, Ann Arbor, MI, USA NASA Investment • $157 million Once the deployment module on the Orbital ATK Pegasus XL launch vehicle reaches its final altitude of approximately 510 kilometers above the Earth’s surface, the eight CYGNSS microsatellite observatories will be deployed in opposing pairs over a two minute period. (Image credit: University of Michigan) Figure 2.2: Deployment of CYGNSS microsatellites. Image credit: University of Michigan. Zenith Satellite System (QZSS), the Chinese BeiDou Navigation Satellite System (BDS), and the DLR-proposed Kepler system [4,57,58]. Thus, Kepler is proposed to function as both GNSS source and receiver. 2.5 Re ectometry with Communication Satellite Signals In July 2010, a NASA-funded airborne experiment over the Chesapeake Bay demonstrated that the ecient encryption and compression of communication signals enables them to be processed like range-coded navigation signals. In the experiment, an upward-looking antenna received Sirius XM Satellite Radio transmissions with LHCP in S band, while a downward-looking antenna received the corresponding re ections from the sea surface. The cross-correlation function of the two signals was sensitive to the distribution of sea surface slopes, allowing the retrieval of wind speed. This demonstration greatly increased the scope of op- portunistic bistatic radar since communication satellite signals have much higher power and a wider diversity of frequencies than their GNSS counterparts [59]. Additional low-altitude experiments with Ku-band DirecTV signals [60] and XM radio signals [61] conrmed the suitability of these SoOps for remote sensing. During the winter of 2015-2016, a proof-of-concept experiment conducted at the Fraser Experimental Forest (FEF) Headquarters, located in the central Rocky Mountains west of Denver, demonstrated the feasibility of observing snow water equivalent (SWE) using P-band SoOps. Tower-based antennas received direct and re ected U.S. Navy Ultra high frequency (UHF) Follow On (UFO) signals transmitted with right-hand circular polarization (RHCP) in the frequency band 10 from 240 to 270 MHz. This demonstration's success means that the concept can be extended to airborne and spaceborne platforms [62]. We note that UFO is being replaced by the Mobile Users Objective System (MUOS), which also operates at UHF [63]. In 2018, the In-space Validation of Earth Science Technologies (InVEST) program element of the NASA Earth Science Technology Oce (ESTO) awarded the SigNals of Opportunity P-band Investigation (SNoOPI) proposed by Purdue University under Research Opportunities in Space and Earth Science (ROSES)- 2017. The SNoOPI CubeSat will be the rst on-orbit demonstration of P-band SoOps. The purpose of the SNoOPI instrument is to measure RZSM and SWE. Also in 2018, JPL researchers presented at the Committee On SPAce Research (COSPAR) conference a multi-SoOp CubeSat concept that would operate in L, K, and Ka bands [64]. In addition to the above-named SoOps of DirecTV in K band, Sirius XM Satellite Radio in S band, various GNSSs in L band, and UFO/MUOS in UHF P band, another potential SoOp is the Orbcomm satellites in very high frequency (VHF) P band at 137:5 MHz [65]. Thus, in all, SoOps suitable for remote sensing span nearly two orders of magnitude of the microwave spectrum. Spatial coverage varies from continental in the case of satellite radio and television to global in the case of MUOS [66]. 2.6 Theoretical Models for Rough Surface Scattering The rst fairly complete treatment of wave scattering from surfaces with large- scale random roughness using the Kirchho approximation (KA) was provided by Isakovich in a 1952 Russian-language paper [67]. English-language descrip- tions of the resulting body of work became available in Rytov, Kravtsov, and Tatarskii [68], Bass and Fuks [69], and Beckmann and Spizzichino [70]. In 1996, NASA researchers Katzberg and Garrison applied this method to derive the power structure of a GPS signal that is re ected from the ocean and processed by cross-correlation with the reference GPS code [34,71]. In 2000, National Oceanic and Atmospheric Administration (NOAA) re- searchers Zavorotny and Voronovich extended the NASA model to include correlation in both time delay and Doppler frequency [72]. The two-dimensional representation of these correlations is known as a delay-Doppler map (DDM). This model, which we shall refer to as ZV2000, represents a bistatic radar equation for SoOps in the geometrical optics (GO) limit of the KA. As such, ZV2000 is valid only for scattering in the strong diuse regime where no coherent component of scattering exists. ZV2000 is the model used by the CYGNSS mission for ocean applications [73]. In Chapter 5, we extend ZV2000 to the modeling of DDMs over land surfaces with topography. 11 In 2017{2018, Voronovich and Zavorotny removed the KA and generalized ZV2000 beyond the strong diuse scattering regime to include weak diuse scattering. The resulting bistatic radar equation, hereafter VZ2018, is valid for any surface whose roughness statistics are spatially homogeneous. VZ2018 accounts for partially coherent scattering from relatively at areas, such as sea ice, lakes, and some land [74,75]. In practice, the coherent component of the scattered signal is frequently modeled in the KA assuming a Gaussian distribution of heights [69, eqs. (24.10) and (20.6); 76, eq. (10); 77; 78, eq. (1); 79, Section 21-10; 80, eq. (12D.10)]. This model forms the basis for the University Corporation for Atmospheric Research (UCAR)/University of Colorado (CU) soil moisture retrieval product [10,11] and other retrievals [6,12] from CYGNSS data. The model is also used in GNSS-R simulators for land applications [81{83]. In Chapter 3, we apply this model to a soil moisture retrieval sensitivity analysis. Although several studies have concluded that a coherent scattering model is valid for airborne GNSS-R data [39, Fig. 4; 4, Fig. 5], we are unaware of any such results for spaceborne data, such as from CYGNSS. Due to the higher altitude of the receiver, spaceborne data have a larger footprint than airborne data. This generally leads to larger surface height variation over the footprint, making coherent scattering less likely. Validity of the coherent model for CYGNSS data over a soil moisture validation site is examined in Chapter 4. Indeed, Ohio State University (OSU) researchers contend that coherent re ections in spaceborne data arise not from land surfaces but rather primarily from inland water bodies [9]. For this reason and in contrast with other soil moisture retrievals, the OSU method is based on Hagfors' law, which is the high frequency limit of the noncoherent physical optics (PO) approximation for an exponentially correlated Gaussian surface [78,84,85]. Since PO is based on the KA, both the OSU method and ZV2000 belong to the category of noncoherent KA methods. Another analytic method used to study wave scattering from random rough surfaces is the SPM [80, 86, 87]. Unlike the KA, which is valid for gently undulating surfaces, the SPM can account for surface roughness on a scale of less than a wavelength. Scattering due to small-scale roughness is sometimes also called Bragg scattering. The SPM requires surface height variation to be small relative to the wavelength. Thus, in contrast to the KA, which is valid in the high frequency limit, the SPM holds in the low frequency limit. Although SPM has been used as part of soil moisture retrievals from traditional backscattering radar [29], we are not aware of any application to GNSS-R at this time. We investigation the application of SPM to GNSS-R in Chapter 4. The advanced integral equation method (AIEM) is a numerical method devel- oped to address the limitations of the KA and the SPM [88]. The AIEM forms the basis for the noncoherent scattering model of the Soil And VEgetation Re ec- tion Simulator (SAVERS) [81]. This simulation also has a coherent component based on the Gaussian KA model of [77]. The coherent component has recently 12 been extended to surfaces with topography, although the topography-specic details of the model are still forthcoming [89,90]. The drawback of AIEM is its computational cost, making it dicult to use for real-time retrievals. Another method developed to address the limitations of the KA and the SPM is the small slope approximation (SSA) [91, 92]. This method has been used to study polarization ratios and features in bistatic scattering from rough surfaces with potential for GNSS-R application [93,94]. The method has also recently been applied to monostatic and bistatic scattering of anisotropic rough soil, such as found in plowed agricultural elds [95]. The method can be used in the framework of VZ2018 to simulate GNSS-R measurements [75], although we are not aware of any such full GNSS-R simulation results for land applications at this time. 13 Chapter 3 Sensitivity Study 3.1 Introduction In this chapter, we investigate the theoretical accuracy of a time series method for surface soil moisture retrieval from GNSS-R measurements, such as from the CYGNSS mission, in the presence of various error sources. Unlike the existing retrievals from CYGNSS data described in Chapter 2, we consider a physics-based method without relying on SMAP or other satellite soil moisture products to calibrate each resolution cell. Thus, we must account for the eects of instrument calibration errors, vegetation cover uncertainty, and surface roughness uncertainty. We focus on exploiting incidence angle diversity within the time series since the eect of surface roughness on scattered power is relatively sensitive to incidence angle and also since CYGNSS lacks other modes of measurement diversication, such as frequency and polarization. In particular, we consider a given resolution cell on the ground through which multiple specular points pass within a given time interval to provide multi-incidence angle re ectometry data. The diversity of incidence angles in the resulting time series is a consequence of the inherent randomness of the corresponding bistatic geometries. At the time this analysis was initiated in early 2018, the CYGNSS calibration team had been working through various calibration problems. One of these problems was inter-observatory statistical biases resulting from a digital-to-analog scaling error, as later documented in [96]. Thus, it was of particular interest to investigate the impact of calibration errors on retrieval performance. Calibration errors are still of interest as track-wise biases continue to be reported [97,98]. We are unaware of any related sensitivity analyses for time series bistatic methods in the literature. The closest reference appears to be [99], which computes the theoretical variance of soil moisture retrieval error for a bare surface of unknown or uncertain roughness using various combinations of monostatic and bistatic airborne radar measurements in various polarizations. We note that (7) of [99] is equivalent to (3.11) below by the Woodbury matrix identity. 14 Various other bistatic sensitivity studies have also been reported [100{103]. An error model for estimates of another biogeophysical parameter from monostatic (i.e., backscattering) radar polarimetry is presented in [104]. However, none of these works has attempted to analyze any time series retrieval methods. This chapter is organized as follows. Section 3.2 denes an observation model and a retrieval cost function. Our goal is not necessarily to develop complex high-delity models but rather only to capture the primary features relevant to sensitivity analysis. From these, we derive a formula for retrieval performance as a function of calibration uncertainty, time series parameters, and prior statistics. Section 3.3 illustrates relationships among various parameters with numerical results. Section 3.4 discusses the results, and Section 3.5 describes limitations of the models and suggests directions for future work. Much of the material in the present chapter is drawn from the author's publication [105]. 3.2 Statistical Modeling 3.2.1 Observation Model Following [10, eq. (S2)], we assume a coherent model for power received from the GNSS source P g = P t G sp t 4(R sp rs +R sp st ) 2 G sp r 2 4 (3.1) where P t is transmit power, G sp t is transmit antenna gain in the direction of the specular point, R sp st is the distance from the transmitter to the specular point on the surface, R sp rs is the distance from the specular point on the surface to the receiver, G sp r is receive antenna gain in the direction of the specular point, and is the polarization-dependent surface power re ectivity. This is the Friis transmission formula for a specular re ection [13, Section 3-3]. We further assume that P g can be obtained from the GNSS-R measurement and that system parametersP t ,G sp t , andG sp r and rangesR sp st andR sp rs are known so that we can solve for . In the case of CYGNSS, can be conveniently obtained from the level 1 (L1) variable brcs, which contains the DDM of bistatic radar cross sections (BRCSs), as follows = (R sp rs +R sp st ) 2 4(R sp rs ) 2 (R sp st ) 2 sp (3.2) where sp is the BRCS in the specular bin of brcs. This result follows directly from (3.1) above and from (7.2) and (7.3) of the CYGNSS Handbook [73]. Thus, we take to be the observation. Following the approach in [106], which in eect combines (11.2) and (11.11) of [13], we model the observation as a function of biogeophysical parameters by (;m v ;ks;) =A vwc (;)A rough (;ks) f (;m v ) (3.3) 15 where A vwc (;) =e (2 sec) (3.4) A rough (;ks) =e (2 cos) 2 (ks) 2 (3.5) =bw c (3.6) and is incidence angle,m v is volumetric soil moisture,k is wavenumber,s is the standard deviation of surface roughness, ks is the surface roughness coecient, is optical depth, A vwc is power attenuation due to vegetation water content (VWC), A rough is attenuation of coherent power due to surface roughness under the KA for Gaussian distribution of heights (references provided in Section 2.6), f is the polarization-dependent Fresnel power re ectivity, and w c is VWC. The vegetation parameter b depends on cover type/structure, wavelength, and polarization. In the following, we takeb = 0:1 based on results from [107, Table 1] near the CYGNSS wavelength of 19 cm. Dening biogeophysical parameter vector x = [m v ks ] T (3.7) and assuming that x is constant over the spatiotemporal extent of the time series and expressing in units of decibels f(;x) = 10 log 10 (;x) = 20 sec ln 10 40 cos 2 ln 10 (ks) 2 + 10 ln 10 ln f (3.8) then we can write the time series of observed surface re ectivities in vector form as y =f(x) + cal (3.9) where y i =f( i ;x) + cal i for observation i = 1;:::;N. Here, we assume that observation noise cal in decibels comes from instrument calibration error and behaves as multiplicative (speckle-like) noise when (3.9) is expressed in linear units. Furthermore, if the multiplicative calibration noise is Gaussian and suciently small (about unity), then the noise remains approximately Gaussian when converted from decibels to linear units, i.e., cal is also Gaussian [108, appendix]. 3.2.2 Cost Function and Error Covariance Matrix We assume that the cost function for the retrieval has the following quadratic form J(x 0 ) = [yf(x 0 )] T C 1 cal [yf(x 0 )] + (x prior x 0 ) T C 1 prior (x prior x 0 ) (3.10) 16 where the prime inx 0 serves to distinguish it from truthx,C cal is the covariance matrix of cal , x prior is the prior mean of x, and C prior is the corresponding covariance matrix. We also assume that cal andx prior x are zero mean and uncorrelated with each other. As an aside, (3.10) remains a valid cost function even if the quantitiesC 1 cal , C 1 prior , andx prior have no statistical meaning attached to them. The following results can be extended to this more general class of retrievals, which includes those based on simplied models [109]. However, for this study we make the statistical assumptions stated in the previous paragraph. To rst order, the error covariance matrix (ECM) for the retrieval is given by C ^ x = E[(^ xx)(^ xx) T ] = (F T C 1 cal F +C 1 prior ) 1 (3.11) where the retrieval ^ x is the minimizer of (3.10) andF is the Jacobian off(x). The rst-order approximation in (3.11) is valid when ^ xx is suciently small relative to higher-order partial derivatives off(x). This, in turn, depends on cal andx prior x being suciently small [109]. If cal andx prior x are jointly Gaussian, then ^ x is the maximum a posteriori probability (MAP) estimator. It is also the minimum mean square error (MMSE) estimator to rst order, and therefore the MMSE is given by the trace of (3.11) to rst order [109]. If no prior information is available, i.e., ifC 1 prior = 0 in (3.10), and if cal is Gaussian, then ^ x is the maximum likelihood (ML) estimator. In this case, ^ x is also the minimum variance unbiased (MVU) estimator to rst order [109]. If prior information is available for some but not all parameters ofx and if these parameters and cal are jointly Gaussian, then ^ x is the mixed ML-MAP estimator [109]. Even if the random variables involved are non-Gaussian, the cost function given by (3.10) may still provide a good retrieval. In any case, (3.11) always holds to rst order under the assumptions of the previous sections regardless of whether the random variables involved are Gaussian [109]. 3.2.3 Partial Derivatives From (3.8), the partial derivatives of the observation model are given by @f @m v = 10 ln 10 1 f @ f @m v (3.12) @f @(ks) = 80 cos 2 ln 10 ks (3.13) @f @ = 20 sec ln 10 : (3.14) Calculation of@ f =@m v in (3.12) requires a soil dielectric model, such as [110,111]. Details are provided in Appendix A. 17 3.2.4 Inner Product Representation We now assume that calibration errors are uncorrelated and homoscedastic with variance 2 cal and that prior uncertainties are uncorrelated so thatD 1 prior =C 1 prior is diagonal. Then, (3.11) becomes C ^ x = 2 cal N ( 1 N F T F + 2 cal N D 1 prior ) 1 : (3.15) Observing that 1 N F T F = 1 N N X i=1 [rf( i )][rf( i )] T (3.16) has the form of an average, we generalize (3.16) by dening a matrix of inner products P = Z [rf()][rf()] T p ()d (3.17) where the elements ofP are given by p lm =hf l ;f m i = Z f l ()f m ()p ()d (3.18) and p () is the probability density function of incidence angles and we have dened the shorthand notation f n () = @f(;x) @x n : (3.19) Here, as in the following, indices l, m, and n will refer to the elements ofx. If we takep () = P N i=1 ( i )=N, then we nd thatP =F T F=N. Thus, we can write (3.15) in the more general form C ^ x = 2 cal N (P + 2 cal N D 1 prior ) 1 (3.20) whereP is given by (3.17) and (3.18). This has the advantage of expressing the ECM in a manner not tied to specic samples of i . 3.2.5 Soil Moisture Retrieval Performance We proceed to derive an expression for the (1; 1)-element ofC ^ x , which represents soil moisture retrieval uncertainty. SinceP is symmetric non-negative denite, we nd the following factorization to be useful P =D 1=2 f RD 1=2 f (3.21) where D f is the diagonal of P with elementshf n ;f n i =kf n k 2 and R is the correlation matrix ofP with elements given by lm = hf l ;f m i p hf l ;f l ihf m ;f m i = hf l ;f m i kf l kkf m k : (3.22) 18 Letting the expression inside the parentheses of (3.20) be e P =P + 2 cal N D 1 prior (3.23) it can be veried by matrix algebra that e P = e D 1=2 f e R e D 1=2 f (3.24) where the elements of diagonal matrix e D f are given by a 2 n kf n k 2 , the diagonal elements of e R are unity, and the o-diagonal elements of e R are given by ~ lm =a l a m lm (3.25) and the attenuation-like quantity is dened by a n = 1 + 2 cal Nkf n k 2 2 n 1=2 (3.26) and n is the uncertainty of the n th prior. We note 0 a n 1. Thus, the priors have the eect of making e P larger and more positive denite (i.e., more diagonal) relative toP . For a matrixM, letM(l;m) denote the submatrix formed by deleting row l and column m fromM. Then, we can factor e P (1; 1) = e D 1=2 f (1; 1) e R(1; 1) e D 1=2 f (1; 1): (3.27) From Cramer's rule and from (3.21) and (3.27), the (1; 1)-element of e P 1 is given by e P 1 1;1 = det h e P (1; 1) i det e P = det h e D 1=2 f (1; 1) i det h e R(1; 1) i det h e D 1=2 f (1; 1) i det e D 1=2 f det e R det e D 1=2 f =a 2 1 kf 1 k 2 det h e R(1; 1) i det e R : (3.28) The desired result follows from (3.20), (3.23), and (3.28) ^ mv = cal p N a mv @f @m v 1 v u u u t det h e R(1; 1) i det e R : (3.29) The determinant factor in (3.29) is given explicitly by s 1 ~ 2 23 1 ~ 2 23 ~ 2 12 ~ 2 13 + 2~ 23 ~ 12 ~ 13 : (3.30) 19 Table 3.1: Soil Parameters Frequency 1:575 42 GHz Water temperature 20 C Salinity of water 4 g salt kg 1 water Bulk density 1:55 g cm 3 Particle density 2:66 g cm 3 Mass fraction of sand 40 % Mass fraction of clay 50 % Mass fraction of silt 10 % Water volume fraction m v 20 % From this, we see that the determinant factor is minimized to unity when ~ 12 = ~ 13 = 0, i.e., when one or more of the following hold @f @m v ; @f @(ks) = @f @m v ; @f @ = 0 (3.31) a ks =a = 0 (3.32) a mv = 0: (3.33) From (3.26), we see that priors are eective when n cal p Nkf n k (3.34) since condition (3.34) drives a n to zero. If instead n cal p Nkf n k (3.35) holds, then the prior has no benet to retrieval performance (3.29) since condition (3.35) drives a n to unity. 3.3 Numerical Results 3.3.1 Parameters Computation off(;x) in (3.8) and its partials in Section 3.2.3 and in Appendix A requires specication of transmit/receive polarization and soil dielectric constant. Following the CYGNSS design, we assume the transmitted signal has RHCP and the receive antenna has LHCP to accommodate polarization change on re ection. We use the soil dielectric model of Peplinski-Ulaby-Dobson [110,111] with soil parameters given in Table 3.1, which is based on [112, Table IV], except that we have increased the value of m v to a more moderate level. 20 10 20 30 40 50 60 70 Incidence Angle (deg) 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 Normalized Sensitivity m v -ks - Figure 3.1: Normalized sensitivities f n ()=kf n k = (@f=@x n )=k@f=@x n k as a function of incidence angle for x n =m v ;ks;. As indicated in the legend, the sign of the sensitivity for ks and has been reversed for plotting convenience. Table 3.2: Sensitivity Norms and Correlations (deg) @f @mv @f @(ks) @f @ mv;ks mv; ks; 10{70 12:5 2:9 13:7 0:893 0:959 0:735 10{40 12:2 3:7 9:8 0:990 0:997 0:977 35{55 12:3 2:7 11:8 0:968 0:992 0:929 40{70 12:7 1:6 16:6 0:912 0:979 0:813 3.3.2 Sensitivities Normalized sensitivities are shown in Fig. 3.1 for 0° 70°. The corresponding norms and correlations (3.22) of the sensitivity functions over various ranges of are given in Table 3.2. We assume has uniform distribution over its corresponding range to compute the underlying inner products (3.18). We note the values ofk@f=@m v k andk@f=@(ks)k depend on the assumed values of m v and ks. Here, m v = 0:20 by Table 3.1 and we take ks = 0:13 based on SMAP ancillary data [106]. The determinant factor (3.30) is calculated in Table 3.3 for the -intervals and-values from Table 3.2. When a parameter indexed byn is listed as \known" in Table 3.3, then n = 0 and thusa n = 0 by (3.26) and the corresponding values of ~ in (3.30) are zero by (3.25). This is equivalent to removing the parameter fromx. 21 Table 3.3: Determinant Factor Known? Incidence Angle Range (degrees) m v ks 10{70 10{40 25{55 40{70 N N N 16:0 236:3 75:9 26:5 N Y N 3:5 13:7 7:9 4:9 N N Y 2:2 7:0 4:0 2:4 N Y Y 1:0 1:0 1:0 1:0 3.3.3 Calibration Accuracy Requirement Solving (3.29) for cal , we obtain the calibration accuracy required to satisfy a specied retrieval accuracy ^ mv as a function of the time-series parameters and priors cal p N ^ mv a mv @f @m v 8 < : det h e R(1; 1) i det e R 9 = ; 1 =2 : (3.36) Fig. 3.2 shows cal as a function of prior surface roughness and prior VWC uncertainties, where we x ^ mv = 4%, N = 4, 0° 70°, and prior mv =1. Here, the prior uncertainties are expressed in decibels of receive power for nadir incidence, i.e., byj@f n ( = 0)j n for n = ks;. This has the advantage of makingkf ks k ks in (3.26) independent of ks. Otherwise, Fig. 3.2 would depend on the value of ks. From (3.13), it follows that a value on the y-axis can be converted to units of ks for an assumed true value ofks by dividing by 34:7ks. Similarly, from (3.14), it follows that a value on the x-axis can be converted to units of by multiplying by 0:115 13. Furthermore, from (3.6), an x-axis value can be converted to VWC in units of kg m 2 by multiplying by 1:1513. For the case where ks !1 in Fig. 3.2, i.e., when surface roughness is a priori unknown, the corresponding required calibration accuracy is shown as the blue line plotted in Fig. 3.3. Three other lines are also plotted in Fig. 3.3 for the three subintervals of incidence angle considered in Tables 3.2 and 3.3. Although these results are only for m v = 0:2, we have found in general that the calibration requirement becomes looser for drier soils and tighter for wetter soils due to m v -dependence of the sensitivity factork@f=@m v k in (3.29). 3.4 Discussion A model for performance of a time series regression for soil moisture retrieval is provided by (3.29). The model is a function of the number of observations in the time series, calibration uncertainty, prior statistics of biogeophysical parameters, and incidence angle distribution. As seen in Fig. 3.1, the normalized sensitivity function for moisture of a typical soil diers only a little from the sensitivity functions for the other two biogeophysical parameters. This is re ected in values 22 0.1 0.1 0.2 0.2 0.2 0.3 0.3 0.3 0.4 0.4 0.5 0.5 0.6 0.6 0.7 0.7 0.8 0.9 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Figure 3.2: Calibration uncertainty (dB) as a function of prior surface roughness uncertainty and prior VWC uncertainty that is required to achieve soil moisture retrieval uncertainty of 4% using 4 observations between 0° and 70° of incidence. 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 (deg) 10-70 10-40 25-55 40-70 Figure 3.3: Calibration uncertainty (dB) as a function of prior VWC uncertainty required to achieve to soil moisture retrieval uncertainty of 4% using 4 observa- tions for various incidence angle intervals when surface roughness is unknown. Multiply x-axis by 1:1513 to convert to kg m 2 . 23 near unity forj mv;ks j andj mv; j in Table 3.2. This, in turn, results in large values for the determinant factor in (3.29), as seen in the rst row of Table 3.3. Thus, incidence angle does a poor job of separating the biogeophysical parameters and, as a result, the regression problem is poorly conditioned. Prior information about the biogeophysical parameters is therefore of signicant interest. We note that, as expected, Table 3.3 indicates that a large range of incidence angles is better than a small range. The results also suggest that larger incidence angles provide more value than smaller angles. A formula for the calibration performance required to meet a specied retrieval accuracy is given by (3.36). The formula is a function of desired soil moisture retrieval accuracy, number of observations in the time series, prior statistics of surface roughness and VWC, and incidence angle distribution. The importance of prior information is highlighted in Fig. 3.2. In this case, calibration accuracy must be less than 0:1 dB to achieve the desired retrieval performance without any prior information. For comparison, the published calibration uncertainty CYGNSS L1 data is 0:39 dB [96]. Thus, the contour line at 0:4 dB approximately denes the region of prior uncertainties required to meet the desired retrieval performance. If we consider surface roughness to be unknown, then Fig. 3.3 suggests that, for the given parameters, soil moisture retrieval with a time series method is feasible for surfaces that are known to be bare or nearly bare. 3.5 Future Directions The validity of (3.29) and of the numerical results in Section 3.3 depends on the many assumptions stated above. The assumption that coherent power can be measured from DDM data is still controversial, as already noted in Chapter 1 and Section 2.6. The simple VWC model (3.4) and (3.6) may need to be developed further to account for -dependent vegetation structure as in [113]. Validity of the KA in (3.5) depends on properties of the rough surface. Although we assumed perfect knowledge of soil dielectric model and soil texture for computing @ f =@m v , our analysis could be extended to include soil texture uncertainty. We also assumed no soil moisture variation with depth. However, such variation may occur and could aect re ected power, e.g., similar to an antire ective coating. The assumption that biogeophysical parameters remain constant over the spatial-temporal extent of the time series needs to be investigated. In the case of CYGNSS, the specular region tends to be long and narrow due to the sensor's relatively long noncoherent integration time. Thus, consecutive measurements within a given cell on the ground may not overlap much [10, Fig. S1], giving rise to spatial noise. Furthermore, soil moisture may change during a time series window due to weather events. These variations could be accommodated, for example, in the framework of Kalman ltering, which allows parameters to vary stochastically in time and space. 24 The cost function (3.10) has desirable optimality properties when the random variables involved are Gaussian as discussed in Section 3.2.2. However, we know that the prior distributions of biogeophysical parameters cannot be perfectly Gaussian since their values are constrained to lie in closed or half-open intervals. Incorporation of these intervals as constraints in the retrieval algorithm may increase its performance. In this case, sensitivity model would need to be updated accordingly. Although we assumed prior estimates of ks and were uncorrelated, if these estimates come from another sensor, such as SMAP, then their uncertainties may be correlated and this would need to be accounted for. Additionally, the assumed model of calibration error as uncorrelated homoscedastic multiplicative noise needs to be examined more closely and developed if needed. In particular, measurement uncertainty may increase for large incidence angles due to increased receive antenna pattern characterization error and decreased SNR. While validity of the ECM (3.11) does not depend on Gaussianity, it does depend on a linear approximation whose validity should be explored. Also, while (3.11) assumes perfect knowledge of all statistics, it can be extended to account for model mismatch in (3.10). It is important to understand the eects of such mismatch since statistics are seldom precisely known. Numerical results in Section 3.3 assume that is uniformly distributed over a given interval, and we see these results are sensitive to interval size and location. Thus, the assumed uniform distribution should be updated to match the actual -distribution from on-orbit data. In the next chapter, we investigate sites for validation of the observation model dened in Section 3.2.1. We also consider the addition of a noncoherent scattering component. 25 Chapter 4 Models for Planar Rough Surfaces 4.1 Introduction In this chapter, the coherent rough surface scattering model of the previous chapter is extended by the addition of a noncoherent part. This approach is similar to SAVERS, which also has separate coherent and noncoherent models. However, unlike SAVERS which uses the AIEM for the noncoherent model, we investigate the use of the SPM to provide improved computational eciency. Background on SAVERS, AIEM, and SPM, including references, is provided in Section 2.6. The modeling of full DDMs based on the coherent scattering model of the previous chapter and the noncoherent part of the SPM is described in Section 4.2. These models are valid for rough surfaces that are planar on average. Our attempts to validate the models against CYGNSS data are reported in Section 4.3. Results are discussed in Section 4.4. The results provide motivation for the topography-based model of the next chapter. 4.2 Theoretical Modeling 4.2.1 Total Power The VZ2018 model, introduced in Section 2.6, decomposes power received in delay-Doppler bin (i;j) into a coherent and a noncoherent component P tot g (i;j) =P coh g (i;j) +P nc g (i;j): (4.1) The subscript g, included for consistency with [73], indicates that power from the GNSS source and rather than from thermal noise in the receiver or antenna. 26 4.2.2 Bistatic Radar Cross Section For CYGNSS data, the Level 1B (L1B) calibration converts the DDM from power received from the re ected GPS waveform in units of watts to surface BRCS in units of square meters. This conversion, given by (7.3) of [73], removes the eects of P t , G t , and G r and of range, atmospheric, and instrument loss factors. The resulting BRCS DDM is provided in the L1 variable brcs. We elect to use this variable as our observation since it benets from the entire CYGNSS calibration. Expressed in the notation of the previous chapter, the L1B-calibrated BRCS is given by (i;j) = (4) 3 (R sp rs ) 2 (R sp st ) 2 2 P t G sp t G sp r L ai (i;j) P L1A g (i;j) (4.2) where L ai (i;j) represents atmospheric losses to and from the surface and known losses introduced by the DDMI and where P L1A g (i;j) is the Level 1A (L1A)- calibrated power received from the GNSS source. But P L1A g (i;j) =L ai (i;j)P tot g (i;j) (4.3) since the total modeled power P tot g (i;j) from (4.1) does not include atmosphere and instrument losses. Thus, from (4.2) and (4.3), a model for BRCS is given by (i;j) = (4) 3 (R sp rs ) 2 (R sp st ) 2 2 P t G sp t G sp r P tot g (i;j): (4.4) 4.2.3 Woodward Ambiguity Function Since we are modeling full DDMs, we will need the point spread function of the coded waveform, known as the Woodward ambiguity function (WAF). Let the coded signal be given by s c (t) =a(t)e i2fct (4.5) wheref c is the carrier frequency,a(t) is the code, andjaj 2 = 1 in a time-averaged sense. In the case of the GPS coarse acquisition (C/A) code, a =1. We dene the WAF to be (;f;t 0 ) = 1 T i Z t0+Ti t0 a(t)a (t)e i2ft dt (4.6) where T i is the coherent integration time of the DDMI. We dene the mean square WAF to be j(;f)j 2 = lim T!1 1 2T Z T T j(;f;t 0 )j 2 dt 0 : (4.7) 27 For the GPS C/A code, which is used by the CYGNSS mission, a good approxi- mation is provided by j(;f)j 2 = 2 ()S 2 (f) (4.8) where () = maxf0; 1jj= c g (4.9) S(f) = sinc(T i f) (4.10) and c is the chip length. This approximation is good to a relative accuracy of 10 4 [75]. 4.2.4 Coherent Component Following VZ2018, we extend the coherent model (3.1) from a scalar power to a full DDM of powers by including a factor of the WAF to account for the point spread function of the waveform P coh g (i;j) = P t G sp t 4(R sp rs +R sp st ) 2 G sp r 2 4 D ( sp i ;f sp j ) 2 E : (4.11) Here, ( sp i ;f sp j ) is the delay and Doppler of of bin (i;j) relative to the delay and Doppler of the specular point, respectively. From (4.1), (4.4), and (4.11), the coherent part of the BRCS DDM is given by coh (i;j) = 4(R sp rs ) 2 (R sp st ) 2 (R sp rs +R sp st ) 2 D ( sp i ;f sp j ) 2 E : (4.12) 4.2.5 Noncoherent Component A model for the noncoherent component of received power is given by [73, eq. (7.1); 75, eq. (13)] P nc g (i;j) = P t 2 (4) 3 Z G t nc 0 G r D j( i ;f j )j 2 E (R rs ) 2 (R st ) 2 d~ r (4.13) where the integration proceeds over the plane of the rough scattering surface, ~ R s is position of the integration variable in this plane where we shall refer to ~ R s as the surface point, ~ R r is receiver position, ~ R t is transmitter position, R rs =j ~ R r ~ R s j is the distance from the receiver to the surface point, R st =j ~ R s ~ R t j is the distance from the surface point to the transmitter, G t is transmitter gain in the direction of the surface point, nc 0 is the noncoherent component of the normalized bistatic radar cross section (NBRCS) of the surface at ~ R s , G r is the receive antenna gain in the direction of the surface point, the delay and Doppler 28 frequency of bin (i;j) relative to the delay of Doppler frequency of ~ R s are dened by i = i (R rs +R st )=c (4.14) f j =f j f D (4.15) and i is the delay of bin i, f j is the Doppler frequency of bin j, the Doppler frequency of ~ R s is given by f D = ( ~ V t ^ R st ~ V r ^ R rs )= (4.16) and ~ V t is transmitter velocity, ~ V r is receiver velocity, and the unit vectors are given by ^ R rs = ~ R r ~ R s j ~ R r ~ R s j and ^ R st = ~ R s ~ R t j ~ R s ~ R t j : (4.17) The SPM formulation for nc 0 is a function of incidence angle i and scattering angles s and s , e.g., see [79, eqs. (21-67) and (21-68)]. Thus, to carry out (4.22), we must compute these angles as a function of ~ R s . Let ^ z be the up-pointing unit vector normal to the plane of the rough surface. We adopt the convention of [75,114] and dene the unit vector ^ y to be the intersection of the rough surface plane and the local incidence plane at ~ R s , with orientation away from the receiver and toward the transmitter. Unit vector ^ x completes the right-handed system. Thus, ^ x, ^ y, ^ z denes a local coordinate system at ~ R s . Then, the incidence and scattering angles are computed by i = arctan( ^ R st ^ y; ^ R st ^ z) (4.18) s = arctan( ^ R rs ^ y; ^ R rs ^ z) (4.19) s = arctan( ^ R rs ^ x; ^ R rs ^ y) (4.20) where the two-argument arctangent is used for numerical stability and is dened such that = arctan(a;b) satises tan =a=b. This is the same convention as atan2 of both Matlab and Python. Following [73], we assume that the support region of the integral in (4.13) is suciently small that several parameters can be approximated by their values at the specular point (sp) and factored out of the integral as follows P nc g (i;j) = 2 P t G sp t G sp r (4) 3 (R sp rs ) 2 (R sp st ) 2 Z nc 0 D j( i ;f j )j 2 E d~ r: (4.21) Thus, from (4.1), (4.4), and (4.21), we obtain the noncoherent component of the BRCS DDM nc (i;j) = Z nc 0 D j( i ;f j )j 2 E d~ r: (4.22) The integral (4.22) can be evaluated numerically using the above formulas if so desired. Alternatively, in the following we propose an approximation to simplify the integration. In particular, we observe that the WAF is much more 29 peaked than nc 0 as a function of the surface integration variable. For a given bin (i;j), the peaks of the WAF are determined by the solutions to i ( ~ R s ) = 0 (4.23) f j ( ~ R s ) = 0 (4.24) where i and f j are given by (4.14) and (4.15), respectively. A second-order Taylor series expansion of (4.14) and (4.15) in (4.23) and (4.24) leads to equations for an ellipse and a hyperbola, respectively. In general, the ellipse and hyperbola will intersect at two dierent points, which we will index by n = 1; 2. Thus, for each bin (i;j), there are two regions R i;j;n , n2f1; 2g, that contribute to the integral (4.22). For bins (i;j) for which no solution to (4.23) and (4.24) exists, the peaks of the WAF can be found approximately by minimizing the distance J( ~ R s ) = 8 < : " i ( ~ R s ) # 2 + " f j ( ~ R s ) f # 2 9 = ; 1 =2 (4.25) where and f D are delay and Doppler bin spacing, respectively. Since the angular extent of each region is relatively small, we assume nc 0 is constant in each of these two regions and that each pair of regions is disjoint R i;j;1 \R i;j;2 =?. We can then write nc (i;j) = X n=1;2 Z Ri;j;n nc 0 ( i ; s ; s ) D j( i ;f j )j 2 E d~ r = X n=1;2 nc 0 [ i (i;j;n); s (i;j;n); s (i;j;n)] A(i;j;n) (4.26) where the arguments (i;j;n) indicate evaluation at the corresponding point of intersection and where we have dened eective surface scattering area by A(i;j;n) = Z Ri;j;n D j( i ;f j )j 2 E d~ r (4.27) for n = 1; 2. We now make the additional approximation that the two areas are equal A(i;j; 1) = A(i;j; 2) (4.28) and we dene A(i;j) = A(i;j; 1) + A(i;j; 2) (4.29) which is consistent with the denition of the eective surface scattering area given by [73, eq. (7.5)]. Thus, from (4.26)-(4.29), we obtain the following approximate model for the noncoherent component of the BRCS DDM nc (i;j) = nc 0 (i;j) A(i;j) (4.30) 30 where nc 0 (i;j) = 1 2 X n=1;2 nc 0 [ i (i;j;n); s (i;j;n); s (i;j;n)]: (4.31) We note that the eective surface scattering area A(i;j) is available in the variable eff scatter of the CYGNSS L1 product. 4.2.6 Change of Polarimetric Basis from Linear to Circular for NBRCS The NBRCS nc 0 in (4.30) and (4.31) depends on the polarization of the trans- mitter and receiver. For CYGNSS, the transmitter is RHCP and the receiver is LHCP. We shall denote the corresponding NBRCS by 0 lr , where we drop the notation \nc" here and in the following. However, the SPM formulation (see Section 2.6 for references) is provided in linear polarization 0 pq for p;q2fh; vg, where h denotes horizontal and v denotes vertical. In the following, we derive an approximation for 0 lr in terms of 0 pq for p;q2fh; vg. From [13, eq. (5.143)], the NBRCS is given by 0 rt = 4 A D ^ p r S ^ p t 2 E (4.32) where r is receive polarization, t is transmit polarization, A is the area of the scatterer described by the 22 scattering matrixS, and ^ p is the two-dimensional antenna polarization unit vector. The angle brackets denote ensemble averaging associated with the random nature of the rough surface. With the convention of ^ k = ^ v ^ h, where ^ k is the unit vector in the direction of wave propagation, ^ v is the vertical unit vector, and ^ h is the horizontal unit vector, we have for left circular polarization ^ p l = 1 p 2 (^ p v i ^ p h ) (4.33) and for right circular polarization ^ p r = 1 p 2 (^ p v + i ^ p) h : (4.34) Thus, ^ p l S ^ p r = 1 2 [(S vv +S hh ) + i(S vh S hv )] (4.35) Hence, ^ p l S ^ p r 2 = 1 4 jS vv +S hh j 2 + 2< [i(S vh S hv )(S vv +S hh ) ] +jS vh S hv j 2 = 1 4 jS vv +S hh j 2 + 2= [(S vh S hv ) (S vv +S hh )] +jS vh S hv j 2 : (4.36) 31 If we were to continue to multiply out the above, we would obtain an expression in terms of the elements of the Mueller matrix or, equivalently, of the 4 4 polarimetric covariance matrix. However, for a re ectometer such as CYGNSS where the azimuth scattering angle is small, we expect the cross-polarized elements of the scattering matrix to be much smaller than the co-polarized elements. Thus, ^ p l S ^ p r 2 1 4 jS vv +S hh j 2 = 1 4 jS vv j 2 +jS hh j 2 + 2<(S hh S vv ) : (4.37) Therefore, 0 lr 4 A 1 4 D jS vv j 2 E + D jS hh j 2 E + 2< D S vv S hh E : (4.38) Dening = hS vv S hh i r D jS vv j 2 ED jS hh j 2 E (4.39) we nd 0 lr 1 4 0 vv + 2<() q 0 vv 0 hh + 0 hh : (4.40) Since1<() 1, we have the approximate bounds p 0 vv p 0 hh 2 ! 2 / 0 lr / p 0 vv + p 0 hh 2 ! 2 : (4.41) 4.2.7 Fresnel Footprint For a given DDM, we would like to be able to visualize the region of active coherent scattering. The question of the location and size of the active scattering region is discussed in Section 2.2 of [70], which concludes that the most important region is generally the rst Fresnel zone. The following two subsections describe how to calculate the rst Fresnel zone and its motion over the noncoherent integration period of a DDM. 4.2.7.1 Elliptical Approximation of First Fresnel Zone In the following, we approximate the surface within the rst Fresnel zone by a plane that is tangent to the specular point. We dene a local coordinate system consistent with [75] and centered at the specular point as follows. The z-axis is normal to the tangent plane. The y-axis is the intersection of the incidence plane and the tangent plane, with +y directed toward the transmitter. The x-axis completes the right-handed system. 32 The total path length is given by f(x;y) =R rs (x;y) +R st (x;y) (4.42) where range notation is the same as previously dened. Our goal is to nd an elliptical approximation to f(x;y) = constant. By the law of the cosines R rs (x;y) = (R sp rs ) 2 + 2R sp rs y sin i +x 2 +y 2 1 =2 (4.43) R st (x;y) = (R sp st ) 2 2R sp st y sin i +x 2 +y 2 1 =2 (4.44) where i is incidence angle. Thus, we nd the derivatives f x (x;y) =x 1 R rs + 1 R st (4.45) f y (x;y) = R sp rs R rs R sp st R st sin i +y 1 R rs + 1 R st (4.46) f xx (x;y) = 1 R rs + 1 R st +x ::: (4.47) f yy (x;y) = 1 R rs + 1 R st (R sp rs ) 2 (R rs ) 3 + (R sp st ) 2 (R st ) 3 sin 2 i +y ::: (4.48) f xy (x;y) =f yx (x;y) =x ::: : (4.49) Hence, f(0; 0) =R sp rs +R sp st (4.50) f x (0; 0) = 0 (4.51) f y (0; 0) = 0 (4.52) f xx (0; 0) = 1 R sp rs + 1 R sp st (4.53) f yy (0; 0) = 1 R sp rs + 1 R sp st cos 2 i (4.54) f xy (0; 0) =f yx (x;y) = 0 (4.55) which leads to the 2nd-order Taylor series expansion R rs +R st =R sp rs +R sp st + 1 2 1 R sp rs + 1 R sp st x 2 + 1 2 1 R sp rs + 1 R sp st y 2 cos 2 i : (4.56) Thus, we nd the elliptical approximation x b 2 + y a 2 = 1 (4.57) 33 Table 4.1: Fresnel zone statistics from day 001 of year 2018. Parameter Value Total number of DDMs 2 758 796 Maximum value of i 72:6° Maximum value of a 1162 m Minimum value of b 218:5 m where b = 2 R sp rs R sp st R sp rs +R sp st R 1 =2 (4.58) a =b sec i (4.59) R =R rs +R st R sp rs R sp st : (4.60) The rst Fresnel zone R=4 is then given by the elliptical region x b 2 + y a 2 1 (4.61) where b = R sp rs R sp st R sp rs +R sp st 2 1 =2 (4.62) a =b sec i : (4.63) We note that equivalent formulas are given in [115] without derivation. Parameter statistics from a single day of CYGNSS data are shown in Table 4.1. 4.2.7.2 Fresnel Footprint of Noncoherent Integration The ellipse of the rst Fresnel zone moves over the surface with specular point velocity ~ V sp for the noncoherent integration period T nc = 1 s, resulting in an elongated footprint. (In mid-2019, a change to the ight software was imple- mented to set T nc = 0:5 s. However, the results below are from 2018 only.) The maximum specular point velocity for all DDMs from day 100 of year 2018 is 6231 m s 1 . Thus, the specular point can move up to T nc 2 max (V sp ) = 3116 m (4.64) away from the center of the noncoherent integration period, which is where the L1 geometry data are sampled. The total Fresnel footprint can be represented as a long, narrow parallelogram with a half-ellipse appended to each end. To nd the corners of the parallelogram, we rst compute the two points where ~ V sp is tangent to the Fresnel ellipse. We 34 parameterize the ellipse in local coordinates by ~ r f () = 0 @ b cos a sin 0 1 A (4.65) where parameterizes the curve. Then, the tangent vector is given by the derivative d~ r f d () = 0 @ b sin a cos 0 1 A : (4.66) Thus, the tangent points~ r f (s tan ) must satisfy d~ r f d ( tan ) ~ V sp = ~ 0: (4.67) Writing the specular point velocity in local coordinates ~ V sp = 0 @ V sp x V sp y 0 1 A (4.68) it follows that bV sp y sin tan +aV sp x cos tan = 0: (4.69) Hence, cos tan = bV sp y r tan (4.70) sin tan = aV sp x r tan (4.71) where we have dened r tan = (aV sp x ) 2 + (bV sp y ) 2 1 =2 : (4.72) Thus, the tangent points are given by ~ r tan f = 0 @ b 2 V sp y =r tan a 2 V sp x =r tan 0 1 A : (4.73) The corners of the parallelogram are then given by ~ r corner = T nc 2 ~ V sp +~ r tan f : (4.74) The ellipses on each end of the parallelogram are given by xV sp x T nc =2 b 2 + yV sp y T nc =2 a 2 1: (4.75) The parallelogram corners ~ r corner and the end ellipses above are sucient to determine if the Fresnel zone passes over a given point during the noncoherent integration period of a given DDM. 35 4.3 Results To validate the rough surface scattering models of the previous section, we begin with the dry lake surface of Kati Thanda-Lake Eyre, Australia. This location has various desirable characteristics but lacks in-situ sensors. We then investigate sites with in-situ sensors and more complicated topographies and land covers. 4.3.1 Kati Thanda-Lake Eyre 4.3.1.1 Dataset Selection For an initial validation, wish to identify sites with the following characteristics: 1. Near sea level. Since the CYGNSS system was designed to operate at sea level, behavior at higher elevations might be more dicult to interpret. 2. No standing water. Standing water produces strong specular scattering not associated with the soil parameters we wish to retrieve. 3. Large, at, level, and bare. The site should be suciently large to ac- commodate at least the rst Fresnel zone over the 1-second CYGNSS noncoherent integration period, suciently at to generate at least some specular re ection, suciently level that the re ection can be geolocated, and suciently bare that attenuation due to vegetation is negligible. 4. In-situ soil data. The site should have instrumentation to provide soil moisture and/or dielectric constant for input to the rough surface scattering models. We are not aware of any sites satisfying all four criteria. Nevertheless, we found news reports that Kati Thanda-Lake Eyre ooded in mid- to late-May of 2018 after a period of drought. Kati Thanda-Lake Eyre is a very large, at, level area with minimum elevation of15 m located in the deserts of central Australia. The oodwaters arrived without any local rain after traveling nearly 1000 km from the outback region of north Queensland where they had originated about two months earlier. Prior to arrival of the oodwaters, conditions in the Kati Thanda-Lake Eyre area were described as \parched" as no local rain had fallen there for eight months. It had been several years since Kati Thanda-Lake Eyre last had a ood of comparable magnitude [116]. We therefore choose to start looking at CYGNSS data over Kati Thanda-Lake Eyre on ordinal (aka Julian) date 100, which corresponds to April 10 and which is well before the arrival of the oodwaters, in the expectation of fullling the rst three criteria. Fig. 4.1 shows samples crossing the Kati Thanda-Lake Eyre area on day 100 of year 2018 from all seven spacecraft available on that day. By convention, samples are indexed from zero. Fig. 4.2 shows samples from spacecraft ight model (FM) 1 only, which tracked the re ection of GPS space vehicle number (SVN) 44 with pseudo-random noise (PRN) code 28 on DDM channel 1 across 36 Figure 4.1: Samples crossing Kati Thanda-Lake Eyre on day 100 of year 2018. the dry lake in a southeasterly direction. Since SVN 44 has block type of IIR, we do not expect it to have the power exing problems associated with block IIF as described in [117]. We select sample number 54 549, which is located in the middle of the east leg of the ground track, for further analysis. This sample has an incidence angle of 31:21°. The DDM SNR from L1 variable ddm snr is 19:15 dB. The only potentially unexpected quality control ag is bit 18, indicating detection of RFI. This ag is set when the kurtosis of the DDM noise oor diers from the ideal Gaussian value of 3:0 by more than 1:0 [117]. The L1 BRCS is shown in Fig. 4.3. 4.3.1.2 Coherent Modeling The coherent component of the BRCS DDM model given by (4.12) and corre- sponding to the CYGNSS data identied in the previous subsection is shown in Fig. 4.4. Here, we have assumed that the surface is smooth (s = 0 m) so as to zero 37 Figure 4.2: Samples crossing Kati Thanda-Lake Eyre on day 100 of year 2018 from FM 1. 38 0 1 2 3 4 5 6 7 8 9 10 Doppler Bin 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 Delay Bin Bistatic RCS from L1 Data 80 85 90 95 100 105 110 115 120 dBsm Figure 4.3: CYGNSS L1 BRCS data from year 2018, day 100, FM 1, channel 1, sample 54 549, over Kati Thanda-Lake Eyre during dry season. 39 0 1 2 3 4 5 6 7 8 9 10 Doppler Bin 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 Delay Bin Bistatic RCS from Forward Model 80 85 90 95 100 105 110 115 120 dBsm Figure 4.4: Forward model for coherent component of BRCS for CYGNSS DDM from year 2018, day 100, FM 1, channel 1, sample 54 549. Dielectric constant is assumed to be 6:27 + 0:627i, and surface roughness assumed to be zero. out the noncoherent component as we expect coherent scattering to dominate for this location. Since no in-situ sensor data are available for computing the soil dielectric constant, we have adjusted the dielectric constant value so as to match the total BRCS of the nine middle bins of the model with the corresponding data. Thus, we obtain a dielectric constant value of r = 6:27 + 0:627i. Comparing the model in Fig. 4.4 and the corresponding L1 BRCS data in Fig. 4.3, we nd good agreement between the model and the data down to at least 20 dB below the peak. One small dierence is that Doppler bins2 have more power in the measurement than in the coherent model. 4.3.1.3 Noncoherent Modeling To assess the validity of the SPM for nc 0 in (4.30) and (4.31), we use the Kati Thanda-Lake Eyre dataset with the dielectric constant value of the previous section to synthesize the noncoherent component of the BRCS DDM. Our objective is not to match the CYGNSS data but rather only to check the 40 First-Order SPM NBRCS Figure 4.5: Delay-Doppler map of nc 0 synthesized using the rst-order SPM and based on the CYGNSS geometry from year 2018, day 100, FM 1, channel 1, sample 54 549 over Kati Thanda-Lake Eyre. reasonableness of results. To this end, we use the rst-order SPM formulation of [79, eqs. (21-67) and (21-68)]. We assume the rough scattering surface is isotropic and exponentially correlated with correlation length ten times greater than the height standard deviation s. To maximize response, we set the surface roughness parameter to ks = 0:3, which represents the upper limit of the roughness validity region for the SPM [13, Section 10-3] and which corresponds to a surface height standard deviation ofs = 9:1 mm at the CYGNSS wavelength. We use the upper bound from (4.41) to convert nc 0 to the circular polarization of the CYGNSS mission. The resulting map of nc 0 (i;j) from (4.31) is shown in Fig. 4.5. We see that the value of nc 0 is nearly constant across the DDM. The slightly nonsymmetric appearance of the DDM is due to the fact that the iso-Doppler lines are not aligned with the plane of incidence for this particular case. 41 Figure 4.6: Eective area A(i;j) from CYGNSS L1 variable eff scatter from year 2018, day 100, FM 1, channel 1, sample 54 549 over Kati Thanda-Lake Eyre. The corresponding map of eective area A(i;j) from the CYGNSS L1 variable eff scatter is shown in Fig. 4.6. We observe two unexpected features. First, eective area in Doppler bin 5 is identically zero for all delays. We instead expect the eective areas in Doppler bin 5 to be similar to those in Doppler bin5. Second, except for the anomalous Doppler bin 5, the eective area does not go to zero with increasingly negative delay. Rather, we see that the eective area goes nearly to zero in delay bin4 and then starts increasing again for more negative delays. Since negative delays represent scattering from a volume located in the atmosphere above the ground (i.e., from noise bins), we expect the corresponding eective areas to be zero within the 4-bin width of the WAF. Thus, the nonzero values from delay bin5 to delay bin8 are unexpected. Nevertheless, the values in the bins near the center of the DDM appear to be reasonable. The center bin has an eective area of 84:4 dBsm, which is equivalent to that of a circle with radius 9:4 km. 42 First-Order SPM BRCS Figure 4.7: Delay-Doppler map of BRCS synthesized using the rst-order SPM and based on the CYGNSS geometry and L1 eective areas from year 2018, day 100, FM 1, channel 1, sample 54 549 over Kati Thanda-Lake Eyre. The corresponding noncoherent BRCS DDM obtained by multiplying nc 0 (i;j) and A(i;j) according to (4.30) is shown in Fig. 4.7. As expected, the eective area artifacts described in a previous paragraph are apparent. Nevertheless, the artifacts are at least 20 dB below the peak noncoherent response. To compare with the coherent response, we rst calculate the coherent attenuation (3.5) due to the assumed surface roughness of ks = 0:3: A rough = exp[(2ks cos i ) 2 ] = 0:768 =1:1437 dB: (4.76) Applying the attenuation to the results of Fig. 4.4, we nd that the synthesized coherent DDM has a peak BRCS of about 118 dBsm. By comparison, the noncoherent peak from Fig. 4.7 is about 83 dBsm, which is a dierence of about 35 dB. Thus, even with surface roughness parameterks set to the greatest value allowed by the SPM, the noncoherent contribution is still negligible relative to the coherent contribution. 43 Figure 4.8: The author, appearing on the left of the photo, digs a hole for soil moisture sensor installation at the SoilSCAPE site of Tonzi Ranch in March 2019. Photo credit: Ruzbeh Akbar. 4.3.2 SoilSCAPE Sites 4.3.2.1 Background The Soil moisture Sensing Controller and oPtimal Estimator (SoilSCAPE) project, which is funded through the NASA ESTO Advanced Information Sys- tems Technology (AIST) program, aims to provide local scale surface-to-depth es- timates of soil moisture with optimal sampling [118,119]. The best-characterized SoilSCAPE site is Tonzi Ranch, where the University of Southern California (USC) Microwave Systems, Sensors, and Imaging Laboratory (MiXIL) research group has already performed detailed soil and vegetation surveys in support of calibration and validation (cal/val) activities for the AirMOSS and SMAP missions. Located east of Sacramento, CA, USA, in the foothills of the Sierra Nevada, Tonzi Ranch is classied as oak savanna woodland. Photographs of eldwork performed at Tonzi Ranch in March 2019, which the author supported, are provided in Figs. 4.8{4.10. In the following subsections, we investigate the application of the scattering models for planar rough surfaces developed earlier in this chapter to the SoilSCAPE site at Tonzi Ranch. 44 Figure 4.9: Installing a cattle fence around a new node at Tonzi Ranch in March 2019. Photo credit: Negar Golestani. Figure 4.10: Sunset on Tonzi Ranch in March 2019. Photo credit: Negar Golestani. 45 Figure 4.11: Fresnel footprint for sample 50420 from channel 3 of FM 1 on day 66 of year 2018 based on location of specular point reported in L1 data. The Fresnel parallelogram is shown in magenta and the end ellipses are plotted in red. Placemarks labeled with sample number designate L1 specular point locations. The cluster of unlabeled placemarks indicates in-situ sensor locations at SoilSCAPE site of Tonzi Ranch. 4.3.2.2 Fresnel Footprint An initial survey of the rst 100 days of year 2018 identied the L1 specular point reported for sample 50420 from channel 3 of FM 1 on day 66 to be located closest to the soil moisture sensors at Tonzi Ranch. The Fresnel footprint of this DDM as calculated by the method of Section 4.2.7 is shown in Fig. 4.11. The parallelogram of the Fresnel footprint is drawn in magenta and the end ellipses in red. Placemarks labeled with sample number designate L1 specular point locations. Unlabeled placemarks indicate the in-situ soil moisture sensor locations. The nearest sensor is located about 1:7 km from the reported specular point. As seen in the gure, the computed Fresnel footprint includes some of the in-situ sensors. An examination of the elevation data from Google Earth reveals that elevation varies more than 10 m across the semimajor axis any given Fresnel ellipse within the Fresnel footprint. By comparison, a surface height standard deviation of only s = 20 cm with i = 60° gives a coherent attenuation factor (3.5) of A rough = exp[(2ks cos i ) 2 ] = 1:22 10 19 =189 dB: (4.77) 46 Thus, the roughness of the land surface at Tonzi Ranch is well outside the regime of measurable coherent scattering. Furthermore, from west to east across the entire Fresnel footprint, we nd that elevation increases from 63 to 193 m. Large-scale undulations such as this can give rise to multiple specular points, known as points of stationary phase, which add in power rather than in amplitude and which can be distributed over a region much larger than the rst Fresnel zone. A model that accommodates large-scale undulations will be investigated in the next chapter. In any case, the topography of the Tonzi Ranch site does not seem to satisfy the assumption of a planar rough surface that is required by the coherent model (4.12) with re ectivity (3.3) and roughness attenuation factor (3.5). We also observe various small inland water bodies, both inside the Fresnel footprint and in the general vicinity. One of these water bodies can be seen in the cell-phone photo of Fig. 4.10. The modeling of inland water bodies in GNSS-R data is an area of ongoing research. As an aside, we also examined L1 specular point locations near the SoilSCAPE site at Lucky Hills, AZ, USA, and we found the topography there to be even more hilly than at Tonzi Ranch on scales of hundreds of meters, although without as much gradient on scales of kilometers. Thus, the land surface roughness at Lucky Hills is also well outside the regime of measurable coherent scattering. The raw counts from the DDM over Tonzi Ranch identied above are shown in Fig. 4.12. No clear signal, such as the WAF point response function seen over Lake Eyre or the typical horseshoe shape associated with diuse ocean scattering, can be seen in this DDM. Indeed, the positive delay bins corresponding to land surface scattering have no more raw counts on average than the thermal noise in the negative delay bins corresponding to points above the Earth's surface. Indeed, we have found that many DDMs over Tonzi Ranch have low SNR. A more comprehensive search to identify useable DDMs over Tonzi Ranch is described next. 4.3.2.3 Tonzi Ranch Dataset Identication To search for useable datasets over Tonzi Ranch, we downloaded version 2.1 of the CYGNSS L1 data for the entire year of 2018. We estimate the total number of DDMs in 2018 to be nearly one billion since, if the system were operating at full capacity, we would have (1 year) (356 day year 1 ) (86 400 s day 1 ) (8 spacecraft s 1 ) (4 channel spacecraft 1 ) (1 DDM channel 1 ) = 1 009 152 000 DDM: (4.78) As shown in Fig. 4.13, we nd a total of 137 DDMs in 2018 whose L1 specular points are located within a radius of 10 km of the SoilSCAPE sensors at Tonzi Ranch. Of these, a total of 6 DDMs have an L1 SNR value above 1:5 dB (shown as green circles) while the remaining 131 DDMs have a value below 1:5 dB (shown as red circles). None of the DDMs have an L1 SNR value above 2:0 dB. For 47 Figure 4.12: Delay-Doppler map of raw counts from sample 50420 from channel 3 of FM 1 on day 66 of year 2018 over SoilSCAPE site of Tonzi Ranch. 48 Figure 4.13: L1s specular point locations of all version 2.1 CYGNSS DDMs from year 2018 located within a radius of 10 km of the SoilSCAPE sensors at Tonzi Ranch. Green circles indicate L1 SNR above 1:5 dB. Red circles indicate L1 SNR below 1:5 dB. The labels on the circles indicate the L1 sample number within the respective day. (Clicking on a circle in Google Earth will highlight the folder with the year, day, FM number, and channel number of the sample, thereby uniquely identifying the corresponding DDM.) The ag indicates the location of the SoilSCAPE sensors. comparison, the Kati Thanda-Lake Eyre DDM has an SNR of 19:15 dB. Also, a recent ESA-funded study used an SNR value of 3 dB as a threshold for data quality control [101]. In spite of the low values of SNR, typical structures can be seen in the six green DDMs as shown in Figs. 4.14 to 4.19. In particular, a specular-like point is evident in Figs. 4.14 to 4.16. Additionally, a horseshoe-like pattern, suggestive of diuse scattering, is apparent in Figs. 4.17 to 4.19. The six DDMs are distributed among three tracks. The days corresponding to the three tracks are 8, 24, and 140. Thus, no two tracks occur within a week of each other. We conclude that the six DDMs could potentially be used for forward model validation. However, they could not realistically support 49 Figure 4.14: BRCS DDM: year 2018, day 008, FM 3, channel 3, sample 77504. validation of a soil moisture retrieval algorithm based on a time series method. We use the track containing the rst DDM as a validation case study for the topography-based scattering model developed in the next chapter. 50 Figure 4.15: BRCS DDM: year 2018, day 008, FM 3, channel 3, sample 77505. Figure 4.16: BRCS DDM: year 2018, day 024, FM 3, channel 2, sample 50372. 51 Figure 4.17: BRCS DDM: year 2018, day 140, FM 4, channel 3, sample 11048. Figure 4.18: BRCS DDM: year 2018, day 140, FM 4, channel 3, sample 11049. 52 Figure 4.19: BRCS DDM: year 2018, day 140, FM 4, channel 3, sample 11050. 4.4 Discussion We have proposed a model for the full DDM of scattering from a bare planar rough surface, where the coherent component (4.12) is based on the KA and the noncoherent component (4.22) is based on the SPM. To evaluate the noncoherent surface integral, we have obtained the approximations (4.30) and (4.31) and the circular polarization bounds (4.41). Finally, we have derived the gure of the rst Fresnel zone (4.61), (4.62), and (4.63) and the Fresnel footprint of noncoherent integration (4.74) and (4.75). As an initial validation of the model, we identied a DDM of a large, at, bare region of Kati Thanda-Lake Eyre during the dry season and found a plausible value of dielectric constant that provided a good match between the model and the data under the assumption of a smooth surface. Without in-situ soil measurements and high-resolution surface roughness data, however, this validation is of limited utility. Using the Kati Thanda-Lake Eyre validation dataset as a basis for simulation, we increased the surface roughness parameter up to the highest value allowed by the SPM and we found that the noncoherent component remained at least 35 dB below the coherent component. This suggests that the proposed model is incapable of representing a mix of coherent and noncoherent scattering. That 53 is, the noncoherent component is always negligible relative to the coherent component. This result casts doubt on the utility of the SPM for modeling noncoherent scattering in re ectometry applications. Moving our attention to validation sites having in-situ soil moisture sensors, analysis in Google Earth of the surface elevation within the rst Fresnel zone and within the Fresnel footprint at the SoilSCAPE sites of Tonzi Ranch and Lucky Hills indicates that surface roughness at these sites is well outside the regime of measurable coherent scattering. That is, the model predicts that the coherent component of scattering is negligible. This result underscores the need for a noncoherent model that can accommodate large-scale undulations and topography. Such a model is investigated in the next chapter. We found that DDMs over the preferred SoilSCAPE validation site of Tonzi Ranch have low SNR. Nevertheless, we identied a total of six Tonzi DDMs from the year 2018 that might have enough signal to be of interest. We also noted the presence of various small inland water bodies, which could potentially complicate the observation of soil moisture. 54 Chapter 5 Models for Rough Surfaces with Topography 5.1 Introduction In this chapter, a model is developed for the DDM of a bare surface for which a DEM is available. The model is derived from rst principles in the GO limit of the KA. This work grew out of eorts described in the previous chapter to dene a model under the assumption that the rough scattering surface is approximately planar with spatially-invariant statistics. We found that the terrain around the SoilSCAPE validation sites of Tonzi Ranch and Walnut Gulch fails to satisfy this assumption. Instead, we found the terrain to be more similar to a typical rough sea surface having large-scale undulations and scattering in the strong diuse regime. Unlike a sea surface, however, the land surface undulations are not entirely random. Rather, they are described at least to some degree by a DEM. As already noted in Section 2.6, DEM-based matchups have already been reported [89,90], although the topography-specic details of the model are still forthcoming. Additionally, a DEM-based method using level rough patches has been pro- posed to account for power ratios that are less than predicted from a strictly planar rough surface but greater than predicted from purely noncoherent scat- tering [120,121]. The method has not yet been validated with GNSS-R data. The derivation in the present chapter closely follows ZV2000 [72] and Rytov, Kravtsov, and Tatarskii [68, Section 5.2], which were introduced in Section 2.6. However, instead of considering the surface height to be a purely random eld as in ZV2000, we decompose the surface height into a deterministic part obtained from a DEM and a random part representing the residual height between the DEM and the rough land surface. This decomposition is similar to that of Bass and Fuks [69, Section 30] for scattering from a body of nite dimension. 55 As noted in Section 2.6, ZV2000 has been generalized in VZ2018 [75] to remove the KA under the assumption of spatially homogeneous roughness statistics. The extension of VZ2018 to land applications with nonhomogeneous statistics is beyond the scope of the present work. The material in this chapter is based on a paper that the author has submitted to the special issue on \Recent Advances in Re ectometry using GNSS and other Signals of Opportunity" of the Journal of Selected Topics in Applied Earth Observations and Remote Sensing (JSTARS) of the Geoscience and Remote Sensing Society (GRSS) of the Institute of Electrical and Electronics Engineers (IEEE). The author presented the results in this chapter in a poster session of the IEEE GRSS Specialist Meeting on Re ectometry using GNSS and Other Signals of Opportunity held in Benevento, Italy, 20-22 May, 2019. This meeting is associated with the JSTARS special issue. The author also gave these results in an oral presentation at the meeting of the CYGNSS science team meeting held in Ann Arbor, MI, USA, 5-7 July, 2019. Section 5.2 derives the general DEM-based model and also provides some CYGNSS-specic details. Section 5.3 presents results of a validation case study with CYGNSS data. Section 5.4 discusses the results, and Section 5.5 identies directions for future work. 5.2 Theoretical Model Subsection 5.2.1 develops the model under the assumption that the surface height is deterministic and known exactly. Subsection 5.2.2 decomposes the height into a deterministic part from the DEM and a random residual height. A DEM-based bistatic radar equation is derived with a DEM-based expression for the normalized BRCS. Subsection 5.2.3 describes CYGNSS-specic modeling details of the BRCS DDM and the GPS C/A-code WAF. 5.2.1 Deterministic Modeling For an electric eld satisfying the scalar Helmholtz equation, the scattered eld at receiver position ~ R r is given by E sca ( ~ R r ) = 1 4 Z E tot ( ~ R s ) @ @n g 0 ( ~ R r ; ~ R s )g 0 ( ~ R r ; ~ R s ) @ @n E tot ( ~ R s ) dS (5.1) where is the scattering surface, E tot is the total electric eld amplitude, ~ R s is the position on of the surface integration variable, @=@n is the surface normal derivative, and g 0 is the free-space Green's function given by g 0 ( ~ R r ; ~ R s ) = exp(ikR rs ) R rs (5.2) 56 1 ⇣ = ⇣ dem + ⇣ res ⇣ dem ⇣ res Receiver at ¯ R Rx Transmitter at ¯ R Tx Rough surface 1 ⇣ = ⇣ dem + ⇣ res ⇣ dem ⇣ res Receiver at ¯ R Rx Transmitter at ¯ R Tx Rough surface 1 ⇣ = ⇣ dem + ⇣ res ⇣ dem ⇣ res Receiver at ¯ R Rx Transmitter at ¯ R Tx Rough surface Earth ellipsoid DEM Rough Surface 1 ⇣ = ⇣ dem + ⇣ res ⇣ dem ⇣ res ¯ V Tx ¯ V Rx ¯ R st = ¯ R s ¯ R t ¯ R rs = ¯ R r ¯ R s ˆ m ˆ n Receiver at ¯ R Rx Transmitter at ¯ R Tx Earth ellipsoid DEM Rough surface 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ¯ V Tx ¯ V Rx ¯ R st = ¯ R s ¯ R t ¯ R rs = ¯ R r ¯ R s ˆ m ˆ n ˆ z R s =¯ ⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ¯ R Rx Transmitter at ¯ R Tx Earth ellipsoid DEM Rough surface 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ¯ V Tx ¯ V Rx ¯ R st = ¯ R s ¯ R t ¯ R rs = ¯ R r ¯ R s ˆ m ˆ n ˆ z ¯ ⇢ ¯ R s =¯ ⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ¯ R Rx Transmitter at ¯ R Tx Earth ellipsoid DEM Rough surface Surface point at 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ¯ V t ¯ V r ¯ R st = ¯ R s ¯ R t ¯ R rs = ¯ R r ¯ R s ˆ m ˆ n ˆ z ¯ ⇢ 2 ⌃ 0 ¯ R s =¯ ⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ¯ R r Transmitter at ¯ R t Integration surface⌃ 0 DEM Rough surface Surface point at 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ¯ V t ¯ V r ¯ R st = ¯ R s ¯ R t ¯ R rs = ¯ R r ¯ R s ˆ m ˆ n ˆ z ¯ ⇢ 2 ⌃ 0 ¯ R s =¯ ⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ¯ R r Transmitter at ¯ R t Integration surface⌃ 0 DEM Rough surface⌃ Surface point at 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ¯ V t ¯ V r ¯ R st = ¯ R s ¯ R t ¯ R rs = ¯ R r ¯ R s ˆ R st ˆ R rs ˆ m ˆ n ˆ z ¯ ⇢ 2 ⌃ 0 ¯ R s =¯ ⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ¯ R r Transmitter at ¯ R t Integration surface⌃ 0 DEM Rough surface⌃ Surface point at 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ¯ V t ¯ V r ¯ R st = ¯ R s ¯ R t ¯ R rs = ¯ R r ¯ R s ˆ R st ˆ R rs ˆ m ˆ n ˆ z ¯ ⇢ 2 ⌃ 0 ¯ R s =¯ ⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ¯ R r Transmitter at ¯ R t Integration surface⌃ 0 DEM Rough surface⌃ Surface point at 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ~ V t ~ V r ~ R st = ~ R s ~ R t ~ R rs = ~ R r ~ R s ˆ R st ˆ R rs ˆ m ˆ n ˆ z ~⇢ 2 ⌃ 0 ~ R s = ~⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ~ R r Transmitter at ~ R t Integration surface⌃ 0 DEM Rough surface⌃ Surface point at 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ~ V t ~ V r ~ R st = ~ R s ~ R t ~ R rs = ~ R r ~ R s ˆ R st ˆ R rs ˆ m ˆ n ˆ z ~⇢ 2 ⌃ 0 ~ R s = ~⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ~ R r Transmitter at ~ R t Integration surface⌃ 0 DEM Rough surface⌃ Surface point at 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ~ V t ~ V r ~ R st = ~ R s ~ R t ~ R rs = ~ R r ~ R s ˆ R st ˆ R rs ˆ m ˆ n ˆ z ~⇢ 2 ⌃ 0 ~ R s = ~⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ~ R r Transmitter at ~ R t Integration surface⌃ 0 DEM Rough surface⌃ Surface point at 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ~ V t ~ V r ~ R st = ~ R s ~ R t ~ R rs = ~ R r ~ R s ˆ R st ˆ R rs ˆ m ˆ n ˆ z ~⇢ 2 ⌃ 0 ~ R s = ~⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ~ R r Transmitter at ~ R t Integration surface⌃ 0 DEM Rough surface⌃ Surface point at 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ~ V t ~ V r ~ R st = ~ R s ~ R t ~ R rs = ~ R r ~ R s ˆ R st ˆ R rs ˆ m ˆ n ˆ z ~⇢ 2 ⌃ 0 ~ R s = ~⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ~ R r Transmitter at ~ R t Integration surface⌃ 0 DEM Rough surface⌃ Surface point at 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ~ V t ~ V r ~ R st = ~ R s ~ R t ~ R rs = ~ R r ~ R s ˆ R st ˆ R rs ˆ m ˆ n ˆ z ~⇢ 2 ⌃ 0 ~ R s = ~⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ~ R r Transmitter at ~ R t Integration surface⌃ 0 DEM Rough surface⌃ Surface point at 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ~ V t ~ V r ~ R st = ~ R s ~ R t ~ R rs = ~ R r ~ R s ˆ R st ˆ R rs ˆ m ˆ n ˆ z ~ r2 ⌃ 0 ~ R s = ~⇢ +(⇣ dem +⇣ res )ˆ z Receiver at ~ R r Transmitter at ~ R t Integration surface⌃ 0 DEM Rough surface⌃ Surface point at 1 ⇣ = ⇣ dem +⇣ res ⇣ dem ⇣ res ~ V t ~ V r ~ R st = ~ R s ~ R t ~ R rs = ~ R r ~ R s ˆ R st ˆ R rs ˆ m ˆ n ˆ z ~ r2 ⌃ 0 ~ R s = ~ r+(⇣ dem +⇣ res )ˆ z Receiver at ~ R r Transmitter at ~ R t Integration surface⌃ 0 DEM Rough surface⌃ Surface point at Figure 5.1: Scattering geometry. where R rs =j ~ R r ~ R s j and k is wavenumber. A diagram of the scattering geometry is provided in Fig. 5.1 for reference throughout this section. Equation (5.1) is known variously as the Helmholtz integral [70, eq. (8) of Ch. 3], the Helmholtz-Kirchho formula [122, eq. (6-17)], a Green's formula [69, eq. (1.4)], and a form of Huygens' principle [123, eq. (8.1.10)]. From the KA, we nd E tot ( ~ R s ) = (1 +R)E inc ( ~ R s ) and (5.3) @E tot ( ~ R s ) @n = (1R) @E inc ( ~ R s ) @n (5.4) where E inc is the incident electric eld amplitude andR is the polarization- sensitive Fresnel re ection coecient. Thus, the scattered electric eld amplitude in the KA is given by [68, eq. (5.34); 69, eq. (19.8); 72, eq. (4)] E sca ( ~ R r ) = 1 4 Z R( ~ R s ) @ @n E inc ( ~ R s ) exp(ikR rs ) R rs dS (5.5) The signal provided by the receive antenna is given by s ant ( ~ R r ) = 1 4 Z A r ( ~ R s )R( ~ R s ) @ @n E inc ( ~ R s ) exp(ikR rs ) R rs dS (5.6) where A r is proportional to the receive antenna pattern expressed as a complex- valued amplitude. In the case that A r is unity, (5.6) reduces to (5.5). We choose to dene A r such thatjs ant j 2 has units of received power in watts. 57 In the special case that the integrand is nonzero over only a small part of the surface so that the angle of arrival of the scattered eld at the antenna is nearly constant, we nd that A r can be moved outside the surface integral. It then follows from (5.5) and (5.6) that received power is given by P r =js ant j 2 =jA r j 2 jE sca j 2 : (5.7) But from antenna theory [13, Sections 2-5 and 5-4], received power is also given by P r =S r A e (5.8) where the average power density at the receive antenna is given by S r = jE sca j 2 =(2 0 ) and antenna eective area is given by A e = 2 4 G r . Here, 0 is the impedance of free space, is wavelength, and G r is receive antenna gain. Equating (5.7) and (5.8), we nd jA r ( ~ R s )j 2 = 1 2 0 2 4 G r ( ~ R s ): (5.9) Here, the dependence of A r and G r on scatterer position ~ R s is shown explicitly. Returning to (5.6) and following ZV2000 [72, eqs. (1)-(3)] for the emitted signal, a model for the amplitude of the electric eld incident on the surface is given by E inc ( ~ R s ;t) = A t ( ~ R s ) R st a(tR st =c) exp(ikR st i2f t D t) (5.10) where t is time relative to the reference time of the coherent integration period T i to be discussed below, A t is proportional to the transmit antenna pattern as a complex amplitude, R st =j ~ R s ~ R t j is the distance from the transmitter position ~ R t to the surface point ~ R s and is assumed to remain xed over T i (i.e., stop-and-hop), a is the modulation function of the transmitted waveform withjaj 2 = 1 in a time-averaged sense, and f t D is the Doppler frequency due to transmitter motion. The time-dependent terms in (5.10) are assumed to change slowly relative to the carrier e i2fct , which is an implicit factor and where f c is the carrier frequency. The average power density incident on the surface is then given by S inc ( ~ R s ) = jE inc ( ~ R s )j 2 2 0 = jA t ( ~ R s )j 2 2 0 R 2 st : (5.11) But from antenna theory we also have S inc ( ~ R s ) = P t G t ( ~ R s ) 4R 2 st (5.12) whereP t is transmitted power and G t is transmit antenna gain. Equating (5.11) and (5.12), we nd jA t ( ~ R s )j 2 = 0 G t ( ~ R s )P t 2 : (5.13) 58 We now substitute (5.10) into (5.6). Since (5.6) was derived for a time- harmonic eld without the modulation factor a and without receiver motion, we must account for the additional propagation delay R rs =c from the surface to the receiver ina and for the additional Doppler due to receiver motion when making this substitution. Thus, we nd s ant ( ~ R r ;t) = 1 4 Z A t ( ~ R s )A r ( ~ R s )R( ~ R s ) (5.14) @ @n a[t (R rs +R st )=c] exp[ik(R rs +R st )] R rs R st e i2fDt dS where the Doppler frequency is given by f D = ( ~ V t ^ R st ~ V r ^ R rs )= (5.15) and ~ V t is transmitter velocity, ~ V r is receiver velocity, and the unit vectors are given by ^ R rs = ~ R r ~ R s j ~ R r ~ R s j and ^ R st = ~ R s ~ R t j ~ R s ~ R t j : (5.16) We note that (5.15) and (5.16) are consistent with (4.16) and (4.17), respectively. Additionally, from (5.9) and (5.13) we nd jA t ( ~ R s )j 2 jA r ( ~ R s )j 2 = 2 (4) 2 P t G t ( ~ R s )G r ( ~ R s ) (5.17) which we shall make use of later. Since exp[ik(R rs +R st )] is the only factor inside the braces of (5.14) that does not vary slowly with ^ R s , the rest of the factors can be moved outside the normal derivative. We then compute @ @n exp[ik(R rs +R st )] = [^ nr s ik(R rs +R st )] exp[ik(R rs +R st )] =i(^ n~ q ) exp[ik(R rs +R st )] (5.18) where ^ n is the surface normal, the gradient is taken with respect to ~ R s , and the scattering vector is dened by ~ q =r s k(R rs +R st ) =k( ^ R rs ^ R st ): (5.19) We now change to integration over a planar or locally planar surface 0 . For example, 0 could represent an Earth ellipsoid or geoid model. Let coordinates x and y lie in 0 and coordinate z be the signed distance from 0 with positive up. Dening d~ r = dx dy and noting that dS = dx dy ^ n ^ z = d~ r ^ n ^ z (5.20) 59 we nd from (5.14), (5.18), and (5.20) that s ant ( ~ R r ;t) = 1 4i Z 0 A t ( ~ R s )A r ( ~ R s )R( ~ R s )a[t (R rs +R st )=c] exp[ik(R rs +R st )] R rs R st e i2fDt ^ n~ q ^ n ^ z d~ r: (5.21) In the KA, the only parts of that scatter are those whose scattering vector ~ q is aligned with the surface normal ^ n. These are the points of stationary phase. Thus, we can replace ^ n by ~ q under the surface integral to obtain s ant ( ~ R r ;t) = Z 0 a[t (R rs +R st )=c]e i2fDt g(~ r ) d~ r (5.22) where g(~ r ) = 1 4i A t ( ~ R s )A r ( ~ R s )R( ~ R s ) exp[ik(R rs +R st )] R rs R st q 2 q z (5.23) and q =j~ qj and q z =~ q ^ z. The DDM instrument cross-correlates the signal s ant received from the antenna with a replica of the code a as follows Y (i;j) = 1 T i Z Ti 0 s ant (t)a (t i )e i2fjt dt (5.24) for each bin (i;j) of the DDM where i is the time delay of delay bin i and f j is the Doppler frequency of Doppler bin j. Substituting (5.22) into (5.24), exchanging the order of integration, and dropping phase factors that are constant in time, we nd Y (i;j) = Z 0 g(~ r )( i ;f j ;t 0 ) d~ r (5.25) where (;f;t 0 ) = 1 T i Z t0+Ti t0 a(t)a (t)e i2ft dt (5.26) is the WAF and where t 0 = (R rs +R st )=c (5.27) and i = i (R rs +R st )=c (5.28) f j =f j f D : (5.29) We note that (5.26), (5.28), and (5.29) are consistent with (4.6), (4.14), and (4.15), respectively. 60 The instrument performs square law detection and noncoherent time averaging of the correlator outputY . Letting the intervals of coherent correlation be indexed by n, we nd the DDM of noncoherently averaged received power is given by P g (i;j) = 1 N N X n=1 P g (i;j;n) (5.30) where P g (i;j;n) =jY (i;j;n)j 2 (5.31) andN is the number of coherent integration intervals in the noncoherent averag- ing. The subscript g, included for consistency with [73], indicates that the power is from the GNSS source and rather than from thermal noise in the receiver and antenna. 5.2.2 Stochastic Modeling We assume that, over the region of support of the integral, the rough surface can be represented as a single-valued function(~ r ) of the vector~ r = (x;y; 0)2 0 with = dem + res (5.32) where dem is a deterministic height representing topography, such as from a DEM, and res is a random residual height representing surface roughness relative to the DEM. From (5.30) and (5.31), the mean power is given by hP g (i;j)i = 1 N N X n=1 hP g (i;j;n)i (5.33) where hP g (i;j;n)i = D jY (i;j;n)j 2 E (5.34) and angle brackets denote ensemble averaging with respect to residual roughness res and, if applicable, the function . Due to topography, the statistics of the surface heights are not necessarily translationally invariant. Thus, for land applications, unlike for sea scattering, the noncoherent averaging operation cannot necessarily be dropped by invoking ergodicity. That is, depending on topography, we may not be able to equate hP g (i;j)i andhP g (i;j;n)i. We now assume that res is suciently small that it can be neglected in all factors of (5.25) except for in the complex sinusoid in the factor g dened by (5.23). By (5.19), a rst-order Taylor series approximation of the phase of the complex sinusoid about = 0 yields k(R rs +R st ) =k(R rs +R st )j =0 q z : (5.35) 61 Thus, from (5.25) and (5.34) we nd that hP g (i;j;n)i = Z 0 hg 0 0 g 00 00 i d~ r 0 d~ r 00 = Z 0 (~ r 0 ;~ r 00 )(~ r 0 ;~ r 00 ) d~ r 0 d~ r 00 (5.36) where (~ r 0 ;~ r 00 ) = 1 16 2 A t A r Rq 2 R rs R st q z 0 A t A r Rq 2 R rs R st q z 00 h 0 00 i (5.37) exp[ik(R 0 rs +R 0 st R 00 rs R 00 st )j =0 iq 0 z 0 dem + iq 00 z 00 dem ] and (~ r 0 ;~ r 00 ) =hexp(iq 0 z 0 res + iq 00 z 00 res )i: (5.38) Here, the notation of prime and double prime indicates evaluation at ~ r 0 and~ r 00 , respectively. We observe that (~ r 0 ;~ r 00 ) has the form of a double point charac- teristic function. We now assume that the support of lies within the region j~ r 00 ~ r 0 j<l (5.39) for some distance l . If l is suciently small, then with change of variables ~ =~ r 0 ~ r 00 (5.40) ~ = (~ r 0 +~ r 00 )=2 (5.41) we can make the approximation~ r 0 =~ r 00 =~ in the integrand factors that change slowly with ~ r 0 and ~ r 00 . Thus, (5.36) can be written in the form of a bistatic radar equation as follows hP g (i;j;n)i = 1 4 Z 0 jA t j 2 jA r j 2 R 2 rs R 2 st j( i ;f j )j 2 0 d~ = P t 2 (4) 3 Z 0 G t 0 G r R 2 rs R 2 st j( i ;f j )j 2 d~ (5.42) where the NBRCS is given by 0 (~ ) = jRj 2 q 4 4q 2 z Z exp[ik(R 0 rs +R 0 st R 00 rs R 00 st )j =0 ] exp[iq 0 z 0 dem + iq 00 z 00 dem ]( ~ ) d ~ : (5.43) The last equality in (5.42) follows from (5.17). The support of 0 is known as the glistening zone. 62 We note that if f is an innitely dierentiable function of~ r, then by Taylor series expansion we have f(~ r 0 )f(~ r 00 ) =f ~ + ~ =2 f ~ ~ =2 =r ? f(~ ) ~ +O( 3 ) (5.44) wherer ? denotes the gradient without the z-component, which is not needed since ~ lies in the local x;y-plane about ~ 2 0 . Applying expansion (5.44) to the exponents of the three factors in the integrand of (5.43) and recalling (5.19), we nd 0 (~ ) = jRj 2 q 4 4q 2 z Z exp(i~ q ? ~ ) exp(iq z r ? dem ~ ) D exp(iq z r ? res ~ ) E d ~ (5.45) where we have made the approximation q 0 z =q 00 z =q z as justied by the limited support of from condition (5.39). We have also assumed that the ranges and the heights in the exponents of (5.43) are innitely dierentiable functions of~ r as justied by being a physical quantity. The approximations leading to (5.45) represent a geometric optics limit since the result is independent of wavelength, as will be seen in (5.46). Letting ~ =q z ~ and moving the expectation outside the integral by Fubini's theorem, we nd 0 (~ ) = jRj 2 q 4 4q 4 z *Z exp i ~ q ? q z +r ? dem +r ? res ~ d~ + = jRj 2 q 4 q 4 z r ? res + ~ q ? q z +r ? dem = jRj 2 q 4 q 4 z p rg ~ q ? q z r ? dem (5.46) wherep rg is the probability density function (pdf) ofr ? res at~ . Here, subscript rg denotes residual gradient, i.e.,r ? res . The application of Fubini's theorem above is justied by R jp rg j = 1<1. We note that, if the DEM gradient term r ? dem is set to zero, then (5.46) above reduces to (34) of ZV2000 [72]. If we now assume that res is a Gaussian eld, thenr ? res is also Gaussian by linearity of the gradient operator and by the assumed dierentiability of res . If we further assume thatr ? res is isotropic with zero mean, then its pdf has form p rg (~ ) = 1 2 2 rg exp 2 2 2 rg (5.47) where 2 rg is the variance of each component ofr ? res . In this case, p rg is characterized by the single parameter rg , which depends both on the surface roughness and on the DEM resolution. 63 5.2.3 CYGNSS-Specic Modeling 5.2.3.1 Surface Bistatic Radar Cross Section As described in Section 4.2.2, the CYGNSS L1B calibration converts the DDM from power received from the re ected GPS signal in units of watts to surface BRCS in units of square meters. The corresponding DDM is provided in the L1 variable brcs. This variable has the advantages of benetting from the entire CYGNSS calibration and representing a physical property of a region on the surface. The DEM-based model for corresponding to brcs is derived as follows. The L1B calibration makes the assumption that gain and range loss factors are constant over the DDM glistening zone and over the noncoherent integration period. Applying the same assumption to (5.42), we nd hP g (i;j;n)i = P t 2 G sp t G sp r (4) 3 (R sp rs ) 2 (R sp st ) 2 Z 0 0 j( i ;f j )j 2 d~ : (5.48) Substituting P g from (5.30) into P tot g in (4.4), we nd (i;j) = (4) 3 (R sp rs ) 2 (R sp st ) 2 2 P t G sp t G sp r 1 N N X n=1 P g (i;j;n): (5.49) Finally, taking the ensemble average of (5.49) and applying (5.48), we obtain h(i;j)i = 1 N N X n=1 Z 0 0 j( i ;f j )j 2 d~ : (5.50) We note that, for simulation purposes, the number of DDMs averaged in (5.50) can potentially be reduced without aecting accuracy. 5.2.3.2 Modeling the Woodward Ambiguity Function As described in Section 4.2.3, a good approximation to the mean square WAF for GPS applications is given by j(;f)j 2 = 2 ()S 2 (f) (5.51) where () = maxf0; 1jj= c g (5.52) S(f) = sinc(T i f) (5.53) and c is the chip length [75]. We note that (5.51), (5.52), and (5.53) are consistent with (4.8), (4.9), and (4.10), respectively. Before we can use this model in (5.50), we need to compute i in (5.28) and f j in (5.29) for each bin (i;j). To compute i from CYGNSS L1 data, we write i = 0 +i = sp ( sp 0 ) +i (5.54) 64 where 0 is the delay of bini = 0, is the bin spacing, and sp is the path delay of the specular point reported in the L1 data. In particular, is obtained from L1 variable delay resolution; sp can be computed from rx to sp range and tx to sp range or from sc pos *, sp pos *, and tx pos *; and the dierence sp 0 can be obtained from brcs ddm sp bin delay row which provides the delay bin number with fractional part for the reported specular point. Similarly, for Doppler we can write f j =f 0 +jf D =f sp (f sp f 0 ) +jf D (5.55) where f 0 is the Doppler frequency of bin j = 0, f D is the Doppler bin spacing, and f sp is the Doppler frequency of the specular point reported in the L1 data. Here, f D is obtained from L1 variable dopp resolution, f sp is computed from * pos * and * vel * and the wavelength, and the dierence f sp f 0 can be obtained from brcs ddm sp bin dopp col, which provides the Doppler bin number with fractional part for the reported specular point. We note that the L1 specular point variables in the above can be used regardless of whether the specular point is correct or whether a single well- dened specular re ection even exists. Of importance is only that the specular point variables provide a tie-in between a point in delay-Doppler space and a scatterer in three-dimensional space. 5.3 Validation Case Study As a validation case study of the DEM-based model, we selected a CYGNSS track over the SoilSCAPE site at Tonzi Ranch near Sacramento, CA, USA. The Tonzi Ranch site is located in the foothills of the Sierra Nevada and surrounded by heterogeneous terrain. To the west lies the Sacramento Valley with relatively at topography. To the east, the Sierra Nevada rises to elevations of over 3000 m with rugged terrain and steep drainages, as seen in Fig. 5.2. Funded through the NASA ESTO AIST Program, SoilSCAPE aims to pro- vide local-scale surface-to-depth estimates of soil moisture with optimal sam- pling [118,119]. The USC MiXIL research group has already performed detailed soil and vegetation surveys at Tonzi Ranch in support of cal/val activities for the AirMOSS and SMAP missions. As noted in Section 4.3.2.3, Tonzi Ranch lies at the northern edge of the CYGNSS coverage zone. Delay-Doppler maps over Tonzi Ranch tend to have large incidence angles and low SNR. Nevertheless, as described in the same section, we identied a total of three tracks over Tonzi Ranch from year 2018 as having enough signal to be potentially useable. We selected the rst of these tracks, which is number 943 from channel 3 of FM 3 on January 8, for the validation case study. The GPS source for this track is SVN 52 from Block IIR-M. The reported specular points from samples 65 121.8 121.4 121.0 120.6 120.2 Longitude 38.0 38.4 38.8 Latitude Sample: 77502 77504 77506 Sacramento Valley Foothills Sierra Nevada SRTM Ellipsoidal Height [km] Reported Specular Point In-situ Sensors 0 1 2 3 (a) Topography. 121.8 121.4 121.0 120.6 120.2 Longitude 38.0 38.4 38.8 Latitude Sample: 77502 77504 77506 Sacramento Valley Foothills Sierra Nevada SRTM Slope Relative to Ellipsoid [deg] Reported Specular Point In-situ Sensors 0 5 10 15 20 25 (b) Slopes. Figure 5.2: Maps derived from SRTM data showing validation track and in-situ sensor locations amid heterogeneous terrain. 77502, 77504, and 77506 (zero-based indexing) of this track from version 2.1 data [96, 124] are shown in Fig. 5.2. The incidence angle for these samples is 67:1°. The DEM heights were obtained from the NASA SRTM global 1-arcsecond DEM [125] and rereferenced to the WGS-84 Earth ellipsoid. Gradients were estimated at each SRTM post using a least squares t to a window of 31 31 samples of ellipsoidal height by a plane. The heights and the slopes of the corresponding gradient estimates are shown in Fig. 5.2. To illustrate the WAF, bin (i;j) = (1; 2) of sample 77502 is shown in Fig. 5.3. If the surface were at, then the delay factor 2 would have an ellipse- like shape centered on the reported specular point. However, for sample 77502, 66 we see that the nominal ellipse is greatly distorted by the mountains to the northeast. Stated dierently, scatterers that would lie on the nominal ellipse in the case of a level surface are displaced into a bin of lower delay by the raised terrain. In contrast, the topography has no perceivable eect on the nominally hyperbola-like shape of the Doppler factor S 2 , which is nearly a straight line for this particular data collection geometry. The Fresnel re ection coecientR in (5.46) was calculated using the circular polarization model of [72]. The dielectric constant for this calculation was obtained from the soil dielectric model of Peplinski-Ulaby-Dobson [110, 111]. Soil moisture content was taken to be the median over all nodes at a depth of 5 cm from the SoilSCAPE daily-average le for Tonzi Ranch [119]. Sand and clay fractions and bulk density were obtained from the Soil SURvey GeOgraphic (SSURGO) database [126]. The residual gradient parameter rg in (5.47) was tuned by a process of trial and error to match the DDM shape to the data. We found that a large value caused the response of the DDM to spread, while a small value resulted in a more focused, specular-like response. We set rg = 0:5° for all three samples. Maps of 0 as computed by (5.46) and (5.47) for the three selected samples are shown in Fig. 5.4. Over the relatively at Sacramento Valley to the southwest of the reported specular point, we see that the glistening zone follows a well- dened Gaussian-like shape whose width is determined by the residual gradient statistic rg . In contrast, over the steep terrain of the Sierra Nevada to the northeast of the reported specular point, the glistening zone spreads out and 0 tends to have much smaller values. Thus, as the reported specular point moves into progressively rougher terrain with increasing sample number, the total cross section of the glistening zone tends to decrease. Plots of the BRCS DDMh(i;j)i as computed by the DEM-based model (5.50) are also shown in Fig. 5.4. To reduce computational cost, the noncoherent av- eraging in (5.50) was approximated by the midpoint DDM at n =N=2. Using this, execution time of surface integral over the region shown in the gures with 1-arcsecond spacing for a CYGNSS DDM of size 17 11 was about 1 minute on Macbook Pro computer. Plots of the corresponding CYGNSS L1 variable brcs are shown for comparison. Similarities between the DEM-based model and the corresponding CYGNSS data are apparent. The responses of both the model and the data decrease as the surface slopes increase into the Sierra Nevada. Additionally, the location and shape of the model and data DDMs appear generally similar. The main dierence between the DEM-based model and the CYGNSS data is that the response of the model is 5 dB higher than that of the data. Additionally, the model response drops o more slowly with increasing sample number than the data. 67 121.8 121.4 121.0 120.6 120.2 Longitude 38.0 38.4 38.8 Latitude Sacramento Valley Foothills Sierra Nevada WAF Delay Factor 2 [dB] Reported Specular Point In-situ Sensors 10 8 6 4 2 0 (a) WAF delay factor for bin i =1 showing distortion due to topog- raphy. 121.8 121.4 121.0 120.6 120.2 Longitude 38.0 38.4 38.8 Latitude Sacramento Valley Foothills Sierra Nevada WAF Doppler Factor S 2 [dB] Reported Specular Point In-situ Sensors 10 8 6 4 2 0 (b) WAF Doppler factor for bin j = 2. 121.8 121.4 121.0 120.6 120.2 Longitude 38.0 38.4 38.8 Latitude Sacramento Valley Foothills Sierra Nevada WAF | | 2 = 2 S 2 [dB] Reported Specular Point In-situ Sensors 10 8 6 4 2 0 (c) Total WAF for delay-Doppler bin (i;j) = (1;2) where indexing is relative to center. Figure 5.3: Illustration of a WAF map from sample 77502. 68 121.8 121.4 121.0 120.6 120.2 Longitude 38.0 38.4 38.8 Latitude Sacramento Valley Foothills Sierra Nevada Modeled Normalized BRCS ( 0 ) [dB] Reported Specular Point In-situ Sensors 20 22 24 26 28 30 5 0 5 Doppler Bin 8 4 0 4 8 Delay Bin Modeled BRCS [dBsm] 106 108 110 112 114 116 5 0 5 Doppler Bin 8 4 0 4 8 Delay Bin SNR: 2.77 dB CYGNSS BRCS [dBsm] 101 103 105 107 109 111 (a) Sample 77502. 121.8 121.4 121.0 120.6 120.2 Longitude 38.0 38.4 38.8 Latitude Sacramento Valley Foothills Sierra Nevada Modeled Normalized BRCS ( 0 ) [dB] Reported Specular Point In-situ Sensors 20 22 24 26 28 30 5 0 5 Doppler Bin 8 4 0 4 8 Delay Bin Modeled BRCS [dBsm] 106 108 110 112 114 116 5 0 5 Doppler Bin 8 4 0 4 8 Delay Bin SNR: 1.61 dB CYGNSS BRCS [dBsm] 101 103 105 107 109 111 (b) Sample 77504. 121.8 121.4 121.0 120.6 120.2 Longitude 38.0 38.4 38.8 Latitude Sacramento Valley Foothills Sierra Nevada Modeled Normalized BRCS ( 0 ) [dB] Reported Specular Point In-situ Sensors 20 22 24 26 28 30 5 0 5 Doppler Bin 8 4 0 4 8 Delay Bin Modeled BRCS [dBsm] 106 108 110 112 114 116 5 0 5 Doppler Bin 8 4 0 4 8 Delay Bin SNR: 0.72 dB CYGNSS BRCS [dBsm] 101 103 105 107 109 111 (c) Sample 77506. Figure 5.4: Each subgure corresponds to a sample in the validation track. The left panel of each subgure shows the map of 0 from the DEM-based model. The positions of the reported specular point and the in-situ sensors within the glistening zone are indicated by red and white dots, respectively. The center panel gives the corresponding DDM from the DEM-based model. The right panel provides the corresponding CYGNSS DDM. Delay bin spacing is 74:8 m, and Doppler bin spacing is 500 Hz. 69 5.4 Discussion A theoretical model for the DDM of the mean power scattered from a bare surface with heterogeneous roughness in the strong diuse regime has been derived from rst principles. The model consists of noncoherent averaging (5.33), a bistatic radar equation (5.42), and a normalized BRCS in the GO limit (5.46). The model makes use of surface heights and gradients derived from a DEM and treats the residual height (i.e., DEM error) as a random eld. If the residual height is locally Gaussian, then the pdf of the residual surface gradient in (5.46) is given by (5.47). For CYGNSS, the bistatic radar equation has been converted from power to BRCS in (5.50) to permit matchups with the calibrated L1 variable brcs. Addi- tionally, a technique to obtain the delay and Doppler coordinates of CYGNSS DDMs over land has been proposed in (5.54) and (5.55). A case study of a CYGNSS track over a soil moisture validation site suggests two mechanism by which DDMs are aected by topography. First, as seen in Fig. 5.3, raised terrain can displace scatterers into bins of more negative delay by reducing the path length R rs +R st in (5.28), thereby requiring a corresponding reduction of i to maintain the same value of i in (5.42) or (5.50). Second, as seen in the mountainous regions in the left panels of Fig. 5.4, when gradientsr ? dem are large relative to the residual gradient statistic rg , then the normalized BRCS 0 tends to decrease sincer ? dem frequently falls outside the support of the pdf in (5.46). Additionally, large gradients tend to spread out the glistening zone since then~ q ? =q z r ? dem sometimes falls inside the support of the pdf for larger values of q ? =q z , which correspond to points farther from the nominal specular point. Thus, the DEM-based model predicts that both raised terrain and steep terrain can aect the location, amplitude, and shape of the DDM. Validation results suggest that the DEM-based model indeed predicts the somewhat asymmetrical structure of the Tonzi Ranch DDMs, as seen in the center and left panels of Fig. 5.4. The model also predicts the decrease in amplitude as the glistening zone moves from the at Sacramento Valley into the steep terrain of the Sierra Nevada. Finally, the location of the modeled DDM in delay-Doppler coordinates agrees with the data, providing validation of the method proposed in 5.2.3.2 using L1 specular point variables. We note that this location method has not yet been tested with high elevation data. As observed in the previous section, the modeled response is signicantly higher than the validation data. Additionally, the model drops o more slowly with increasing roughness than the data. These observations suggest a number of areas for further investigation, as discussed in the following section. 70 5.5 Future Directions 5.5.1 Land Cover Measurements from GNSS-R are sensitive not only to soil moisture content and surface roughness but also to land cover, such as vegetation and inland water bodies. For validation sites with land cover, corresponding models need to be included. In particular, vegetation attenuation may help explain why the CYGNSS validation data over Tonzi Ranch are lower than the DEM-based bare-surface model. Attenuation is expected to be higher for the large incidence angles of the Tonzi Ranch data since the optical depth of a vegetation layer increases with incidence angle [13, eq. (11.3)]. Inland water bodies can also have a signicant contribution to DDMs [9] and require special treatment since they violate the assumption of strong diuse scattering used in the derivation of the DEM-based model. In particular, the area surrounding Tonzi Ranch contains multiple ponds and small lakes, as seen in Fig. 4.13. The modeling of such inland water bodies in DDMs is an area of ongoing research. 5.5.2 Calibration and Site Selection The DDMs used for the validation study are subject to several potential calibra- tion errors. First, the large incidence angles corresponding to the high latitude of the validation site put the data at the edge of the receive antenna pattern where attenuation and uncertainty are greatest. Second, quality control ags for the validation samples indicated that attitude error was larger than normal, which could cause a gain error by coupling through the receive antenna pattern. Third, the GPS-IIRM source is subject to power exing. To reduce calibration uncertainty, additional sites with less extreme incidence angles and higher SNR need to be identied for further analysis and validation of the DEM-based model. These could include sites currently being identied by the CYGNSS science team for land application cal/val. Additionally, we anticipate improved calibration of GPS source power and antenna patterns with the release of version 3 of CYGNSS science data. 5.5.3 Geometrical Optics Limit The validity of the GO limit needs to be investigated since land surfaces can have small-scale roughness that is not accounted for in the GO limit of (5.46). Attenuation of quasi-specular scattering due to small-scale roughness may help explain why the CYGNSS validation data are lower than predicted by the DEM-based model. For surfaces with homogeneous roughness, VZ2018 provides more general framework based on the theory of scattering amplitude (SA) [75]. In VZ2018, the normalized BRCS 0 in the bistatic radar equation is not constrained by 71 the GO approximation. Rather, models for 0 that also account for small-scale roughness, such as those based on the SSA, can be applied. The extension of VZ2018 to surfaces with heterogeneous roughness is therefore a topic of interest. 5.5.4 Fluctuations The DEM-based model provides the mean power or the mean BRCS. However, the actual power or actual BRCS will uctuate about the mean as the glistening zone moves with respect to surface irregularities. It is therefore of interest to model and characterize these uctuations and their spatial correlations. Sections 25, 26, and 31 of [69] are potentially applicable in this regard. 5.5.5 Heterogeneities of Non-DEM Parameters To obtain the results of Section 5.3, a number of factors were approximated as constant over the surface integral for convenience, including rg ,R, G r , G t , and range loss factors. Additionally, the medium was approximated as homogeneous with depth (i.e., not layered). The validity of these approximations needs to be investigated. For example, if the DEM gradients are less accurate over rough terrain than over smooth, which seems likely, then we could expect rg to be higher over mountainous regions, such as the Sierra Nevada, than over at regions, such as the Sacramento Valley. This eect could potentially help explain why the data fall o more quickly than the modeled response for the validation track. Additionally, techniques have been developed for upscaling in-situ soil mois- ture measurements [5]. This could improve the accuracy ofR in the model since the glistening zone in Fig. 5.4 extends far to the west of the sensors. If the surface slopes are suciently large, then the glistening zones observed by consecutive DDMs will overlap one another signicantly. This provides the opportunity to apply DDM assimilation [127] to retrieve heterogeneous elds of biogeophysical parameters. 5.5.6 Gradient and Gradient Error Estimation The method of estimating gradients from the DEM still needs to be optimized for model performance. We have found that if we make the window size too small, then the map of 0 becomes noisy due to gradient error, resulting in excessive roughness. On the other hand, if we make the window size too large, then the map of 0 loses resolution. Additionally, an automated method needs to be developed to estimate the parameter rg from the DDM structure. For example, rg could be optimized based on a statistic to which it is sensitive, such as the leading or trailing edge slope. 72 5.5.7 Computational Eciency Optimization of the spatial grid size for surface integral evaluation still needs to be investigated. In particular, computational costs could be reduced by pre-averaging the map of 0 , which is independent of delay-Doppler bin. The surface integral could then be carried out over a grid with resolution coarser than the DEM for each bin. The amount of pre-averaging that could be performed without aecting accuracy may depend on the topography itself. Alternatively, the surface integral could potentially be carried out in a more computationally ecient manner by transforming the integral into delay-Doppler space, expressing the result as a convolution, and evaluating the convolution with fast Fourier transforms. This method is well known for ocean applications, where the transformation to delay-Doppler space is performed using a Jacobian [114]. For applications with complex topography, techniques more sophisticated than the Jacobian may be required. 73 Appendix A Partial Derivatives of Fresnel Re ectivity The purpose of this appendix is to compute the sensitivity@ f =@m v that appears in (3.12) using the soil dielectric model of Peplinski-Ulaby-Dobson [110, 111]. Since we may eventually wish to expand the sensitivity analysis to include other soil parameters in addition to m v , we let a represent any soil parameter and we proceed to calculate @ f =@a. We then replace the generic soil parameter a by m v at the appropriate point in the calculation. The Fresnel power re ectivity for an RHCP source and an LHCP receiver is given by f =jR lr j 2 (A.1) where R lr = R vv +R hh 2 (A.2) R vv = r() cos r() + cos (A.3) R hh = cosr() cos +r() (A.4) r() = p sin 2 (A.5) Here,R rt is the Fresnel re ection coecient of a t-polarized transmitter and an r-polarized receiver. The polarizations of right-hand circular, left-hand circular, vertical, and horizontal are denoted by subscripts r, l, v, and h, respectively. Additionally, is incidence angle and is relative dielectric constant. The re ection coecient convention is consistent with that of [73, eq. (8.18)]. 74 Dierentiating, we nd @ f @a = @ @a jR lr j 2 = 2< R lr @R lr @a (A.6) @R lr @a = 1 2 @R vv @a + @R hh @a (A.7) @R vv @a = (cos)[=r() 2r()] [r() + cos] 2 @ @a (A.8) @R hh @a = cos r()[cos +r()] 2 @ @a (A.9) where<(z) denotes the real part of the complex number z. By the Peplinski-Ulaby-Dobson semiempirical dielectric mixing model [110, 111], the relative complex-valued dielectric constant of a soil is given by = 0 + i 00 m (A.10) where the linear correction for the real part is given by 0 = 1:15 0 m 0:68 (A.11) the real and imaginary parts of the dielectric mixing model are given by 0 m = 1 + b s ( s 1) +m 0 v 0 fw m v 1 = (A.12) 00 m = m 00 v 00 fw 1 = =m 00 = v 00 fw (A.13) respectively, m v is water volume fraction (or volumetric moisture content) of the mixture, b is bulk density in grams per cubic centimeter, s is the specic density of the solid soil particles, is an empirically determined constant, and 0 and 00 are empirically determined soil-type dependent constants given by 0 = 1:2748 0:519S 0:152C (A.14) 00 = 1:337 97 0:603S 0:166C (A.15) and S and C represent the mass fractions of sand and clay, respectively, with 0S;C 1, the dielectric constant of the soil solids is given by s = (1:01 + 0:44 s ) 2 0:062 (A.16) the real and imaginary parts of the dielectric constant of free water are given by 0 fw = w1 + w0 w1 1 + (2f w ) 2 (A.17) 00 fw = 2f w ( w0 w1 ) 1 + (2f w ) 2 + e 2 0 f s b b m v (A.18) 75 respectively, 0 is permittivity of free space, w is relaxation time for water, f is frequency in hertz, w0 is the static dielectric constant for water, w1 is the high-frequency limit of 0 fw , and eective conductivity is given by e = 0:0467 + 0:2204 b 0:4111S + 0:6614C: (A.19) Dierentiating, we nd @ @a = @ 0 @a + i @ 00 m @a = 1:15 @ 0 m @a + i @ 00 m @a (A.20) Finally, taking a =m v , we nd @ 0 m @m v = 0(1) m 0 m 0 1 v 0 fw 1 (A.21) @ 00 m @m v = 00 m 00 1 v 00 fw +m 00 = v @ 00 fw @m v (A.22) = 00 m 00 m v + 1 00 fw @ 00 fw @m v (A.23) @ 00 fw @m v = e 2 0 f s b s m 2 v (A.24) 76 Acronyms ABoVE Arctic Boreal Vulnerability Experiment AIEM advanced integral equation method AirMOSS Airborne Microwave Observatory of Subcanopy and Subsurface AIST Advanced Information Systems Technology ALT active layer thickness AMSR2 Advanced Microwave Scanning Radiometer 2 AMSR-E Advanced Microwave Scanning Radiometer for EOS ASCAT Advanced SCATterometer ASP Advanced Study Program BDS BeiDou Navigation Satellite System BRCS bistatic radar cross section cal/val calibration and validation CCI Climate Change Initiative COSPAR Committee On SPAce Research CU University of Colorado CYGNSS Cyclone Global Navigation Satellite System DDM delay-Doppler map DDMI delay-Doppler mapping instrument DEM digital elevation model DLL delay-locked loop DLR Deutsches Zentrum f ur Luft- und Raumfahrt 77 DOI digital object identier ECM error covariance matrix ECV essential climate variable EOS Earth Observing System ESA European Space Agency ESAS 2017 2017-2027 Decadal Survey for Earth Science and Applications from Space ESTO Earth Science Technology Oce EUMETSAT European Organization for the Exploitation of Meteorological Satellites EV-1 Earth Venture 1 EVM-1 Earth Venture Mission 1 FEF Fraser Experimental Forest FM ight model FMPL Flexible Microwave Payload GAGAN GPS Aided Geo Augmented Navigation GCOM-W1 Global Change Observation Mission-Water 1 GCOS Global Climate Observing System GEROS-ISS GNSS rE ectometry, Radio Occultation, and Scatterometry on board the International Space Station GLONASS GLObal NAvigation Satellite System GNSS global navigation satellite system GNSS-R global navigation satellite system re ectometry GO geometrical optics GPS Global Positioning System C/A coarse acquisition GRSS Geoscience and Remote Sensing Society IEEE Institute of Electrical and Electronics Engineers IF intermediate frequency 78 InVEST In-space Validation of Earth Science Technologies IRNSS Indian Regional Navigation Satellite System JAXA Japanese Aerospace Exploration Agency JPL Jet Propulsion Laboratory JSTARS Journal of Selected Topics in Applied Earth Observations and Remote Sensing KA Kirchho approximation L1 level 1 L1A Level 1A L1B Level 1B LEO low Earth orbit LHCP left-hand circular polarization MAP maximum a posteriori probability MEO medium Earth orbit MetOp Meteorological Operational MiXIL Microwave Systems, Sensors, and Imaging Laboratory ML maximum likelihood MMSE minimum mean square error MUOS Mobile Users Objective System MVU minimum variance unbiased NASA National Aeronautics and Space Administration NBRCS normalized bistatic radar cross section NCAR National Center for Atmospheric Research NOAA National Oceanic and Atmospheric Administration NOC National Oceanography Centre NSF National Science Foundation OSU Ohio State University pdf probability density function 79 PO physical optics PRETTY Passive RE ecTometry and DosimeTrY PRN pseudo-random noise QZSS Quasi-Zenith Satellite System RFI radio frequency interference RHCP right-hand circular polarization RO radio occultation ROSES Research Opportunities in Space and Earth Science RZSM root zone soil moisture SA scattering amplitude SAR synthetic aperture radar SAS Space and Airborne Systems SAVERS Soil And VEgetation Re ection Simulator SGR-ReSI Space GNSS Receiver Remote Sensing Instrument SIR-C Spaceborne Imaging Radar C-Band SMAP Soil Moisture Active Passive SMEX02 Soil Moisture Experiment 2002 SMOS Soil Moisture and Ocean Salinity SNoOPI SigNals of Opportunity P-band Investigation SNR signal-to-noise ratio SoilSCAPE Soil moisture Sensing Controller and oPtimal Estimator SoOp signal of opportunity SoOp-R re ectometry from signals of opportunity SPM small perturbation method SSA small slope approximation SSTL Surrey Satellite Technology Limited SSURGO Soil SURvey GeOgraphic SVN space vehicle number 80 SWE snow water equivalent SwRI Southwest Research Institute TDS-1 TechDemoSat-1 UCAR University Corporation for Atmospheric Research UFO UHF Follow On UHF ultra high frequency UK-DMC United Kingdom Disaster Monitoring Constellation UKSA United Kingdom Space Agency UN United Nations UPC Universidad Polit ecnica Catalunya USC University of Southern California VHF very high frequency VWC vegetation water content WAF Woodward ambiguity function WMO World Meteorological Organization 81 References [1] National Academies of Sciences, Engineering, and Medicine, Thriving on OurChangingPlanet: ADecadalStrategyforEarthObservationfromSpace. Washington, DC: The National Academies Press, 2018, doi:10.17226/24938. [2] E. S. Krueger, T. E. Ochsner, S. M. Quiring, D. M. Engle, J. Carlson, D. Twidwell, and S. D. Fuhlendorf, \Measured soil moisture is a bet- ter predictor of large growing-season wildres than the Keetch{Byram drought index," Soil Sci. Soc. Am. J., vol. 81, no. 3, p. 490, 2017, doi:10.2136/sssaj2017.01.0003. [3] V. U. Zavorotny, S. Gleason, E. Cardellach, and A. Camps, \Tu- torial on remote sensing using GNSS bistatic radar of opportunity," IEEE Geosci. Remote Sens. Mag., vol. 2, no. 4, pp. 8{45, Dec. 2014, doi:10.1109/MGRS.2014.2374220. [4] S. V. Nghiem, C. Zuada, R. Shah, C. C. Chew, S. T. Lowe, A. J. Mannucci, E. Cardellach, G. R. Brakenridge, G. Geller, and A. Rosenqvist, \Wetland monitoring with global navigation satellite system re ectometry," Earth Space Sci., vol. 4, no. 1, pp. 16{39, Jan. 2017, doi:10.1002/2016ea000194. [5] D. Clewley, J. B. Whitcomb, R. Akbar, A. R. Silva, A. Berg, J. R. Adams, T. Caldwell, D. Entekhabi, and M. Moghaddam, \A method for upscaling in situ soil moisture measurements to satellite footprint scale using Random Forests," IEEE J. Sel. Topics Appl. Earth Observ. Remote Sens., vol. 10, no. 6, pp. 2663{2673, Jun. 2017, doi:10.1109/jstars.2017.2690220. [6] H. Kim and V. Lakshmi, \Use of Cyclone Global Navigation Satellite Sys- tem (CyGNSS) observations for estimation of soil moisture," Geophys. Res. Lett., vol. 45, no. 16, pp. 8272{8282, Aug. 2018, doi:10.1029/2018gl078923. [7] C. C. Chew, E. E. Small, and L. Read, \How much closer can GNSS-R hydrology retrievals from CYGNSS get us to `the truth'," presented at the IEEE Specialist Meeting on Re ectometry Using GNSS and Other Signals of Opportunity held in Benevento, Italy, on May 21, 2019. 82 [8] A. Azemati and M. Moghaddam, \Circular-linear polarization sig- natures in bistatic scattering models applied to signals of opportu- nity," in 2017 International Conference on Electromagnetics in Ad- vanced Applications (ICEAA), Vorona, Sep. 2017, pp. 1825{1827, doi:10.1109/ICEAA.2017.8065654. [9] M. M. Al-Khaldi, J. T. Johnson, A. J. O'Brien, A. Balenzano, and F. Mat- tia, \Time-series retrieval of soil moisture using CYGNSS," IEEE Trans. Geosci. Remote Sens., pp. 1{10, 2019, doi:10.1109/tgrs.2018.2890646. [10] C. C. Chew and E. E. Small, \Soil moisture sensing using spaceborne GNSS re ections: Comparison of CYGNSS re ectivity to SMAP soil moisture," Geophys. Res. Lett., vol. 45, Apr. 2018, doi:10.1029/2018gl077905. [11] \UCAR/CU CYGNSS soil moisture product." [Online]. Available: https://data.cosmic.ucar.edu/gnss-r/soilMoisture/cygnss/ [12] M. P. Clarizia, N. Pierdicca, F. Costantini, and N. Floury, \Analysis of CYGNSS data for soil moisture retrieval," IEEE J. Sel. Topics Appl. Earth Observ. Remote Sens., pp. 1{9, 2019, doi:10.1109/jstars.2019.2895510. [13] F. T. Ulaby, D. G. Long, W. J. Blackwell, C. Elachi, A. K. Fung, C. Ruf, K. Sarabandi, H. A. Zebker, and J. V. Zyl, Microwave Radar and Radio- metric Remote Sensing. University of Michigan Press, 2014. [14] A. Belward, M. Bourassa, M. Dowell, S. Briggs, H. Dolman, K. Holmlund, R. Husband, S. Quegan, R. Saunders, A. Simmons, B. Sloyan, D. Stammer, K. Steen, T. Tanhua, M. Verstraete, C. Richter, S. Eggleston, K. Hill, and C. Tassone, The Global Observing System for Climate: Implementation Needs. World Meteorological Organization, 2016, GCOS-195. [Online]. Available: https://public.wmo.int/en/programmes/global-climate-observing-system [15] D. Entekhabi, E. G. Njoku, P. E. O'Neill, K. H. Kellogg, W. T. Crow, W. N. Edelstein, J. K. Entin, S. D. Goodman, T. J. Jackson, J. Johnson, J. Kimball, J. R. Piepmeier, R. D. Koster, N. Martin, K. C. McDonald, M. Moghaddam, S. Moran, R. Reichle, J. C. Shi, M. W. Spencer, S. W. Thurman, L. Tsang, and J. V. Zyl, \The Soil Moisture Active Passive (SMAP) mission," Proc. IEEE, vol. 98, no. 5, pp. 704{716, May 2010, doi:10.1109/jproc.2010.2043918. [16] Y. H. Kerr, P. Waldteufel, J.-P. Wigneron, S. Delwart, F. Cabot, J. Boutin, M.-J. Escorihuela, J. Font, N. Reul, C. Gruhier, S. E. Ju- glea, M. R. Drinkwater, A. Hahne, M. Mart n-Neira, and S. Mecklen- burg, \The SMOS mission: New tool for monitoring key elements of the global water cycle," Proc. IEEE, vol. 98, no. 5, pp. 666{687, May 2010, doi:10.1109/JPROC.2010.2043032. 83 [17] J. Font, A. Camps, A. Borges, M. Mart n-Neira, J. Boutin, N. Reul, Y. H. Kerr, A. Hahne, and S. Mecklenburg, \SMOS: The challenging sea surface salinity measurement from space," Proc. IEEE, vol. 98, no. 5, pp. 649{665, May 2010, doi:10.1109/jproc.2009.2033096. [18] J. Figa-Salda~ na, J. J. Wilson, E. Attema, R. Gelsthorpe, M. R. Drinkwater, and A. Stoelen, \The advanced scatterometer (ASCAT) on the meteo- rological operational (MetOp) platform: A follow on for European wind scatterometers," Canadian Journal of Remote Sensing, vol. 28, no. 3, pp. 404{412, 2002, doi:10.5589/m02-035. [19] L. Brocca, S. Hasenauer, T. Lacava, F. Melone, T. Moramarco, W. Wagner, W. Dorigo, P. Matgen, J. Mart nez-Fern andez, P. Llorens, J. Latron, C. Martin, and M. Bittelli, \Soil moisture estimation through ASCAT and AMSR-E sensors: An intercomparison and validation study across Europe," Remote Sens. Environ., vol. 115, no. 12, pp. 3390{3408, Dec. 2011, doi:10.1016/j.rse.2011.08.003. [20] W. Wagner, S. Hahn, R. Kidd, T. Melzer, Z. Bartalis, S. Hasenauer, J. Figa-Salda~ na, P. de Rosnay, A. Jann, S. Schneider, J. Komma, G. Kubu, K. Brugger, C. Aubrecht, J. Z uger, U. Gangkofner, S. Kienberger, L. Brocca, Y. Wang, G. Bl oschl, J. Eitzinger, and K. Steinnocher, \The ASCAT soil moisture product: A review of its specications, validation results, and emerging applications," Meteorol. Z., vol. 22, no. 1, pp. 5{33, Feb. 2013, doi:10.1127/0941-2948/2013/0399. [21] L. Brocca, W. T. Crow, L. Ciabatta, C. Massari, P. de Rosnay, M. Enenkel, S. Hahn, G. Amarnath, S. Camici, A. Tarpanelli, and W. Wagner, \A review of the applications of ASCAT soil moisture products," IEEE J. Sel. Topics Appl. Earth Observ. Remote Sens., vol. 10, no. 5, pp. 2285{2306, May 2017, doi:10.1109/jstars.2017.2651140. [22] M. Kachi, K. Naoki, M. Hori, and K. Imaoka, \AMSR2 validation results," in 2013 IEEE International Geoscience and Remote Sensing Symposium (IGARSS), Jul. 2013, doi:10.1109/igarss.2013.6721287. [23] K. Imaoka, M. Kachi, H. Fujii, H. Murakami, M. Hori, A. Ono, T. Igarashi, K. Nakagawa, T. Oki, Y. Honda, and H. Shimoda, \Global Change Ob- servation Mission (GCOM) for monitoring carbon, water cycles, and climate change," Proc. IEEE, vol. 98, no. 5, pp. 717{734, May 2010, doi:10.1109/jproc.2009.2036869. [24] H. Kim, R. Parinussa, A. G. Konings, W. Wagner, M. H. Cosh, V. Lak- shmi, M. Zohaib, and M. Choi, \Global-scale assessment and combina- tion of SMAP with ASCAT (active) and AMSR2 (passive) soil mois- ture products," Remote Sens. Environ., vol. 204, pp. 260{275, Jan. 2018, doi:10.1016/j.rse.2017.10.026. 84 [25] D. Geudtner, \Sentinel-1 system overview and performance," in Proc. SPIE, vol. 8528, Nov. 2012, doi:10.1117/12.977485. [26] H. Lievens, R. H. Reichle, Q. Liu, G. J. M. D. Lannoy, R. S. Dunbar, S. B. Kim, N. N. Das, M. Cosh, J. P. Walker, and W. Wagner, \Joint Sentinel-1 and SMAP data assimilation to improve soil moisture esti- mates," Geophys. Res. Lett., vol. 44, no. 12, pp. 6145{6153, Jun. 2017, doi:10.1002/2017gl073904. [27] W. Wagner and W. Dorigo, \ESA CCI SM at a glance," Climate Change Initiative Soil Moisture Newsletter, Mar. 2018. [Online]. Available: http://www.esa-soilmoisture-cci.org [28] E. Chapin, A. Chau, J. Chen, B. Heavey, S. Hensley, Y. Lou, R. Machuzak, and M. Moghaddam, \AirMOSS: An airborne P-band SAR to measure root-zone soil moisture," in 2012 IEEE Radar Conference, May 2012, doi:10.1109/radar.2012.6212227. [29] A. Tabatabaeenejad, M. Burgin, X. Duan, and M. Moghaddam, \P-band radar retrieval of subsurface soil moisture prole as a second-order polyno- mial: First AirMOSS results," IEEE Trans. Geosci. Remote Sens., vol. 53, no. 2, pp. 645{658, Feb. 2015, doi:10.1109/tgrs.2014.2326839. [30] R. H. Chen, A. Tabatabaeenejad, and M. Moghaddam, \A time-series active layer thickness retrieval algorithm using P- and L-band SAR observations," in 2016 IEEE International Geoscience and Remote Sensing Symposium (IGARSS), Jul. 2016, pp. 3672{3675. [31] C. D. Hall and R. A. Cordey, \Multistatic scatterometry," in International Geoscience and Remote Sensing Symposium (IGARSS), Remote Sensing: Moving Toward the 21st Century, 1988, doi:10.1109/igarss.1988.570200. [32] J.-C. Auber, A. Bibaut, and J.-M. Rigal, \Characterization of multipath on land and sea at GPS frequencies," in Proceedings of the 7th International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GPS 1994), Salt Lake City, Utah, Sep. 1994, pp. 1155{1171. [33] M. Mart n-Neira, \A passive re ectometry and interferometry system (PARIS): Application to ocean altimetry," ESA Journal, vol. 17, no. 4, pp. 331{355, 1993. [34] J. L. Garrison, S. J. Katzberg, and M. I. Hill, \Eect of sea roughness on bistatically scattered range coded signals from the Global Positioning System," Geophys. Res. Lett., vol. 25, no. 13, pp. 2257{2260, Jul. 1998, doi:10.1029/98gl51615. [35] B. Lin, S. J. Katzberg, J. L. Garrison, and B. A. Wielicki, \Rela- tionship between GPS signals re ected from sea surfaces and surface winds: Modeling results and comparisons with aircraft measurements," J. Geophys. Res. Oceans, vol. 104, no. C9, pp. 20 713{20 727, Sep. 1999, doi:10.1029/1999jc900176. 85 [36] S. T. Lowe, J. L. LaBrecque, C. Zuada, L. J. Romans, L. E. Young, and G. A. Hajj, \First spaceborne observation of an Earth-re ected GPS signal," Radio Sci., vol. 37, no. 1, pp. 7-1{7-28, Jan. 2002, doi:10.1029/2000rs002539. [37] NSF/NCAR Earth Observing Laboratory C-130 Hercules aircraft, doi:10.5065/D6WM1BG0. [38] D. Masters, \Initial results of land-re ected GPS bistatic radar measure- ments in SMEX02," Remote Sens. Environ., vol. 92, no. 4, pp. 507{520, Sep. 2004, doi:10.1016/j.rse.2004.05.016. [39] S. J. Katzberg, O. Torres, M. S. Grant, and D. Masters, \Utilizing cal- ibrated GPS re ected signals to estimate soil re ectivity and dielectric constant: Results from SMEX02," Remote Sens. Environ., vol. 100, no. 1, pp. 17{28, 2005, doi:https://doi.org/10.1016/j.rse.2005.09.015. [40] S. Gleason, S. Hodgart, Y. Sun, C. Gommenginger, S. Mackin, M. Adjrad, and M. Unwin, \Detection and processing of bistatically re ected GPS signals from low Earth orbit for the purpose of ocean remote sensing," IEEE Trans. Geosci. Remote Sens., vol. 43, no. 6, pp. 1229{1241, Jun. 2005, doi:10.1109/tgrs.2005.845643. [41] S. Gleason, \Remote sensing of ocean, ice and land surfaces using bistati- cally scattered GNSS signals from low Earth orbit," Ph.D. dissertation, Univ. Surrey, Guildford, U.K., 2006. [42] M. P. Clarizia, C. P. Gommenginger, S. T. Gleason, M. A. Srokosz, C. Galdi, and M. D. Bisceglie, \Analysis of GNSS-R delay-Doppler maps from the UK-DMC satellite over the ocean," Geophys. Res. Lett., vol. 36, no. 2, p. L02608 (5 pages), Jan. 2009, doi:10.1029/2008gl036292. [43] S. Gleason, S. Lowe, and V. U. Zavorotny, GNSS Applications and Methods. Artech House, 2009, ch. Remote Sensing with Bistatic GNSS Re ections, pp. 399{436. [44] S. Gleason, \Towards sea ice remote sensing with space detected GPS signals: Demonstration of technical feasibility and initial consistency check using low resolution sea ice information," Remote Sensing, vol. 2, no. 8, pp. 2017{2039, Aug. 2010, doi:10.3390/rs2082017. [45] ||, \Space-based GNSS scatterometry: Ocean wind sensing using an empirically calibrated model," IEEE Trans. Geosci. Remote Sens., vol. 51, no. 9, pp. 4853{4863, Sep. 2013, doi:10.1109/tgrs.2012.2230401. [46] M. Unwin, S. Duncan, P. Jales, P. Blunt, and M. Brenchley, \Implementing GNSS re ectometry in space on the TechDemoSat-1 mission," Proceedings of the Institute of Navigation, pp. 1222{1235, 2014. 86 [47] M. Unwin, P. Jales, J. Tye, C. Gommenginger, G. Foti, and J. Rosello, \Spaceborne GNSS-re ectometry on TechDemoSat-1: Early mission operations and exploitation," IEEE J. Sel. Topics Appl. Earth Observ. Remote Sens., vol. 9, no. 10, pp. 4525{4539, Oct. 2016, doi:10.1109/jstars.2016.2603846. [48] G. Foti, C. Gommenginger, and M. Srokosz, \First spaceborne GNSS- re ectometry observations of hurricanes from the UK TechDemoSat-1 mission," Geophys. Res. Lett., vol. 44, no. 24, pp. 12 358{12 366, Dec. 2017, doi:10.1002/2017gl076166. [49] C. C. Chew, R. Shah, C. Zuada, G. Hajj, D. Masters, and A. J. Mannucci, \Demonstrating soil moisture remote sensing with observations from the UK TechDemoSat-1 satellite mission," Geophys. Res. Lett., vol. 43, no. 7, pp. 3317{3324, Apr. 2016, doi:10.1002/2016gl068189. [50] M. Unwin, J. Rawlinson, L. King, C. Gommenginger, G. Foti, M. Ham- mond, and J. Rosello, \Towards a GNSS re ectometry service from UK TDS-1 satellite," presented at the IEEE Specialist Meeting on Re ectom- etry Using GNSS and Other Signals of Opportunity held in Benevento, Italy, on May 20, 2019. [51] C. S. Ruf, S. Gleason, Z. Jelenak, S. Katzberg, A. Ridley, R. Rose, J. Scher- rer, and V. U. Zavorotny, \The CYGNSS nanosatellite constellation hurri- cane mission," in 2012 IEEE International Geoscience and Remote Sensing Symposium (IGARSS), Jul. 2012, doi:10.1109/igarss.2012.6351600. [52] C. S. Ruf, R. Atlas, P. S. Chang, M. P. Clarizia, J. L. Garrison, S. Gleason, S. J. Katzberg, Z. Jelenak, J. T. Johnson, S. J. Majumdar, A. O'Brien, D. J. Posselt, A. J. Ridley, R. J. Rose, and V. U. Zavorotny, \New ocean winds satellite mission to probe hurricanes and tropical convection," Bull. Amer. Meteor. Soc., vol. 97, no. 3, pp. 385{395, Mar. 2016, doi:10.1175/bams-d- 14-00218.1. [53] D. Masters, S. Esterhuizen, P. Jales, E. Ibrahimi, J. Spark, J. Cappaert, O. Nogu es-Correig, T. Yuasa, V. Nguyen, V. Irisov, T. Duly, L. Tan, and V. Zavorotny, \Design and planning for the rst Spire GNSS-R missions of 2019," presented at the IEEE Specialist Meeting on Re ectometry Using GNSS and Other Signals of Opportunity held in Benevento, Italy, on May 20, 2019. [54] J. A. Ruiz-de-Azua, J. F. Mu~ noz, L. Fern andez, M. Badia, D. Llaver a, C. Diez, A. Aguilella, A. P erez, O. Milian, M. Sobrino, A. Navarro, H. Lle o, O. Balcells, M. Sureda, M. Soria, and A. Camps, \The Flexible Microwave Payload-1: A combined GNSS-re ectometer and L-band radiometer for the 3 Cat-4 CubeSat mission," presented at the IEEE Specialist Meeting on Re ectometry Using GNSS and Other Signals of Opportunity held in Benevento, Italy, on May 20, 2019. 87 [55] J. F. Mu~ noz-Mart n, L. Fern andez, J. Ruiz-de-Azua, and A. Camps, \The Flexible Microwave Payload-2: A SDR-based GNSS-R instrument for CubeSats," presented at the IEEE Specialist Meeting on Re ectometry Using GNSS and Other Signals of Opportunity held in Benevento, Italy, on May 20, 2019. [56] J. Wickert, O. Andersen, J. Bandeiras, L. Bertino, E. Cardellach, A. Camps, N. Catarino, B. Chapron, G. Foti, C. Gommenginger, J. Hatton, P. Hoeg, A. Jaggi, M. Kern, T. Lee, M. Martin-Neira, H. Park, N. Pierdicca, J. Rosello, M. Semmling, C. Shum, C. Zuada, F. Soulat, A. Sousa, and J. Xi, \Innovative sea surface monitoring with GNSS-re ectometry aboard ISS: Overview and recent results from GEROS-ISS," in 2016 IEEE International Geoscience and Remote Sensing Symposium (IGARSS), Jul. 2016, doi:10.1109/igarss.2016.7730465. [57] K. Balidakis, M. Semmling, F. Zus, M. Asgarimehr, G. Michalak, G. Bey- erle, J. Wickert, and H. Schuh, \Kepler observing Earth | A re ectometry concept for ocean altimetry," presented at the IEEE Specialist Meeting on Re ectometry Using GNSS and Other Signals of Opportunity held in Benevento, Italy, on May 20, 2019. [58] G. Giorgi, T. D. Schmidt, C. Trainotti, R. Mata-Calvo, C. Fuchs, M. M. Hoque, J. Berdermann, J. Furthner, C. G unther, T. Schuldt, J. Sanjuan, M. Gohlke, M. Oswald, C. Braxmaier, K. Balidakis, G. Dick, F. Flechtner, M. Ge, S. Glaser, R. K onig, G. Michalak, M. Murb ock, M. Semmling, and H. Schuh, \Advanced technologies for satellite navigation and geodesy," Advances in Space Research, Jun. 2019, doi:10.1016/j.asr.2019.06.010. [59] R. Shah, J. L. Garrison, and M. S. Grant, \Demonstration of bistatic radar for ocean remote sensing using communication satellite signals," IEEE Geosci. Remote Sens. Lett., vol. 9, no. 4, pp. 619{623, Jul. 2012, doi:10.1109/lgrs.2011.2177061. [60] R. Shah and J. L. Garrison, \Correlation properties of direct broad- cast signals for bistatic remote sensing," in 2012 IEEE International Geoscience and Remote Sensing Symposium (IGARSS), Jul. 2012, doi:10.1109/igarss.2012.6351895. [61] ||, \Application of the ICF coherence time method for ocean remote sensing using digital communication satellite signals," IEEE J. Sel. Topics Appl. Earth Observ. Remote Sens., vol. 7, no. 5, pp. 1584{1591, May 2014, doi:10.1109/jstars.2014.2314531. [62] R. Shah, X. Xu, S. Yueh, C. S. Chae, K. Elder, B. Starr, and Y. Kim, \Remote sensing of snow water equivalent using P-band coherent re ection," IEEE Geosci. Remote Sens. Lett., vol. 14, no. 3, pp. 309{313, Mar. 2017, doi:10.1109/lgrs.2016.2636664. 88 [63] D. Wadsworth, \Military communications satellite system multiplies UHF channel capacity for mobile users," in 1999 IEEE Military Communications (MILCOM) Conference Proceedings (Cat. No.99CH36341), vol. 2, Oct. 1999, pp. 1145{1152, doi:10.1109/MILCOM.1999.821382. [64] T. Freeman and R. Shah, \Designing a spacecraft to ex- ploit signals-of-opportunity for Earth system science," Jul. 2018, presentation at COSPAR in Pasadena, CA. [Online]. Avail- able: https://www.researchgate.net/publication/326723972 Designing a spacecraft to exploit Signals-of-Opportunity for Earth System Science [65] S. H. Yueh, X. Xu, R. Shah, Y. Kim, J. L. Garrison, A. Komanduru, and K. Elder, \Remote sensing of snow water equivalent using coherent re ection from satellite signals of opportunity: Theoretical modeling," IEEE J. Sel. Topics Appl. Earth Observ. Remote Sens., vol. 10, no. 12, pp. 5529{5540, Dec. 2017, doi:10.1109/jstars.2017.2743172. [66] K. McCaney, \Researchers take high-bandwidth communications to the South Pole," Aug 24, 2015. [Online]. Available: https://defensesystems. com/articles/2015/08/24/afrl-navy-lockheed-muos-antarctica-tests.aspx [67] M. A. Isakovich, \The scattering of waves from a statistically rough surface," Zhurn. Eksp. Teor. Fiz., vol. 23, pp. 305{314, 1952, in Russian. [68] S. M. Rytov, Y. A. Kravtsov, and V. I. Tatarskii, Principles of Statistical Radiophysics 4: Wave Propagation Through Random Media. Springer- Verlag, 1989. [69] F. G. Bass and I. M. Fuks, Wave Scattering from Statistically Rough Surfaces. Pergamon Press, 1979, doi:10.1016/C2013-0-05724-6. [70] P. Beckmann and A. Spizzichino, The Scattering of Electromagnetic Waves from Rough Surfaces. New York: Pergamon Press, 1963. [71] S. J. Katzberg and J. L. Garrison Jr., \Utilizing GPS to determine iono- spheric delay over the ocean," NASA Langley Research Center, Hampton, VA 23681-0001, Technical Memorandum 4750, December 1996. [72] V. U. Zavorotny and A. G. Voronovich, \Scattering of GPS signals from the ocean with wind remote sensing application," IEEE Trans. Geosci. Remote Sens., vol. 38, no. 2, pp. 951{964, Mar. 2000, doi:10.1109/36.841977. [73] C. Ruf, P. S. Chang, M.-P. Clarizia, S. Gleason, Z. Jelenak, S. Majumdar, M. Morris, J. Murray, S. Musko, D. Posselt, D. Provost, D. Starkenburg, and V. U. Zavorotny, CYGNSS Handbook. Michigan Publishing, Apr. 2016. [Online]. Available: http://clasp-research.engin.umich.edu/missions/ cygnss/reference/cygnss-mission/CYGNSS Handbook April2016.pdf 89 [74] A. G. Voronovich and V. U. Zavorotny, \The transition from weak to strong diuse radar bistatic scattering from rough ocean surface," IEEE Trans. Antennas Propag., vol. 65, no. 11, pp. 6029{6034, Nov. 2017, doi:10.1109/tap.2017.2752219. [75] ||, \Bistatic radar equation for signals of opportunity revisited," IEEE Trans. Geosci. Remote Sens., vol. 56, no. 4, pp. 1959{1968, Apr. 2018, doi:10.1109/tgrs.2017.2771253. [76] B. J. Choudhury, T. J. Schmugge, A. Chang, and R. W. Newton, \Eect of surface roughness on the microwave emission from soils," J. Geophys. Res., vol. 84, no. C9, p. 5699, 1979, doi:10.1029/jc084ic09p05699. [77] A. Fung and H. Eom, \Coherent scattering of a spherical wave from an irregular surface," IEEE Trans. Antennas Propag., vol. 31, no. 1, pp. 68{72, Jan. 1983, doi:10.1109/TAP.1983.1142979. [78] R. D. De Roo and F. T. Ulaby, \Bistatic specular scattering from rough dielectric surfaces," IEEE Trans. Antennas Propag., vol. 42, no. 2, pp. 220{231, 1994, doi:10.1109/8.277216. [79] A. Ishimaru, Wave propagation and scattering in random media. Academic press New York, 1978, vol. 2. [80] F. T. Ulaby, R. K. Moore, and A. K. Fung, Microwave Remote Sensing: Active and Passive, Volume II. Artech House, 1982. [81] N. Pierdicca, L. Guerriero, R. Giusto, M. Brogioni, and A. Egido, \SAVERS: A simulator of GNSS re ections from bare and vegetated soils," IEEE Trans. Geosci. Remote Sens., vol. 52, no. 10, pp. 6542{6554, Oct. 2014, doi:10.1109/TGRS.2013.2297572. [82] L. Thirion-Lefevre, E. Colin-Koeniguer, and C. Dahon, \Bistatic scatter- ing from forest components. Part I: coherent polarimetric modelling and analysis of simulated results," Waves Random Complex Medium, vol. 20, no. 1, pp. 36{61, Feb. 2010, doi:10.1080/17455030903499680. [83] M. Kurum, M. Deshpande, A. T. Joseph, P. E. O'Neill, R. H. Lang, and O. Eroglu, \SCoBi-Veg: A generalized bistatic scattering model of re ectometry from vegetation for signals of opportunity applications," IEEE Trans. Geosci. Remote Sens., vol. 57, no. 2, pp. 1049{1068, Feb. 2019, doi:10.1109/tgrs.2018.2864631. [84] J. Johnson, K. Warnick, and P. Xu, \On the geometrical optics (Hagfors' law) and physical optics approximations for scattering from exponentially correlated surfaces," IEEE Trans. Geosci. Remote Sens., vol. 45, no. 8, pp. 2619{2629, Aug. 2007, doi:10.1109/tgrs.2007.900682. 90 [85] M. Saillard, R. D. De Roo, and F. T. Ulaby, \Comments on `Bistatic specu- lar scattering from rough dielectric surfaces' [with reply]," IEEE Trans. An- tennas Propag., vol. 43, no. 2, pp. 224{226, Feb. 1995, doi:10.1109/8.366388. [86] S. O. Rice, \Re ection of electromagnetic waves from slightly rough sur- faces," Commun. Pure Appl. Math., vol. 4, no. 2-3, pp. 351{378, Aug. 1951, doi:10.1002/cpa.3160040206. [87] L. Tsang, J. A. Kong, and R. T. Shin, Theory of Microwave Remote Sensing. Wiley, 1985. [88] T.-D. Wu and K.-S. Chen, \A reappraisal of the validity of the IEM model for backscattering from rough surfaces," IEEE Trans. Geosci. Remote Sens., vol. 42, no. 4, pp. 743{753, Apr. 2004, doi:10.1109/tgrs.2003.815405. [89] N. Pierdicca, L. Guerriero, D. Comite, F. Costantini, and L. Dente, \Ex- perimental and modelled GNSS re ectometry response over land," in 2018 International Conference on Electromagnetics in Advanced Applications (ICEAA), Sep. 2018, doi:10.1109/iceaa.2018.8520343. [90] L. Guerriero, L. Dente, D. Comite, and N. Pierdicca, \Simulations of spaceborne GNSS-R signal over land," presented at the IEEE Specialist Meeting on Re ectometry Using GNSS and Other Signals of Opportunity held in Benevento, Italy, on May 20, 2019. [91] A. G. Voronovich, \Small-slope approximation for electromagnetic wave scattering at a rough interface of two dielectric half-spaces," Waves in Random Media, vol. 4, no. 3, pp. 337{367, Jul. 1994, doi:10.1088/0959- 7174/4/3/008. [92] ||, Wave Scattering from Rough Surfaces, 2nd ed. Springer, 1999, doi:10.1007/978-3-642-59936-1. [93] V. U. Zavorotny and A. G. Voronovich, \Bistatic GPS signal re ections at various polarizations from rough land surface with moisture content," in IEEE International Geoscience and Remote Sensing Symposium Proceed- ings, vol. 7, 2000, pp. 2852{2854, doi:10.1109/IGARSS.2000.860269. [94] J. T. Johnson and J. D. Ouellette, \Polarization features in bistatic scat- tering from rough surfaces," IEEE Trans. Geosci. Remote Sens., vol. 52, no. 3, pp. 1616{1626, Mar. 2014, doi:10.1109/tgrs.2013.2252909. [95] D. Comite and N. Pierdicca, \Monostatic and bistatic scattering modeling of the anisotropic rough soil," IEEE Trans. Geosci. Remote Sens., vol. 57, no. 5, pp. 2543{2556, May 2019, doi:10.1109/tgrs.2018.2874540. [96] S. Gleason, C. S. Ruf, A. J. O'Brien, and D. S. McKague, \The CYGNSS level-1 calibration algorithm and error analysis based on on-orbit measure- ments," IEEE J. Sel. Topics Appl. Earth Observ. Remote Sens., vol. 12, no. 1, pp. 37{49, Jan. 2019, doi:10.1109/jstars.2018.2832981. 91 [97] F. Said, Z. Jelenak, J. Park, S. Soisuvarn, and P. S. Chang, \Ocean calibration of CYGNSS measurements using wind and wave dependent GMF," presented at the CYGNSS science team meeting held in Ann Arbor, MI, USA, on Jun. 5, 2019. [98] Z. Jelenak, \Implementing track-wise debiasing for CYGNSS L2 data," presented to the CYGNSS science team by videoconference on Jul. 10, 2019. [99] N. Pierdicca, L. Pulvirenti, F. Ticconi, and M. Brogioni, \Radar bistatic congurations for soil moisture retrieval: A simulation study," IEEE Trans. Geosci. Remote Sens., vol. 46, no. 10, pp. 3252{3264, Oct. 2008, doi:10.1109/tgrs.2008.921495. [100] P. Ferrazzoli, L. Guerriero, N. Pierdicca, and R. Rahmoune, \Forest biomass monitoring with GNSS-R: Theoretical simulations," Adv. Space Res., vol. 47, no. 10, pp. 1823{1832, May 2011, doi:10.1016/j.asr.2010.04.025. [101] H. Carre~ no-Luengo, G. Luzi, and M. Crosetto, \Sensitivity of CyGNSS bistatic re ectivity and SMAP microwave radiometry brightness tempera- ture to geophysical parameters over land surfaces," IEEE J. Sel. Topics Appl. Earth Observ. Remote Sens., vol. 12, no. 1, pp. 107{122, Jan. 2019, doi:10.1109/jstars.2018.2856588. [102] A. Camps, M. Vallllossera, H. Park, G. Portal, and L. Rossato, \Sensi- tivity of TDS-1 GNSS-R re ectivity to soil moisture: Global and regional dierences and impact of dierent spatial scales," Remote Sensing, vol. 10, no. 11, p. 1856, Nov. 2018, doi:10.3390/rs10111856. [103] A. Camps, H. Park, M. Pablos, G. Foti, C. P. Gommenginger, P.- W. Liu, and J. Judge, \Sensitivity of GNSS-R spaceborne observa- tions to soil moisture and vegetation," IEEE J. Sel. Topics Appl. Earth Observ. Remote Sens., vol. 9, no. 10, pp. 4730{4742, Oct. 2016, doi:10.1109/jstars.2016.2588467. [104] S. Hensley, S. Oveisgharan, S. Saatchi, M. Simard, R. Ahmed, and Z. Had- dad, \An error model for biomass estimates derived from polarimetric radar backscatter," IEEE Trans. Geosci. Remote Sens., vol. 52, no. 7, pp. 4065{4082, Jul. 2014, doi:10.1109/tgrs.2013.2279400. [105] J. D. Campbell and M. Moghaddam, \GNSS-R parameter sensitiv- ities for soil moisture retrieval," in 2018 International Conference on Electromagnetics in Advanced Applications (ICEAA), Sep. 2018, doi:10.1109/iceaa.2018.8520512. [106] S. Yueh, A. Colliander, R. Shah, and X. Xu, \CYGNSS soil moisture algorithm and validation," presented to the CYGNSS science team by videoconference on Apr. 11, 2018. 92 [107] T. J. Jackson and T. J. Schmugge, \Vegetation eects on the microwave emission of soils," Remote Sens. Environ., vol. 36, no. 3, pp. 203{212, Jun. 1991, doi:10.1016/0034-4257(91)90057-d. [108] S.-B. Kim, L. Tsang, J. T. Johnson, S. Huang, J. J. van Zyl, and E. G. Njoku, \Soil moisture retrieval using time-series radar observations over bare surfaces," IEEE Trans. Geosci. Remote Sens., vol. 50, no. 5, pp. 1853{1863, May 2012, doi:10.1109/tgrs.2011.2169454. [109] P. Swerling, \Topics in generalized least squares signal estimation," SIAM J. Appl. Math., vol. 14, no. 5, pp. 998{1031, 1966, doi:10.1137/0114082. [110] N. R. Peplinski, F. T. Ulaby, and M. C. Dobson, \Dielectric properties of soils in the 0.3{1.3-GHz range," IEEE Trans. Geosci. Remote Sens., vol. 33, no. 3, pp. 803{807, May 1995, doi:10.1109/36.387598. [111] ||, \Corrections to dielectric properties of soils in the 0.3{1.3-GHz range," IEEE Trans. Geosci. Remote Sens., vol. 33, no. 6, p. 1340, Nov. 1995, doi:10.1109/tgrs.1995.477193. [112] M. Burgin, D. Clewley, R. M. Lucas, and M. Moghaddam, \A generalized radar backscattering model based on wave theory for multilayer multi- species vegetation," IEEE Trans. Geosci. Remote Sens., vol. 49, no. 12, pp. 4832{4845, Dec. 2011, doi:10.1109/tgrs.2011.2172949. [113] A. Azemati and M. Moghaddam, \Modeling and analysis of bistatic scattering from forests in support of soil moisture retrieval," in 2017 IEEE International Symposium on Antennas and Propagation and USNC/URSI National Radio Science Meeting, Jul. 2017, pp. 1833{1834, doi:10.1109/apusncursinrsm.2017.8072959. [114] J. March an-Hern andez, A. Camps, N. Rodr guez- Alvarez, E. Valencia, X. Bosch-Lluis, and I. Ramos-P erez, \An ecient algorithm to the simula- tion of delay-Doppler maps of re ected global navigation satellite system signals," IEEE Trans. Geosci. Remote Sens., vol. 47, no. 8, pp. 2733{2740, Aug. 2009, doi:10.1109/tgrs.2009.2014465. [115] G. A. Hajj and C. Zuada, \Theoretical description of a bistatic system for ocean altimetry using the GPS signal," Radio Sci., vol. 38, no. 5, pp. 10-1{10-19, Oct. 2003, doi:10.1029/2002rs002787. [116] K. Sullivan and D. O'Callaghan, \Tourists and feral animals ock to Lake Eyre as oodwaters end 1,000 km journey," May 2018. [Online]. Avail- able: https://www.abc.net.au/news/rural/programs/sa-country-hour/ 2018-05-24/lake-eyre-lls-water-tourists-feral-animals/9792352 [117] \Algorithm theoretical basis document level 1B DDM calibration," Aug. 2018, University of Michigan document number 148-0137. [Online]. Avail- able: http://clasp-research.engin.umich.edu/missions/cygnss/reference/ 148-0137 ATBD L1B DDMCalibration Rev2 Aug2018 release.pdf 93 [118] M. Moghaddam, D. Entekhabi, Y. Goykhman, K. Li, M. Liu, A. Mahajan, A. Nayyar, D. Shuman, and D. Teneketzis, \A wireless soil moisture smart sensor web using physics-based optimal control: Concept and initial demonstrations," IEEE J. Sel. Topics Appl. Earth Observ. Remote Sens., vol. 3, no. 4, pp. 522{535, Dec. 2010, doi:10.1109/jstars.2010.2052918. [119] M. Moghaddam, A. Silva, D. Clewley, R. Akbar, S. A. Hussaini, J. Whit- comb, R. Devarakonda, R. Shrestha, R. B. Cook, G. Prakash, S. K. Santhana Vannan, and A. G. Boyer, \Soil moisture proles and tempera- ture data from SoilSCAPE sites, USA," 2016, ORNL DAAC, Oak Ridge, Tennessee, USA, doi:10.3334/ornldaac/1339. [120] W. Gu, H. Xu, and L. Tsang, \A numerical Kirchho simulator for GNSS- R land applications," Progress In Electromagnetics Research, vol. 164, pp. 119{133, 2019, doi:10.2528/pier18121803. [121] J. Zhu, L. Tsang, and H. Xu, \A physical patch model for GNSS-R land applications," Progress In Electromagnetics Research, vol. 164, pp. 93{105, 2019. [122] A. Ishimaru, Electromagnetic Wave Propagation, Radiation, and Scattering. Wiley, May 2017, doi:10.1002/9781119079699. [123] W. C. Chew, Waves and Fields in Inhomogeneous Media. IEEE Press, 1995. [124] CYGNSS, \CYGNSS level 1 science data record version 2.1," 2018, NASA Physical Oceanography DAAC, doi:10.5067/cygns-l1x21. [125] Earth Resources Observation and Science Center, United States Geological Survey, \Shuttle Radar Topography Mission (SRTM) 1 arcsecond global," 2017, doi:10.5066/f7pr7tft. [126] Soil Survey Sta, Natural Resources Conservation Service, United States Department of Agriculture, \Soil Survey Geographic Database (SSURGO)," 2015, doi:10.15482/usda.adc/1242479. [127] F. Huang, J. L. Garrison, M. Leidner, B. Annane, and R. Homan, \A GNSS-R forward model for delay-Doppler map assimilation," in 2018 IEEE International Geoscience and Remote Sensing Symposium (IGARSS), Jul. 2018, doi:10.1109/igarss.2018.8518987. 94

Abstract (if available)

Abstract Global land surface observations of soil moisture from space are key to understanding water, energy, and carbon cycles. As such, soil moisture observations support many applications of human interest, including weather and climate forecasting, flood and landslide prediction, aquifer modeling, drought analysis, crop productivity evaluation, wildfire danger assessment, and human health. Current soil moisture observation missions rely on large, expensive Earth-observing instruments, such as radiometers and scatterometers, which tend to lack spatial and temporal resolution. ❧ An alternative approach for Earth observation is the use of reflectometry from signals of opportunity (SoOp-R). In particular, SoOps from global navigation satellite systems and digital communication satellites are continually reflected from Earth's surface, picking up information about soil moisture and other biogeophysical parameters in the process. Due to their low weight and power, SoOp-R receivers can fly on small satellites with the potential for improved resolution and revisit rates at a lower cost. The aim of the present research is to develop physics-based models for satellite remote sensing of land surfaces using SoOp-R. ❧ We begin by investigating the theoretical accuracy of a time series method for soil moisture retrieval without relying on any non-SoOp-R soil moisture products to calibrate each resolution cell. The method exploits the incidence angle diversity that is inherent in SoOp-R data. The method assumes coherent scattering of electromagnetic waves in the presence of instrument calibration errors, vegetation cover uncertainty, and surface roughness uncertainty. We conclude that the SoOp-R-only retrieval problem for current missions, which have only a single frequency and a single polarization, is poorly conditioned in general but may be feasible in the case of a bare rough surface. ❧ Next, we add a noncoherent scattering model based on the small perturbation method (SPM) and we attempt to validate the models with satellite data. Initial results from a large dry lakebed using the coherent model are promising. However, we find that contribution of the noncoherent SPM is always negligible relative to that of the coherent model within their common regime of validity. Furthermore, for validation sites having in-situ sensors, we find that both models fail due to topography within the regions of active scattering. ❧ Finally, we develop a new model for the scattering of electromagnetic waves from land surfaces with topography. The model makes use of a digital elevation model and represents an extension of a well-known existing model for scattering from a rough sea surface. Initial validation results are encouraging, and areas for future work are identified.

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Physics-based bistatic radar scattering model for vegetated terrains in support of soil moisture retrieval from signals of opportunity

PDF

Physics-based models and microwave sensors for mapping soil carbon in arctic permafrost using long-wavelength radars

PDF

Multi-constrained inversion algorithms for microwave imaging

Asset Metadata

Creator Campbell, James David (author)

Core Title Electromagnetic scattering models for satellite remote sensing of soil moisture using reflectometry from microwave signals of opportunity

School Viterbi School of Engineering

Degree Doctor of Philosophy

Degree Program Electrical Engineering

Publication Date 10/28/2019

Defense Date 09/04/2019

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag Cyclone Global Navigation Satellite System (CYGNSS),delay-Doppler map (DDM),digital elevation model (DEM),electromagnetic reflection,electromagnetic scattering by rough surfaces,global navigation satellite system-reflectometry (GNSS-R),Global Positioning System (GPS),OAI-PMH Harvest,reflectometry,remote sensing,signals of opportunity-reflectometry (SoOp-R),soil measurements,soil moisture,Soil moisture Sensing Controller and oPtimal Estimator (SoilSCAPE)

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Moghaddam, Mahta (committee chair), Sanders, Kelly (committee member), Tabatabaeenejad, Alireza (committee member)

Creator Email james@jamesdcampbell.net,jamesdca@usc.edu

Permanent Link (DOI) https://doi.org/10.25549/usctheses-c89-227944

Unique identifier UC11673097

Identifier etd-CampbellJa-7887.pdf (filename),usctheses-c89-227944 (legacy record id)

Legacy Identifier etd-CampbellJa-7887.pdf

Dmrecord 227944

Document Type Dissertation

Rights Campbell, James David

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the a...

Repository Name University of Southern California Digital Library

Repository Location USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA