Qi Gao
email: qi.gao@isardsat.cat

Mehrez Zribi
email: mehrez.zribi@cesbio.cnes.fr

Maria Jose Escorihuela
email: mj.escorihuela@isardsat.cat

Nicolas Baghdadi
email: nicolas.baghdadi@teledetection.fr

Synergetic Use of Sentinel-1 and Sentinel-2 Data for Soil Moisture Mapping at 100 m Resolution

Keywords: Water resources, remote sensing, synthetic aperture radar (SAR), soil moisture, altimetry, Sentinel Recursos hídrics, teledetecció, radar d'obertura sintètica (SAR), humitat del sòl, altimetria

The estimation of the water resources of the continental surfaces at a regional and global scale is fundamental for good water resources management. This estimation covers a wide range of topics and fields, including the characterisation of soils and water resources at the basin scale, hydrological modelling and flood prediction and mapping. In this context, the characterisation of the states of the continental surface, to obtain better input parameters for hydrological models, is essential to improve the precision in the simulation of flows, droughts, and floods. The estimation of the water content in the system, including the different water bodies and the free water in the soil, is especially necessary for a precise description of the hydrological processes and, in general, of the water cycle on the continental surfaces. To better characterise hydrological processes, human interventions cannot be neglected. Humans influence the water cycle, mainly through irrigation and the construction of reservoirs, which must be correctly quantified.

The objective of the thesis is the improvement of the remote estimation of water resources, including the quantification of human factors, using several sensors recently launched, taking advantage of recent developments in remote sensing technology. With the arrival of the Sentinel constellations (Sentinel-1, 2, 3), we have better tools to estimate water resources, including human impacts, with greater precision and coverage.

This thesis consists mainly of two parts where human interventions in the water cycle are considered: irrigation cartography (as an application of soil moisture), and the forcing of reservoirs in hydrological simulations (as an application of altimetry).

Firstly, soil moisture is estimated from the statistical analysis of Sentinel-1 SAR data. Two methodologies are developed to obtain soil moisture with a spatial resolution of 100 m based on the interpretation of Sentinel-1 data collected with the VV polarization (vertical-vertical), which is combined with optical data of Sentinel-2 for the analysis of the effects of vegetation.

Secondly, irrigation is mapped under various meteorological conditions, including high spatial and temporal resolution. A methodology for irrigation mapping is proposed using SAR data obtained in VV (vertical-vertical) and VH (vertical-horizontal) polarizations. With Sentinel-1 time series, different statistics and metrics are analysed, including the mean value, the variance of the signal, the correlation length and the fractal dimension, based on which the classification of irrigated trees, irrigated crops, and non-irrigated crops are derived.

Finally, the level of the reservoirs is estimated from the Sentinel-3 altimetry data, with the SAR altimeter (SRAL), based on different algorithms to improve the accuracy. This study presents three specialised algorithms or retrackers designed to obtain the level of the surface of the studied inland bodies of water, minimising the contamination of the waveforms due to the surrounding soil. The performance of the selection method of the proposed wave portion is compared with three retrackers, that is, the centre of gravity retracker (OCOG) and the two-step physical-based retracker. Temporal series of the water level of reservoirs located in the basin of the Ebro River (Spain) are obtained.

As an application, the level series of the reservoirs obtained are used to force the reservoirs in hydrological simulations.

all my PhD study. I would specially thank Maria Jose, who offered me an internship opportunity after my master study followed by the PhD opportunity and who offered generous support and encouragement in my study. Also, I especially thank Mehrez, who dedicated a lot of time in guiding my study with invaluable ideas and showed me the way to do scientific research. I would like to thank Pere as well for his support and patience to take me as his first PhD student.

Furthermore, I would like to acknowledge the jury members and the members of my thesis committee, who offered valuable suggestions and help for my PhD study.

I would also like to acknowledge with much appreciation the isardSAT team with whom I have been working with more than three years. I am very thankful for the opportunity isardSAT offered me to be a part of the company and the financial support which makes my PhD study possible. A special thanks to Eduard Makhoul who generously offered his knowledge to help me in the altimetry research. Another special thanks to Mercè Natividad who saved me from countless administrative processes.

I also thank the CESBIO staff, who offered me a very friendly and helpful working environment during my stay in Toulouse and the Observatori de l'Ebre team and Universitat Ramon Llull team who offered great support during my study, particularly Anaïs Barella, who helped me with the SASER model.

I also sincerely appreciate the funding source I received apart from isardSAT offered, including the grant DI-15-08105 from the Spanish Education Ministry (MICINN), with the support of the Secretariat of Universities and Research of the Business and Knowledge Department of the RESUMEN La estimación de los recursos hídricos de las superficies continentales a escala regional y global es fundamental para una buena gestión de los recursos hídricos. Esta estimación cubre una amplia gama de temas y campos, incluyendo la caracterización de los suelos y de los recursos hídricos a escala de cuenca, la modelización hidrológica y la predicción y la cartografía de inundaciones. En este contexto, la caracterización de los estados de la superficie continental, para obtener mejores parámetros de entrada para los modelos hidrológicos, es esencial para mejorar la precisión en la simulación de caudales, sequías e inundaciones. La estimación del contenido de agua en el sistema, incluidas las diferentes masas de agua y el agua libre en el suelo, es especialmente necesaria para una descripción precisa de los procesos hidrológicos y, en general, del ciclo del agua en las superficies continentales. Una caracterización precisa de los procesos hidrológicos requiere no descuidar las intervenciones humanas. El hombre influye en el ciclo del agua, principalmente mediante el riego y la construcción de embalses, lo que se debe cuantificar correctamente.

El objetivo de la tesis es la mejora de la estimación remota de los recursos hídricos, incluyendo la cuantificación de los factores humanos, mediante el uso de varios sensores lanzados recientemente, aprovechando recientes desarrollos en la tecnología de teledetección. Con la llegada de las constelaciones Sentinel (Sentinel-1, 2, 3), disponemos de mejores herramientas para estimar los recursos hídricos, incluyendo los impactos humanos, con una mayor precisión y cobertura.

Este trabajo de tesis consta principalmente en dos ejes de investigación donde se estiman las intervenciones humanas en el ciclo hidrológico: la cartografía del riego (como aplicación en humedad del suelo), y el forzamiento de embalses en simulaciones hidrológicas (como aplicación de la altimetría).

En relación al primer eje, se estima la humedad del suelo a partir del análisis estadístico de los datos SAR de Sentinel-1. Se desarrollan dos metodologías para obtener la humedad del suelo con una resolución espacial de 100 m basándose en la interpretación de los datos de Sentinel-1 obtenidas con la polarización VV (vertical-vertical), que se combina con datos ópticas Sentinel-2 para el análisis de los efectos de la vegetación.

Como aplicación de la humedad del suelo, se cartografía el riego en diversas condiciones meteorológicas, y con una alta resolución espacial y temporal. Se propone una metodología para la cartografía del riego mediante datos SAR obtenidos en polarizaciones VV (verticalvertical) y VH (vertical-horizontal). A partir de la serie temporal Sentinel-1, se analizan diferentes estadísticas y métricas, incluyendo el valor medio, la varianza de la señal, la longitud de la correlación y la dimensión fractal, a partir de los cuales se clasifican los árboles irrigados, los cultivos irrigados y los cultivos no irrigados.

En el segundo eje, se estima el nivel de los embalses a partir de los datos de altimetría de Sentinel-3, con el altímetro SAR (SRAL), basándose en diferentes algoritmos para mejorar la precisión. Este estudio presenta tres algoritmos especializados o retrackers destinados a obtener el nivel de la superficie de los cuerpos de agua estudiados, minimizando la contaminación de las formas de onda debido al suelo que los rodea. Se compara el rendimiento del método propuesto de selección de la porción de onda con tres retrackers, es decir, un retracker de umbral, el retracker xi xii RESUMEN del centro de gravedad (OCOG) y un retracker de base física de dos pasos. Se obtienen series temporales del nivel de la lámina de agua de embalses situados en la cuenca del río Ebro (España).

Como aplicación, las series de nivel de los embalses obtenidas se utilizan para forzar los embalses en simulaciones hidrológicas.

Palabras clave: Recursos hídricos, teledetección, radar de apertura sintética (SAR), humedad del suelo, altimetría

Context

Water is one of the most common substances on Earth [START_REF] Thompson | ENVIRONMENTAL HEALTH PRACTI-TIONER MANUAL: A RESOURCE MANUAL FOR ENVIRONMENTAL HEALTH PRAC-TITIONERS WORKING WITH ABORIGINAL AND TORRES[END_REF]. It exists in the oceans, on the surface, below the ground, and in the atmosphere. Water resources are defined as "Water available or capable of being made available for use, in sufficient quantity and adequate quality" (WMO, 2012). Therefore, water resources normally refer to freshwater that mostly exists in rivers, lakes, surfaces, and underground.

Although water needs are generally associated with access to drinking water for the population, it is also essential for many industrial and agricultural sectors. Agriculture, which uses water for crop irrigation, accounts for 3/4 of the world's current demand [START_REF] Postel | Water scarcity[END_REF][START_REF] Shiklomanov | Appraisal and assessment of world water resources[END_REF][START_REF] Cai | Global water demand and supply projections part1: A modeling approach[END_REF][START_REF] Wisser | Global irrigation water demand: Variability and uncertainties arising from agricultural and climate data sets[END_REF][START_REF] Sauer | Agriculture and resource availability in a changing world: The role of irrigation[END_REF][START_REF] Siebert | Development and validation of the global map of irrigation areas[END_REF]. This high demand for water is due to the expansion and intensification of agriculture to cope with the increase in the world's population and thus the need for food. Food and Agriculture Organization (FAO) of the United Nations estimated that the global demand for food, feed and fibre will grow by 70 percent in the first half of this century (FAO, 2009). In addition to this growing demand, climate change, which is already being felt in different parts of the world, is impacting water resources [START_REF] Vörösmarty | Global water resources: Vulnerability from climate change and population growth[END_REF].

Climate change and interaction with water resources

Climate change, as a scientific consensus, is undergoing and is primarily due to human activities [START_REF] Cook | Quantifying the consensus on anthropogenic global warming in the scientific literature[END_REF]. The Earth is globally warmer by 0.85 • C on average since the beginning of the industrial era, according to the 5th report of the Intergovernmental Panel on Climate Change (IPCC) [START_REF] Pachauri | Intergovernmental Panel on Climate Change[END_REF]. Also, the global surface temperature is projected to rise by the end of the 21st century (2081-2100) to about 0.3 • C to 1.7 • C under RCP2.6 (Representative Concentration Pathway with radiative forcing of 2.6 W/m 2 ) and 2.6 • C to 4.8 • C under RCP8.5, with reference to 1986-2005 as shown in Figure 1.1 (a) [START_REF] Pachauri | Intergovernmental Panel on Climate Change[END_REF]. It is very likely that heat waves will occur more often and last longer, and extreme precipitation events, as well as drought events, will become more intense and frequent in many regions [START_REF] Easterling | Climate extremes: Observations, modeling, and impacts[END_REF]. The high latitudes and the equatorial Pacific are likely to experience an increase in annual mean precipitation, while many mid-latitude and dry subtropical regions will likely suffer from drought as shown in (left) and RCP8.5 (right) scenarios. RCP is short for Representative Concentration Pathway, which is a greenhouse gas concentration (not emissions) trajectory adopted by the IPCC for its fifth Assessment Report in 2014 [START_REF] Pachauri | Intergovernmental Panel on Climate Change[END_REF].

The impact of climate change is still uncertain in many areas of the world. To have knowledge of the climate change impact over different regions, a Regional Climate Change Index (RCCI) is developed based on regional mean precipitation change, mean surface air temperature change 1.1. CONTEXT Figure 1.2: Regional Climate Change Index (RCCI) over 26 land regions of the world calculated from 20 coupled Atmospheric-Ocean General Circulation Models and 3 IPCC emission scenarios (A1B, A2, B1) [START_REF] Giorgi | Climate change hot-spots[END_REF]PAGES , 2018). and change in precipitation and temperature interannual variability by [START_REF] Giorgi | Climate change hot-spots[END_REF]. The RCCI identifies the most responsive regions to climate change or Hot-Spots. As shown in Figure 1.2, the Mediterranean and North Eastern European regions emerge as the primary Hot-Spots, followed by high latitude northern hemisphere regions and by Central America [START_REF] Giorgi | Climate change hot-spots[END_REF].

Global warming will lead to increased evaporation, leading to an intensification (or acceleration) of the water cycle. The intensification of the global water cycle could affect water availability and may increase the frequency and intensity of tropical storms, floods, and droughts [START_REF] Huntington | Evidence for intensification of the global water cycle: Review and synthesis[END_REF]. In humid climates, increased evaporation will likely result in more precipitations and increase the risk of flooding. Over arid areas, increased evaporation will reduce water availability, prolong dry periods and aggravate droughts. In some regions, mean annual precipitation will increase, but the number of events will be reduced, which could lead to longer dry seasons and a higher risk of floods.

Furthermore, the increase evapotranspiration (ET) will decrease the effective precipitation in many areas of the world even with the precipitation increases, affecting water resources and agriculture. In a regional perspective of view, the uncertainty of the regional climate, the diversity of the hydrology, and the human activities could all affect the regional water resources. Climate change is projected to reduce renewable surface water and groundwater resources in most dry semitropical regions [START_REF] Pachauri | Intergovernmental Panel on Climate Change[END_REF]. Extreme challenges for the sustainable management of water resources exist at all levels from the local to the global scale [START_REF] Trenberth | Challenges and Opportunities in Water Cycle Research: WCRP Contributions[END_REF]. Therefore, the water resources context and problems need to be known both in the global and regional scale.

Water resources in global scale

Freshwater is a renewable resource, yet large regions of the world faces the problem of water scarcity with its patterns varying in time and space. From the perspective of the entire water cycle, water resources are significantly affected by global climate change and human activities. Climate change or global warming increases the uncertainty associated with the future availability and variability of freshwater resources, and may even lead to permanent desertification of certain regions of the world (UNESCO, 2011). Extreme weather happens more often than before including droughts and floods, bringing more challenges in water resources management. More than two billion people are affected by water scarcity worldwide [START_REF] Oki | Global hydrological cycles and world water resources[END_REF], and a total of 2.7 billion have water scarcity for at least one month of the year [START_REF] Hoekstra | Global monthly water scarcity: Blue water footprints versus blue water availability[END_REF]. In less than 25 years, two-thirds of the world's population will be living in water-stressed countries (UNESCO, 2011).

Water scarcity is exacerbating in many parts of the world due to climate change. In the meanwhile, the demand for water resources is increasing due to population growth and economic development.

Population density and per-capita resource use have increased dramatically over the past few decades. The United Nations predictions estimate that the world population will reach 9 billion people in 2050. This exponential growth in population -a major driver of energy consumption and anthropogenic climate change -is also the key driver behind the hydrologic change and its impacts (UNESCO, 2011).

Urbanisation with a large number of people moving from the countryside to cities aggravates the impact of population growth. More than half of the world's population currently lives in cities, and this number will reach 70 percent by 2050 (UNESCO, 2011). With the increase of the population density in cities, the water demand is significantly increasing making water stress more severe for cities.

The demand for agricultural products from the food we eat, to the cotton we wear, is increasing with the rapid population growth. More than 40% of global food and agricultural commodity output [START_REF] Shiklomanov | Assessment of water resources and water availabilit in the world: Scientific and technical report[END_REF]; WMO, 1997; Vörösmarty and Sahagian, 2000) depends on irrigation, which requires a large amount of water. Up to 70% of water is currently withdrawn for agricultural use [START_REF] Hamdy | Coping with water scarcity: water saving and increasing water productivity[END_REF][START_REF] Wolfe | [END_REF]), yet both irrigated and rainfed agriculture are exposed to water risk. Intensification of agriculture is presently contributing to deforestation and desertification, which in turn reduces the holding water capacity of the soil and increases the water scarcity [START_REF] Johnson | Governing Global Desertification: Linking Environmental Degradation, Poverty and Participation[END_REF][START_REF] Keys | Global change and the intensification of agriculture in the tropics[END_REF].

As the global population continues to increase, authorities are finding ways to secure water supplies through dam construction and water transfer. The number of dams and water transfer systems has risen dramatically in the last century to deal with the uneven distribution of water in space and season. According to [START_REF] Lehner | High-resolution mapping of the world's reservoirs and dams for sustainable river-flow management[END_REF] [START_REF] Lehner | High-resolution mapping of the world's reservoirs and dams for sustainable river-flow management[END_REF], there are approximately 40 000 registered large dams in the world, with most of them being used for irrigation, followed by the use for hydropower, water supply, and flood control.

Arid and semi-arid areas are facing the most significant challenges of water stress, while many densely populated parts of the humid tropics and temperate zones are also affected by water stress. Water stress is transcending national boundaries and is increasingly leading to competition for water at local, regional, and international levels.

On the other hand, floods which occur periodically in many parts of the world put a significant threat to crops, settlements, infrastructure, and most importantly life. From 1980 to 2009, floods caused more than 500 000 deaths and affected more than 2.8 billion people worldwide. [START_REF] Doocy | The human impact of floods: a historical review of events 1980-2009 and systematic literature review[END_REF]. It is essential to monitor and manage the water resources in order to provide early-warning alerts for flood predictions.

1.1. CONTEXT 5 Water risk, which comes both from the water stress and the risk of flooding, is a global issue as shown in Figure 1.3 [START_REF] Gassert | Aqueduct global maps 2.1: Constructing decision-relevant global water risk indicators[END_REF]. The impacts of floods and droughts in many areas will have to be managed more frequently than before, and few areas of the world are not affected [START_REF] Kundzewicz | Coping with variability and change: Floods and droughts[END_REF][START_REF] Sadoff | Better water resources management: Greater resilience today, more effective adaptation tomorrow[END_REF]. As problems related to water spatiotemporal variability manifest themselves at smaller scales, a significant challenge in freshwater assessments is how to handle this on different geographical scales (Sophocleous, 2004a). Water resources in regional scale are the priority concern for local government to make policy decisions and for scientists to learn the mechanism of the water cycle.

Water resources in regional scale

At the regional scale, water resources are often considered within a river basin, within which the water budget can be closed (groundwater permitting). Humans directly change the dynamics of the water cycle through dams that are constructed for water storage or hydropower, and through water withdrawals for industrial, agricultural, or domestic purposes [START_REF] Haddeland | Global water resources affected by human interventions and climate change[END_REF]. Everything we do affects the river basin, while in reverse, water resources in the river basin affect humans by determining the food yield, the diversity of animal species, the speed of economic growth, and the number of people and livestock that can be sustainably supported. Among all the human activities, agriculture is the dominant component of human water use. Hence, many of the solutions to water-related food and environmental security problems come from the increasing of the efficiency of water use for the productivity of agriculture (Sophocleous, 2004b).

River basins in arid and semi-arid areas face a more significant threat to water shortage and a bigger challenge of water resources management when encountering the need to support a large amount of population. In a semi-arid environment, the water cycle has particular characteristics. Precipitation and many planetary boundary layer processes are affected by soil moisture variability and, precisely in these areas, soil moisture is very variable in space and time, making the areas very vulnerable to climate change with consequences that may have serious social and environmental effects (UNESCO, 2011). The water resources in semi-arid regions, with their vulnerability to climate change, need to be analysed, and the impact of climate change needs to be assessed.

The Mediterranean, as mentioned before, is one of the primary Hot-Spots of climate change. The precipitation is predicted to decrease as shown in Figure 1.1 (b), and droughts will occur more frequently. The Mediterranean region is often cited as an example of the extreme tension on water resources: 180 million inhabitants have less than 1000 m3/year per capita of renewable water resources (425 in Tunisia and 750 in Morocco), and 80 million will have less than 500 m3/year in 2025 (UNEP/MAP , 2009; EASAC , 2010). In most Mediterranean countries, overexploitation of water resources is already observed, and this situation is bound to deteriorate further. The water management is then, with its corollary of food security, the critical issue of the Mediterranean in the 21st century.

Water cycle and anthropogenic interventions

Water is not a stationary substance on Earth. It is moving all the time. To study the water resources, which form water cycle, their movements need to be known first. The water cycle is the movement of water between the oceans, the atmosphere, the surface, and the fresh water. It consists of several processes including precipitation, evaporation, transpiration, infiltration, and runoff, as shown in Figure 1.4 [START_REF] Green | Global energy and water cycles[END_REF]. Humans influence the water cycle mainly by influencing evapotranspiration (evaporation + transpiration) through irrigation for agriculture or land use change and runoff through damming for water storage, hydro powers or flood control. The other processes are subsequently influenced. Anthropogenic influences have affected the global water cycle since the beginning of 1.2. WATER RESOURCES ESTIMATION agriculture, and the impact has increased dramatically since the 1960s, as dams were built extensively [START_REF] Pachauri | Intergovernmental Panel on Climate Change[END_REF].

With the construction of dams, reservoirs are usually formed from the water accumulating above them. Proper management of dams is of great importance since it relates not only to food production but also to floods, droughts, and energy. Moreover, the severe mismanaging of the dams may cause droughts in the downstream of the river or increase the risk of catastrophic floods. Several studies [START_REF] Li | Reduced runoff due to anthropogenic intervention in the loess plateau, china[END_REF][START_REF] Haddeland | Global water resources affected by human interventions and climate change[END_REF][START_REF] Gao | Impact of climate change and anthropogenic activities on stream flow and sediment discharge in the wei river basin, china[END_REF][START_REF] Guo | Effects of the three gorges dam on the yangtze river flow and river interaction with poyang lake, china: 2003-2008[END_REF] found that the construction of dams reduces river flows and changes the stream flow seasonally, which influences the water cycle significantly. The impacts of reservoirs are most profound in Asia, Europe, and Africa, where in some months the total flux of fresh water into the ocean is 10% less compared to a naturalised situation [START_REF] Biemans | Impact of reservoirs on river discharge and irrigation water supply during the 20th century[END_REF]. Therefore, in addition to the satisfaction of the irrigation or energy purpose, the impact and role of dams in the water cycle need to be accounted for.

Equally important as dam construction, but still related to the dams, irrigation is another anthropogenic intervention of the water cycle. To fulfil the increasing global food demand, irrigation is increasing, making more contributions to evaporation as an important step in the water cycle. Irrigation creates an anthropogenic change of the soil moisture, and can ultimately influence the planetary boundary layer, local circulation, clouds and precipitation through evaporation from soil and transpiration from thriving vegetation. Several studies (Tuinenburg et [START_REF] Lo | Irrigation in california's central valley strengthens the southwestern u.s. water cycle[END_REF] identified that irrigation in the Central Valley (California, United States) increased summer precipitation by 15%, causing correspondingly the increase in Colorado River stream flow of about 30%. Boucher et al., 2004 (Boucher et al., 2004) estimated a global mean radiative forcing in the range of 0.03 to +0.1 W/m 2 due to the increase in water vapour from irrigation, and a massive surface cooling of up to 0.8 K over irrigated land areas. Climate change with rising temperature globally increases atmospheric evaporation demand, and also, evaporation caused by irrigation influences the climate. The water cycle involves the exchange of energy, which leads to temperature changes. The heat exchanges from water evaporation and condensation influence climate and the impact need to be studied and involved in the water cycle process. Also, as the primary purpose of reservoirs, irrigation has been found to decrease the river discharge to the oceans both on a continental and global scale [START_REF] Gerten | Causes of change in 20th century global river discharge[END_REF][START_REF] Biemans | Impact of reservoirs on river discharge and irrigation water supply during the 20th century[END_REF].

Land use changes also influence the water cycle by changing the evapotranspiration. Land use changes may have many situations including agriculture to buildings, agriculture type changes, etc. With the complexity, it is difficult to quantify the influences of land use change, and thus, out of the scope of this thesis.

Water resources estimation

To deal with droughts, water shortages, and scarcity in many areas of the world, effective water management is needed to provide some of society's most basic needs from agricultural crops to electricity. Monitoring the state of water resources and their change over time is integral to the development of effective strategies for sustainable water resources management. Water shortage or flood events forecast is essential for minimising the effects of disasters and planning the water allocation for irrigation and industrial usage, yet the forecast depends on the knowledge of the current status of the water resources and the historical changing patterns.

Among all the water resources in a basin scale, the surface water that exists in the reservoirs, CHAPTER 1. INTRODUCTION lakes and rivers is the most in need of monitoring and management since it is often the direct source of irrigation, hydropower and flood. The reservoirs and lakes work as a buffer to the river flow and storage to water which can be further used for irrigation. Soil moisture is also a valuable water resource for agriculture, especially for rainfed fields. It is the direct water sources for the crops and decides the amount of water needed for irrigation.

Hydrological models and land surface models

In response to water-related challenges, many hydrological models have been developed to analyse, understand, and explore solutions for sustainable water management in the basin scale, to support decision makers and operational water managers. Hydrological models are mainly used for understanding hydrological processes, estimating the water resources, and predicting system behaviour. They have evolved from empirical rainfall-runoff models [START_REF] Mulvany | On the use of self-registering rain and flood gauges in making observations of the relations of rainfall and of flood discharges in a catchment[END_REF][START_REF] Sherman | Stream flow from rainfall by the unit graph method[END_REF][START_REF] Horton | Surface Runoff Phenomena. Part 1: Analysis of the Hydrograph[END_REF] to conceptual models [START_REF] Linsley | Computation of a synthetic streamflow record on a digital computer[END_REF] [START_REF] Todini | Hydrological catchment modelling: past, present and future[END_REF]. The recent physical models use parameters to characterise the unique aspects of the system being studied [START_REF] Nepal | Assessment of spatial transferability of process-based hydrological model parameters in two neighbouring catchments in the himalayan region[END_REF]. However, the parameters may differ in different catchments due to physical and hydrological differences [START_REF] Blöschl | On hydrological predictability[END_REF].

In addition to hydrological models, Land Surface Models (LSMs) are equally important to learn the regional land surface status with more complexity since it is physical with the simulation of the energy balance. The LSMs are the components of Global Climate Models (GCMs) that simulate land surface processes, such as the absorption and partitioning of radiation, moisture, and carbon [START_REF] Abramowitz | Evaluating the performance of land surface models[END_REF]. The LSMs model the surface energy balance, surface water balance and carbon balance. They are typically provided with meteorological conditions as inputs (from a boundary layer atmospheric model or reanalysis data) and produce outputs that include the surface runoff, latent and sensible heat fluxes, CO2 fluxes, deep soil drainage, soil moisture, etc. LSMs are essential to assess the impact of climate change.

A model is a simplified representation of real-world system [START_REF] Sharma | Hydrological Modelling in Arid and Semi-Arid Areas[END_REF], with the best one giving results close to reality with the use of least parameters and model complexity [START_REF] Devia | A review on hydrological models[END_REF]. However, most of the current models do not incorporate (or very little) human intervention in the water cycle, bringing more uncertainties into the assessment, analysis and forecasting. The uncertainties are even more significant in semi-arid areas since models are not accurately developed, and the impact of irrigation and reservoir operation may be amplified [START_REF] Zhou | The contribution of reservoirs to global land surface water storage variations[END_REF].

Incorporating anthropogenic intervention in the models helps to understand the hydrological processes and the climate change better, yet long-term and consistent records that focus on water management and human interventions on the water cycle are not widely available [START_REF] Zhou | The contribution of reservoirs to global land surface water storage variations[END_REF]. The modelling requires better accuracy, both in time and spatial scale, for which researchers are still working on.

To solve water-related challenges and to assess the impact of climate change, better characterised hydrological models and LSMs are needed. For that, the anthropogenic intervention needs to be introduced in the water cycle and models. Involving anthropogenic intervention in the water cycle needs to be done by introducing the irrigation status and the water availability in rivers and reservoirs considering the effects of dams.

LSMs can benefit from reliable soil moisture data for model validation and data assimilation. Soil moisture, which is directly linked with evaporation and decides the irrigation strategy,

WATER RESOURCES ESTIMATION

needs be to mapped with better efficiency to obtain a higher spatial and temporal resolution. Also, LSMs need a good mapping of irrigation, in order to have realistic covers and explicitly simulate irrigation, which is needed for cropland irrigation management (Ambika et al., 2016). The amount of irrigation surfaces, which is of great importance to water resources management, is still unclear. Therefore, to optimise the models, surfaces of irrigated and non-irrigated fields need to be mapped.

Furthermore, now that LSMs have started to include dams, they need dam level/volume data, which is often difficult to obtain. Temporal continuous water availability in rivers and reservoirs needs to be known to assess the effects of dams in hydrological processes. Precise water levels over reservoirs should be monitored to evaluate the available water amount for irrigations. The buffering effects of dams should be included and assessed in the models for optimisation.

Measurements and remote sensing techniques

The irrigation status is traditionally measured by field campaign, which is time-consuming and more importantly impossible to be done for all areas. Soil moisture is traditionally measured with in situ measurements either by taking soil samples or by the soil moisture probe. However, in situ networks represent single point locations and usually cover relatively short periods of observation (Pratola et al., 2015).

Considering the aspect of water levels, the traditional way of monitoring is using the gauging station or stream gauge, which is a location used by hydrologists to monitor the terrestrial bodies of water. Unfortunately, in situ gauging stations are not always available in many parts of the world, or they are otherwise not publicly available and are maintained by local authorities. Furthermore, over the last few decades, in situ observations of hydrologic variables have generally been in decline [START_REF] Trenberth | WCRP grand challenge: Changes in water availability[END_REF].

The access to both data sources is critical, and moreover, continuous temporal data is needed for the monitoring, assessing, modelling and predicting. Yet data are sparse in developing countries both in spatial and in temporal dimensions. Even where such data have been collected, they are rarely shared across ministries or institutions, while thousands of freshwater bodies are the shared resource of several countries.

At the same time, new observation methods that use satellite sensors have thrived over the past few decades. Thanks to the remote sensing technology, we now possess a powerful tool for monitoring the water resources both within inland water systems where no in situ data are available and within transboundary river basins. Wherever satellite data are available, the surface soil moisture, the status of irrigation, and the surface water levels can be assessed and monitored at both regional and global scales.

Challenges

The remote sensing technology used for hydrology has been developed for a few decades to overcome the limitation of lacking data on the time evolution and spatial structure, yet challenges remain to achieve the spatial and temporal resolution needed for hydrological studies.

Satellite systems, both passive and active, have already demonstrated their capability to provide reliable soil moisture measurements [START_REF] Lacava | A multi-sensor (smos, amsr-e and ascat) satellite-based soil moisture products intercomparison[END_REF] the Aqua satellite [START_REF] Njoku | Soil moisture retrieval from amsr-e[END_REF] and the Advanced SCATterometer (ASCAT) onboard the MetOp (Meteorological Operational) satellite [START_REF] Wagner | The ascat soil moisture product: A review of its specifications, validation results, and emerging applications[END_REF]. However, the soil moisture retrieved from these instruments has a low spatial resolution. SMOS and SMAP data have a spatial resolution of about 40 km (El Hajj et al., 2018), while AMSR-E and ASCAT data have a spatial resolution of about 60 km and 50 km, respectively.

Currently, ASCAT, which is a real-aperture radar instrument, provides the only active global soil moisture dataset. It was designed to observe wind speed and direction over the oceans, but has been shown to be useful to measure large-scale soil moisture [START_REF] Wagner | The ascat soil moisture product: A review of its specifications, validation results, and emerging applications[END_REF]. Although backscatter data from active sensors have the potential to monitor soil moisture, there is currently no operational soil moisture product from SAR active microwave (Escorihuela and Quintana-Seguí, 2016). SAR imagery has been shown to be advantageous for the estimation of soil surface characteristics, in particular surface roughness and soil moisture (Aubert et . However, almost five years after Sentinel-1 satellite (carries a C-band SAR instrument with a spatial resolution about 10 m) was launched, there is still no operational soil moisture product available. Soil moisture retrieval under dense vegetation cover (irrigated crop fields) is still needs to be studied further.

The similar problem exists for irrigation mapping. Studies using remote sensing to map irrigated fields remain relatively rare (Ozdogan et al., 2010). There are only a few studies which identify irrigated areas globally (Salmon et For the studying of water levels, satellite altimeters were designed to monitor homogeneous surfaces such as oceans or ice sheets, resulting in poor performance over small inland water bodies due to the contribution from land contamination in the returned waveforms. The advent of SAR altimetry with its improved along-track spatial resolution has enabled the measurement of inland water levels with better accuracy and an increased spatial resolution. However, many previous studies focused on relatively large water bodies, and most studies used low-resolution mode (LRM) altimeters onboard Envisat and ERS. The performance of SAR altimeter over middle-sized and small-sized water bodies still needs to be explored, especially over water bodies with challenging environments.

About the thesis 1.3.1 Objectives

Within the context of water scarcity and anthropological influence, the objective of this PhD study is trying to better estimate the water resources by involving the human interventions in the hydrological process by using remote sensing data. Considering the two main human interventions, irrigation and dam construction, this PhD study tries to contribute to both parts to provide parameters with higher accuracies for LSMs or hydrological models.

First, the irrigation is studied by mapping irrigated areas to better support water resources management, droughts and floods forecasts and agricultural development. In the meantime, the soil moisture will be monitored to evaluate the status of the field along with the health of the crops. Soil moisture in itself is an essential part of the regional water resources, and needs to be quantified more precisely. Therefore, the objective of this part is to provide soil moisture product and irrigation maps with satellite data at better spatial accuracy. SAR data provides the possibility of monitoring land surface under any weather conditions and with vegetation cover. Therefore, this PhD study used the recent SAR product from Sentinel-1 mission to develop methodologies of soil moisture retrieval and irrigation mapping at a field scale.

The second aspect considered is the influence of dam constructions for the river basin. Dams constructed aim not just for hydropower, but more importantly for irrigation purpose. To study the influence of dams, a river flow simulation model is considered with which the buffering effects of dams should be included. The water volumes within reservoirs, which is directly linked to the water levels, will be estimated both for analysing the dam effects in the river model and for irrigation activities, whose water withdraw depends on the storage of water in reservoirs. To retrieve the water levels of the reservoirs, the recent SAR altimetry product from the Sentinel-3 mission is considered, since it provides possibilities of water level retrieval over middle-sized and small-sized water bodies.

To conclude, this PhD thesis is to improve the regional water resource estimation by improving both soil moisture and water level products. It will involve the human interventions in each part for a better assessment, including the irrigation study and involving of the dam levels in river flow simulation model.

Structure of the thesis

As a first step, soil moisture retrieval has been done by two change detection methodologies, without field calibrations, and validated over an area characterised by a dense vegetation cover (irrigated crop fields). The soil moisture is retrieved from the synergetic interpretation of Sentinel-1 and Sentinel-2 data. The mean soil moisture is computed at a resolution of 100 m, which is compatible with agricultural applications. Before the 100 m study, soil moisture was also retrieved at 1 km resolution as a test of the methodologies. Despite Sentinel-1 spatial resolution being around 10 m, the soil moisture is estimated a lower resolution (100 m) in order to decrease uncertainties caused by different types of heterogeneities in agricultural fields such as local changes in roughness, heterogeneities in vegetation cover, etc.

For the second step, the irrigation mapping method has been developed using SAR data time series at a field scale. Four metrics, including the mean temporal value, the signal variance, the correlation length, and the fractal dimension were analysed for the mapping. The classification was done with the Support Vector Machine (SVM).

In the third step, water level retrieval has been improved by three specialised algorithms or retrackers to retrieve water levels from SAR altimeter data over inland water bodies dedicated to minimising land contamination from the waveforms. Three different retrackers were modelled, namely, the threshold retracker, the Offset Centre of Gravity (OCOG) retracker, and the 2step physical-based retracker; the parameters of these retrackers were adjusted for inland water bodies. Additionally, the digital elevation model (DEM) information was included to limit land contamination within the received waveforms and to improve the accuracy for the monitoring of middle-sized and small-sized water bodies.

In the last step, the satellite-derived water levels were converted to water volumes and involved in a river flow model to see the buffering effect of dams in river flow simulations. The observed levels were used to tell the model to try to be close to those levels, which is a technique to include observed information in the simulation.

The following manuscript is structured in 6 parts:

• Chapter 2, offering a detailed description of the remote sensing technology related to soil moisture and water level monitoring.

• Chapter 3, presenting the study area and database

Introduction

According to Caloz and Collet (Caloz and Collet, 2001), remote sensing is defined as the technique of remote observation by the measurement and the analysis of the electromagnetic radiation emitted or reflected by the studied object to interpret the information about its nature, its properties and status. It is based on the use of sensors to record electromagnetic radiation and converting it into a digital output signal.

Satellite-based remote sensing has been rapidly developed after the first launch of an artificial satellite named Sputnik I in 1957 [START_REF] Tatem | Fifty years of earth observation satellites: Views from above have lead to countless advances on the ground in both scientific knowledge and daily life[END_REF]. It is a valuable tool in many study areas including weather forecast, climate, hydrology, oceanography, and land use.

There are two types of remote sensing, namely passive and active (Gupta et al., 2013). Passive remote sensing needs an external energy source for objects to reflect or radiate, which in most cases is the sun. Active remote sensing system has its own energy source (Campbell and Wynne, 2011). The sensor onboard the satellite or aircraft emits a signal towards the target surface, and receives its reflection by the target object for further processing.

Electromagnetic waves

Electromagnetic (EM) waves are energy transported through space as a result of periodic disturbances of an electric field and a magnetic field [START_REF] Rammer | Quantum field-theoretical methods in transport theory of metals[END_REF]. Most sensors record information about the Earth's surface by measuring the transmission of energy from the surface in different portions of the EM spectrum (Figure 2.1) [START_REF] Aggarwal | Principles of remote sensing[END_REF]. The EM 16 CHAPTER 2. REMOTE SENSING PRINCIPLES spectrum can be divided into seven different regions: gamma rays (highest frequency and shortest wavelength), X-rays, ultraviolet, visible light, infrared, microwaves and radio waves (lowest frequency and longest wavelength) [START_REF] Wait | Electromagnetic Wave Theory[END_REF].

Remote sensing technology uses different bands of electromagnetic spectrum for different purposes. The most common bands include visible spectrum, infrared and microwave with microwave domain most sensitive to soil moisture and can be used to measure distance, which is the principle of altimetry. Data collected over a large number of wavelength bands is called multispectral or hyperspectral data.

Characteristics of electromagnetic waves

EM waves travel through space at the same speed, c = 2.99792458 * 10 8 m/s, commonly known as the speed of light. Besides this, EM waves have certain characteristics, including amplitude, wavelength/frequency, phase and polarization. Amplitude Amplitude is the distance from the maximum vertical displacement of the wave to the middle of the wave. It measures the magnitude or intensity of the oscillation of a wave, with larger values corresponding to higher energy and lower values corresponding to lower energy.

Wavelength Wavelength is the distance of one full cycle of the wave oscillation. Wavelength and frequency are related to the speed of light by the equation:

c = f * λ (2.1)
where c is the speed of light, f is the frequency, and λ is the wavelength. Longer wavelength corresponds to lower frequency, which is proportional to energy, whilst shorter wavelength corresponds to the higher frequency and higher energy. Phase Phase is the position of a point in time on a waveform cycle. Phase information can be used to correct Doppler effects, calculate time delay, and is an important parameter in interferometry for deformation monitoring.

Polarization Polarization is a property of the transverse waves that specifies the geometrical orientation of the oscillations [START_REF] Shipman | An Introduction to Physical Science[END_REF]. The EM wave consists of two fluctuating fields which are the electric field (E) and magnetic field (B). The two vectors are orthogonal to one another, and both are perpendicular to the direction of travel (Figure 2.2). The polarization of the electromagnetic waves refers to the orientation of E [START_REF] Smith | Electromagnetic wave propagation in media with indefinite permittivity and permeability tensors[END_REF]. If both E and B remain in their respective planes, the radiation is called plane or linearly polarized [START_REF] Sandwell | Satellite remote sensing sio 135/sio 236[END_REF]. There are two types of polarization including vertically polarization with E parallel to the plane of incidence, and horizontal polarization with E perpendicular to the plane of incidence. The plane of incidence is the plane defined by the vertical and the direction of propagation. Linear polarization is normally used by active remote sensing systems.

Microwave remote sensing

In particular, microwave remote sensing uses electromagnetic radiation with a wavelength between 1 cm and 1 m (300 MHz to 30 GHz commonly referred to as microwaves) as a measurement tool.

While electromagnetic waves in the visible and infrared portions of the spectrum primarily characterised by wavelength, microwave portions of the spectrum are often referenced according to both wavelength and frequency (Natural Resources Canada, 2015). The microwave region of the spectrum is quite large, relative to the visible and infrared, and there are several wavelength ranges or bands commonly used which given code letters during World War II, and remained to today as listed in Table 2.1 (du Preez and Sinha, 2016).

Microwave remote sensing has many advantages compared to optical remote sensing. It has a day-and-night imaging capability, and ability to penetrate could cover, haze, fog, and to some extent, rain as shown in Figure 2.3. For lower frequencies at L-and P-band, the transmission rate is almost 100%. Moreover, microwave remote sensing is more sensitive to the dielectric constant of the medium and structural characteristics of the surface, such as roughness or canopy structure compared to optical sensors (

Radar

Similar to remote sensing systems, microwave remote sensing systems are also classified into two groups: passive and active. Passive systems collect the radiation in microwave bands that is naturally emitted by the target surface, and they are generally characterised by relatively low spatial resolutions.

On the contrary, active microwave systems emit their own waves in the microwave bandwidth from the transmitter towards the observed surface. The power of the backscattered signal can be used to discriminate different targets within the scene, while the time between the sent and the received signal is used to measure the distance of the target. A system that operates in this way is called radar (RAdio Detection And Ranging), and also can obtain a microwave image of the observed scene.

The fundamental relation between the characteristics of the radar, the target, and the received signal is called the radar equation (Ulaby et al., 1986b). The geometry of scattering from an isolated radar target (scatterer) is shown in Figure 2.4, along with the parameters that are involved in the radar equation. When a power P t is transmitted by an antenna with gain G t , the power per unit solid angle in the direction of the scatterer is P t * G t . The power density at the scatterer S s is:

S s = (P t G t )( 1 4πR 2 t )A rs (2.2)
where A rs is the effective receiving area of the scatterer, and the spreading loss

( 1 4πR 2 t
) is the reduction in power density associated with spreading of the power over a sphere of radius R surrounding the antenna.

Some of the power received by the scatterer is absorbed in losses in the scatterer, unless it is a perfect conductor or a perfect isolator; the rest is re-radiated in various directions (Ulaby et al., 1982). The fraction absorbed is F a , so the fraction reradiated is 1 -F a .

With the spreading factor for the re-radiation

( 1 4πR 2 r
), and the effective aperture of the receiving antenna A r , the gain at the receiver can be expressed as:

P r = (P t G t )( 1 4πR 2 t )A rs (1 -f a )G ts ( 1 4πR 2 r )A r = ( P t G t A r (4π) 2 R 2 t R 2 r )[A ts (1 -f a )G ts ] (2.3)
The factors associated with the scatterer are combined in the square brackets, which are difficult to measure individually, and their relative contributions are irrelevant for the size of the received radar signal. Hence they are normally combined into one factor, the radar scattering cross section:

σ = A ts (1 -f a )G ts (2.4)
Therefore, the final form of the radar equation is obtained by rewriting equation 2.3:

P r = ( P t G t A r (4π) 2 R 2 t R 2 r )σ (2.5)
Most radar sensors are designed to transmit microwave radiation either horizontally polarized (H) or vertically polarized (V). Similarly, the antenna receives either the horizontally or vertically polarized backscattered energy, or both. Thus, there can be four combinations considering the two polarization states for transmission and reception: CHAPTER 2. REMOTE SENSING PRINCIPLES HH -for horizontal transmit and horizontal receive, VV -for vertical transmit and vertical receive, HV -for horizontal transmit and vertical receive, VH -for vertical transmit and horizontal receive. Radar can measure amplitude, which is the strength of the reflected signal, and phase, which is the position of a point in time on a waveform cycle [START_REF] Podest | Basics of synthetic aperture radar[END_REF]. Amplitude depends on the target properties such as structure and dielectric properties, while the phase is a function of the distance between the sensor and the target as well as target properties. Therefore, radar sensors are generally divided into two distinct categories: imaging and non-imaging.

Synthetic Aperture Radar

Similar to conventional radar, Synthetic Aperture Radar (SAR) sequentially transmits electromagnetic waves towards the Earth surface. The microwave energy scattered back to the antenna is measured. The consecutive time of transmission and reception converts into different locations due to the platform movement [START_REF] Moreira | A tutorial on synthetic aperture radar[END_REF]). An appropriate coherent combination of the received signals makes it possible to construct a virtual aperture much longer than the physical antenna length [START_REF] Skolnik | Radar Handbook, Third Edition, Electronics electrical engineering[END_REF]. The slant-range resolution δ r is inversely proportional to the system bandwidth with the twoway path from transmission to reception, being able to detect the minimum distance between two targets:

δ r = c 0 2B r (2.6)
where c 0 is the speed of light.

The azimuth resolution δ a is provided by the synthetic aperture, which is the path length during which the radar receives echo signals from a point target [START_REF] Chan | An introduction to synthetic aperture radar (sar)[END_REF]. The beamwidth of an antenna with length d a can be approximated with wavelength λ:

θ a = λ d a (2.7)
The synthetic aperture length is given by:

L sa = θ a * r 0 = λr 0 d a (2.8)
where r 0 is the distance between the antenna and the target object. Therefore, the synthetic beam width can be modelled as:

θ sa = λ 2L sa = d a 2r 0 (2.9)
The azimuth resolution is half of the real radar aperture calculated as follows:

δ a = r 0 * θ sa = d a 2 (2.10)

Imaging systems

The imaging radars map variations in microwave backscatter at fine spatial scales (1 to 50 m) and can be used to measure variations in surface roughness and surface moisture.

The received echo signal data is stored in a two-dimensional data matrix in forms of complex samples, with each sample given by its real and imaginary part, thus representing an amplitude and phase value [START_REF] Moreira | A tutorial on synthetic aperture radar[END_REF]. The first dimension corresponds to the range direction or fast time, while the second dimension corresponds to the azimuth direction or slow time.

In contrast to optical sensors, visualising raw SAR data does not give any useful information on the target [START_REF] Lillesand | Remote Sensing and Image Interpretation[END_REF]. The interpretable image is only available after signal processing, which steps are shown in Figure 2.6. The whole process can be comprehended as two separate matched filter operations along the range and azimuth dimensions. In the first step, the transmitted chirp signals are compressed into a short pulse, and a convolution is performed in the time domain instead of the frequency domain due to the much lower computational load [START_REF] Moreira | A tutorial on synthetic aperture radar[END_REF]. Thus, a range of the compressed image is obtained with each range line multiplied in the frequency domain by the complex conjugate of the spectrum of the transmitted chirp. Until now, only information about the relative distance between the radar and possible targets is known. In the second step, the compression in azimuth direction is solved by convolving with the reference function, which is the complex conjugation of the response expected from a target position on the ground.

Scattering mechanisms

SAR images represent an estimate of the radar backscatter for the target area on the ground. The strength of the reflected signal is called backscattering coefficient (σ 0 ), and is expressed in decibels (dB). Higher values of the backscatter coefficient mean that an expansive part of the radar energy was reflected back to the radar, shown as bright features, while low values imply that very little energy was reflected, shown as dark features. The values of the backscatter coefficient over a target area with a particular wavelength varies a lot for different conditions and different characteristics of the targets, including surface roughness and dielectric properties.

Surface roughness: For a fixed looking angle of the satellite, the terrain geometry changes the incidence angle, therefore, changing the backscattered signal received by the sensors. For flat terrain, specular reflection results in very low energy of the returned signal. For a rough surface, on the other hand, incidence energy is scattered in all directions and return a significant portion of the incident energy back to the antenna (ESA, 2007), as shown in Figure 2.7. Vegetation is typically moderately rough compared with the scale of most radar wavelengths and returns moderate scale of energy, which results in grey or light grey colour in the radar image. A particularly bright response comes from a corner reflector, which produces a double bounce. City streets and buildings often perform as corner reflectors and, subsequently, appear to be very bright in radar images.

Dielectric properties:

Radar backscatter also relies upon the dielectric properties of the target: for metal and water the dielectric constant is high (80), while for most other materials it is relatively low: in dry conditions, the dielectric constant ranges from 3 to 8. The wetness of soils can deliver a prominent increase in radar signal reflectivity. Based on this phenomenon, SAR systems are also used to retrieve the soil moisture content (primarily of bare soils), with wetter objects appearing bright, and drier targets appearing dark (Massonet and K.L.Feigl , 1998). The exception to this is a smooth water body, which acts as a flat surface. The incoming pulses of the sensor are reflected away by the flat water bodies, which appear to be dark in radar images.

Topographic effects

Topographic effects are image distortions due to the side-looking geometry of the radar sensor, particularly in areas with pronounced topography [START_REF] Bayer | Terrain influences in sar backscatter and attempts to their correction[END_REF] as shown in Figure 2.8. Foreshortening is a dominant effect in SAR images of mountainous areas, with the time delay between the radar echoes received from two different points determining their distance in the image (ESA, 2018). In case of a very steep slope or a very big incidence angle, layovers will occur. The targets in the valley will have a larger slant range than related mountain tops, showing a "reversed" effects in the slant range image.

These distortions are of a geometric and radiometric nature and need be reduced by terrain correction with the use of a digital elevation model, e.g. SRTM 1 ArcSec (Jarvis et al., 2008).

Speckle effect

A particular effect existing in SAR images is the so-called speckle effect, which is caused by the presence of many elemental scatterers with a random distribution within a resolution cell [START_REF] Moreira | A tutorial on synthetic aperture radar[END_REF]. The coherent sum of their amplitudes and phases results in strong fluctuations of the backscatters from cell to cell of the image.

A certain process needs to be applied to minimise the speckle effect in order to get a cleaner image for interpretation. One methodology is using multi-look processing. Multi-look method divides the radar beam into several, narrower sub-beams with each one representing one "look" at the scene. By summing and averaging the different "looks" together, the amount of speckle can be reduced in the final output image [START_REF] Podest | Basics of synthetic aperture radar[END_REF]. The other speckle reduction methodology commonly used is spatial filtering using a moving window over the whole image. A mathematical calculation of all pixels within the window is performed to replace the value of the central pixel. Speckle filtering reduces the visual appearance of speckle and applies a smoothing effect of the image.

Soil moisture principle

Compared with other components of the hydrologic cycle, the volume of soil moisture is small. However, it is essential for hydrology, meteorology and agriculture.

In soil moisture study, surface soil moisture products have been produced both based on the active and passive microwave, i.e. radars and radiometers separately. Two types of radiometers can be identified: L-band (at a frequency of 1.4 GHz) and C-band (at a frequency of 7 GHz). Notable L-band radiometers include SMOS [START_REF] Kerr | The SMOS mission : new tool for monitoring key elements of the global water cycle[END_REF] and SMAP (Entekhabi et al., 2010). Because of the vast size required for the antenna in order to have sufficient sensitivity to capture signals, combined with a fixed weight for the payload, the resolution of L-band data is limited to 30-60 km [START_REF] Kerr | The SMOS mission : new tool for monitoring key elements of the global water cycle[END_REF][START_REF] Njoku | Passive microwave remote sensing of soil moisture[END_REF][START_REF] Schmugge | Observations of coherent emissions from soils[END_REF]. However, this resolution is too coarse for the majority of hydrological and agricultural applications. Thanks to the radar sensors, the resolution of the image could be improved to several tens of meters. The high resolution achieved by active sensors (ERS, ALOS, Sentinel-1) is due to the intensity of the signal backscattered is much higher than the energy emitted by land surfaces measured by radiometers.

The radar signals backscattered by the surface can be modelled as the sum of the radar signals scattered by the bare soil (soil moisture and soil roughness) and attenuated by vegetation effects, and also the signals scattered by the vegetation cover. These two contributions can be expressed as:

σ 0 cover = σ 0 veg + λ 2 (θ)σ 0 soil (2.11)
where λ 2 (θ) = exp[-2τ /cos(θ)] is the two-way vegetation canopy transmissivity, θ is the incidence angle, and τ is the optical thickness parameter that depends on the type of geometrical structure and vegetation water content of the canopy (Attema and Ulaby, 1978;Gupta et al., 2013).

The irrigation status can be assessed by analysing the temporal changes and characteristics of the backscatter coefficient from the SAR signal since it largely relates to the soil moisture of the field. Details of the methodology are presented in Chapter 5.

Ranging system and principle of altimetry

Non-imaging radars, including altimeters, measure the elevation of the earth's surface, and scatterometers are used to estimate ocean wind speed. In this section, the altimeter, which is the ranging system, is discussed.

Space-borne radar altimeters are essential tools for monitoring the oceans, which has been performed for over 20 years (Benveniste, 2011). Satellite altimeters have also proven to be valuable tools for monitoring the water levels within inland water systems, including lakes and rivers (Alsdorf Space-borne radar altimeters transmit a short microwave pulse in the nadir direction, and the signal reflected by the surface is received by the instruments. As shown in Figure 2.9, the altimeter measures the time delay taken for a radar pulse to travel to the surface and back again. The elapsed time corresponds to the range between the satellite and the Earth's surface (Fernandes et al., 2014) by the following equation:

R = c * T /2 (2.12)
where c is the speed of light and T is the time delay of the signal. The shape of the returned echo received by the altimeter is what we call the waveform. The water level from an altimeter is derived from the measured range (R), which is subtracted from the altitude of the satellite (H alt ), after which geophysical corrections are applied as shown in the equation [START_REF] Calman | Monitoring continental surface waters by satellite altimetry[END_REF]AVISO+, 2017).

H water level = H alt -R -(C dry tropo + C wet tropo + C iono + C solid earth tide + C pole tide + C ocean tide ) -C geoid (2.13)
Corrections applied including the wet troposphere, dry troposphere, ionosphere, solid earth tide, geocentric pole tide and ocean loading tide corrections, and the geoid correction.

Trackers and Retrackers

To determine the altimeter range accurately, which is related to the water level measurement accuracy, trackers and retrackers are needed to locate the waveform window and to determine the accurate tracking point on the waveform.

Waveforms are recorded by the tracking system placed onboard the satellite. The purpose of the on-board tracker is to keep the reflected signal from the Earth's surface to be within the altimeter analysis window [START_REF] Vignudelli | Coastal Altimetry[END_REF]. There are two tracking modes available, namely the closed-loop tracking and open-loop tracking. Closed-loop tracking is the traditional mode to be used where the altimeter range window is autonomously positioned based on onboard Near Real Time (NRT) analysis of previous waveforms. Alternatively, the open-loop tracking mode is available where the altimeter range window is positioned using a priori knowledge of the surface height stored onboard the instrument in a one-dimensional along track Digital Elevation Model (DEM).

In order to obtain the highest possible accuracy on range measurements, the precise tracking point located on the leading edge (Figure 2.10) needs to be found out, which is called "waveform retracking" (Deng and Featherstone, 2006). The ground-based retracking aims to fit a model or functional form to the measured waveforms, and retrieve geophysical parameters such as the range and echo power. Functional forms can be purely empirical, or physical as in the case of the Brown ocean retracker. Details about both empirical and physical-based retrackers are presented in chapter 6.

Pulse-limited altimetry

Conventional radar altimeters are pulse-limited altimetry. The effective footprint of a pulselimited altimeter is related with the pulse duration and the width of the waveform window [START_REF] Chelton | Pulse compression and sea level tracking in satellite altimetry[END_REF][START_REF] Chelton | Chapter 1 satellite altimetry[END_REF]. The altimeter sends a signal in the nadir direction after the pulse hits the water surface, and the area contributing to the reflected power is expanding as a circle, which, after reaching its maximum size, keeps on spreading into an expanding ring with increasing diameter but a constant surface area (Figure 2.11) [START_REF] Vignudelli | Coastal Altimetry[END_REF]. The diameter of the pulse-limited footprint is generally in the range of 2 and 7 km, depending on the significant wave height, which is the average of the highest one-third of waves [START_REF] Alcorn | Chapter 17 -wave energy[END_REF]. On the contrary, when the returned pulse is determined by the width of the beam, is called beam-limited altimetry, which is the principle of laser altimeters (e.g. GLAS on ICESAT).

SAR/Delay-Doppler altimetry

Until CryoSat-2 was launched, the satellite altimeters were pulse limited in low-resolution mode (LRM). CryoSat-2 is the first altimeter with synthetic aperture radar (SAR) mode available, followed by Sentinel-3. Unlike classical pulse-limited altimeters, SAR altimeters exploit coherent processing of groups of transmitted pulses to make the most efficient use of the power reflected from the surface (Radar Altimetry Tutorial and Toolbox , 2017). Therefore, the alongtrack resolution improves significantly. SAR altimetry was first described as Delay-Doppler altimetry in 1998 (Raney, 1998) since it uses Doppler effects caused by the satellite movement in the along-track direction. Its key innovation is echo delay compensation, analogous to range cell migration correction in a burstmode SAR (Raney, 1998).

The antenna altimeter emits spherical microwave radiation towards the surface in the nadir direction. These pulses are frequency linearly modulated signals and are emitted at regular intervals defined by the Pulse Repetition Frequency (PRF) (Radar Altimetry Tutorial and Toolbox , 2017). 64 coherent pulses are emitted in a burst, and each echo of the pulse contains the response of the water surface over the original footprint as in conventional altimetry. A Fast Fourier Transform (FFT) is performed in the along-track dimension (corresponding to time) over the echoes, generating corresponding Doppler beams [START_REF] Boy | Cryosat-2 sar-mode over oceans: Processing methods, global assessment, and benefits[END_REF]. The echo power is sampled and recorded in the waveform tracking window (normally with 64 or 128 gates or waveforms samples), forming the echo waveform. The footprint of each Doppler beam on the ground is referred to as a Doppler cell (Figure 2.9).

The Delay-Doppler technique not only improves the spatial resolution in the along-track direction (down to approximately 300 m for CryoSat-2 and Sentinel-3), but additionally reduces speckle noise for a given spatial resolution cell due to a higher number of independent samples (looks) to be averaged [START_REF] Boy | Cryosat-2 sar-mode over oceans: Processing methods, global assessment, and benefits[END_REF].

Conclusion

This chapter presented the remote sensing principles and the principles related to soil moisture, irrigation and altimetry, which are the main objectives of this thesis. With the SAR technique, the resolution of the imaging system and the ranging system improves dramatically, offers more opportunities to acquire soil moisture map, irrigation map and water level products with a much better resolution. Additionally, microwave remote sensing offers the tool which is functional under various weather conditions. Therefore, the temporal resolution can be improved mostly over cloudy areas. The combination use of multi-sensor microwave remote sensing is able to enhance the potential of water resources estimation more accurately and efficiently. The Ebro river basin located in Iberian peninsula, Europe, which is one of the primary Hot-Spots of climate change, is chosen as the study area for this PhD research. Ebro River basin is located in the Iberian peninsula as shown in Figure 3.1. It is one of the most important rivers on the Iberian Peninsula that flows into the Mediterranean, with a length of approximately 928 km and a drainage basin with an area of approximately 85 550 km2 [START_REF] Barceló | The Ebro River Basin[END_REF]. With the Pyrenees located in the north-east and the Cantabrian Range in the north-west, the Atlantic wind can hardly get through the basin. Therefore, the topography of the Ebro River basin determines the semi-arid Mediterranean, with a continental influence. In the central part of the basin, aridity is the main climate characteristic with rainfall decreased dramatically [START_REF] Barceló | The Ebro River Basin[END_REF].

The river flow is irregular throughout the year, with low levels at the end of summer and high levels during the spring due to melting run-off in the Pyrenees, leading to a danger of flooding. The Ebro River is of great importance for agriculture in the summer, during which drought often 29 CHAPTER 3. STUDY AREA AND DATASETS occurs due to the continental Mediterranean climate. Nevertheless, the mean annual flow has decreased by approximately 29% during the 20th century due to many causes: the construction of dams, the increasing demands for irrigation, and the evaporation (which is higher than the precipitation due to low rainfall amounts, high sunshine intensities and strong, dry winds) from reservoirs within the river basin [START_REF] Barceló | The Ebro River Basin[END_REF].

The population in the basin is about 2.8 million people, with a density of 33 inhabitants per km 2 . However, the population is heterogeneously distributed with nearly half of the population concentrated in the cities located in the centre of the Ebro valley [START_REF] Barceló | The Ebro River Basin[END_REF]. Land use in the Ebro River basin has been traditionally based on agricultural crops, such as vineyards, orchards, and maize. Up to 783 948 hectares are dedicated to agriculture where, mainly in the mid and lower Ebro sections, they are irrigated. Nowadays, with the development of the economy, this industry is also gaining importance in most cities. Hydroelectric energy produces about 8 297 m 3 /s in 340 hydroelectric plants at the Ebro River basin. Water of the Ebro River is also used for cooling nuclear and thermic plants ( [START_REF] Barceló | The Ebro River Basin[END_REF].

The frequency and intensity of extreme events, which will be increasing with climate change, affect the regional water availability, and also influences human activities and government decisions from agricultural irrigation and energy production to flood control. According to the Intergovernmental Panel on Climate Change [START_REF] Pachauri | Intergovernmental Panel on Climate Change[END_REF], almost the entire of Europe will be negatively affected by climate change. A higher risk of droughts and floods are faced by the whole of Europe, with southern Europe more likely to suffer from droughts with the decrease of summer flows. Moreover, water withdrawals are also expected to increase in south-3.1. STUDY AREA ern Europe, amplifying the risks associated with climate change. The Mediterranean regions including the Ebro River basin are more exposed to drought risk.

In the Ebro river basin, due to the semi-arid climate, irrigation has been intensified there since the mid-20th century, though it is still limited to the main floodplains in the middle (Rodriguez , 2018). Many reservoirs (138, total water storage capacity of 6 837 hm3) and canals were constructed during this century for agricultural irrigation, electricity production and domestic use [START_REF] Comín | Management of the ebro river basin: Past, present and future[END_REF].

The models applied in the Ebro river basin contain many uncertainties. Moreover, current models do not incorporate (or very little) human intervention in the water cycle, which makes the uncertainty even larger. Within this context, the primary goal of the water resources estimation over the Ebro river basin is to estimate the water levels and soil moisture, and to involve the anthropogenic influences including the irrigation study and the dams.

Urgell site

For the soil moisture estimation and irrigation study, the Urgell site within the Ebro river basin is selected. The study area covers a 60 km by 60 km square located in Urgell, Catalonia, as shown in Figure 3.2. The 60 km by 60 km area is used for the soil moisture study at the 1km spatial resolution, while the 20 km by 20 km area, which is shown as the small red square inside the 60 km area, is specifically used for soil moisture and irrigation study at a better spatial resolution. The Urgell climate is typically Mediterranean, with continental influence; it is mild in winter and warm in summer, with a very dry season in summer and two rainy seasons in autumn and spring (Escorihuela and Quintana-Seguí, 2016). The average annual precipitation is around 376 mm (347 mm in 2015, 385 mm in 2016, and 397 mm in 2017).

CHAPTER 3. STUDY AREA AND DATASETS

The study area is mainly irrigated using two irrigation methods: inundation, which is the old irrigation method, and sprinkler or drip, which are the new irrigation systems. The status of the irrigation is mostly static from year to year. Fields that have access to water are always irrigated, and fields that do not have access to water are not. In the old irrigated district, the open channel leads water into the agricultural fields, which causes vegetation flourishing in this area (recognized as being greener areas within the region). Surrounding the irrigated area, the land is much drier without irrigation. Currently, a new irrigation Segarra-Garrigues (SG) system is being carried to convert most of the current dry-land fields into irrigated fields, surrounding the old irrigated system, shown as the blue polygon in Figure 3.2. SG system uses the sprinkler or dripping irrigation, however, most farmers have not yet installed this irrigation system, which means that the SG systems can still be regarded as dry land [START_REF] Fontanet | The value of satellite remote sensing soil moisture data and the dispatch algorithm in irrigation fields[END_REF]. Only scatter irrigation fields are located in the middle and north of the study area. The irrigation period mainly occurs in summer, from May to September, and the frequency depends on the irrigation district. In the old district, the irrigation frequency is mainly every two weeks, and in the new area it may be daily, but it depends on the farmers. In the fields with the new irrigation method, the drippers are typically used to irrigate fruit trees, and the sprinklers are typically used to irrigate crops. The croplands dominate most parts of the area, including wheat, corn, and alfalfa. The croplands include fields irrigated using inundation, sprinkler, or subsurface drippers, and rainfed fields. The productive trees are mainly fruit trees, olive trees, or vineyards. Forest and wild trees are spreading throughout the whole area, mostly in the mountainous places in the northwestern part.

Remote sensing dataset

Different spatial datasets are used in this thesis including Sentinel-1 SAR data, Sentinel-2 optical data, Sentinel-3 altimetry data, and MODIS (Moderate Resolution Imaging Spectroradiometer) optical data. Digital elevation data and SMOS (Soil Moisture and Ocean Salinity) data are also described below, which are used as auxiliary data.

Sentinel data

ESA is developing a new family of missions called Sentinels, specifically for the operational needs of the Copernicus programme. Each Sentinel mission is based on a constellation of two satellites to fulfil the revisit and coverage requirements, providing robust datasets for Copernicus Services (ESA, 2017).

These missions carry a range of technologies, such as radar and multi-spectral imaging instruments for the land, ocean and atmospheric monitoring. The data used in this thesis are Sentinel products, including Sentinel-1, 2, and 3.

• Sentinel-1

The Sentinel-1 mission provides data from a dual-polarization C-band Synthetic Aperture Radar (SAR) instrument. The constellation consists of the Sentinel-1A and Sentinel-1B satellites, with a temporal resolution of 12 days for each one. Detailed parameters about the satellites are listed in Table 3 The primary operational imaging mode, which is the IW, combined with the Wave mode, satisfies most currently service requirements, and preserves revisit performance with robust and high reliability. For continuity reasons and emerging user requirements, the SM and EW modes are also provided when needed. The SAR instrument measurement modes and their characteristics are listed in Table 3.2.

Except for the Wave mode, which is a single-polarisation mode (selectable between HH and VV), the SAR instrument supports operation in dual polarisation (selectable between HH+HV and VV+VH). The products distributed by ESA include three levels: Level-0 raw data, which is only for specific usage, Level-1 Single Look Complex (SLC) or Ground Range Detected (GRD) data, and Level-2 Ocean (OCN) data.

The Level-1 SLC data is systematically distributed but limited to specific relevant areas, while Level-1 GRD data are systematically distributed globally with multi-looked intensity. Level-2 OCN data is only for retrieving geophysical parameters of the ocean. Therefore, in this PhD study, the Level-1 GRD data is used for soil moisture study (see Chapter 4) and irrigation mapping (see Chapter 5), having a resolution of 20 m by 5 m. Both VV and VH polarization are available with Sentinel-1 data. The investigation period of Sentinel-1A is from June 2015 until September 2017, while Sentinel-1B is from September 2016 until September 2017.

The Sentinel-1 data were downloaded from Google Earth Engine and were pre-processed by the following three steps: CHAPTER 3. STUDY AREA AND DATASETS -Thermal noise removal, which improves the signal-to-noise ratio (SNR) by removing the thermal noise as a result of the random movement of electrons in the circuitry.

-Radiometric calibration, which calibrates images acquired by satellite sensors to radiance or reflectance values, allowing a more accurate assessment of ground surface properties.

-Terrain correction using SRTM DEM at 30 m, which reduces topographic effects in SAR images.

• Sentinel-2 Sentinel-2 is a polar-orbiting, multispectral high-resolution imaging mission for land monitoring to provide, for example, the imagery of vegetation, soil and water cover, inland waterways and coastal areas. Sentinel-2 can also deliver information for emergency services. Sentinel-2A was launched on 23 June 2015 and Sentinel-2B followed on 7 March 2017. Detailed parameters are listed in Table 3.3. The Sentinel-2 Multispectral Instrument (MSI) samples 13 spectral bands: four bands at 10 metres, six bands at 20 metres and three bands at 60 metres spatial resolution (details shown in Table 3.4). Only the Level-1C and Level-2A products are released to users with all data are systematically processed to Level-1C by the Payload Data Ground Segment (PDGS).

In this PhD study, Sentinel-2 Level-1C product is used to calculate Normalized Difference Vegetation Index (NDVI) which provide the vegetation information for soil moisture retrieval (see Chapter 4). Band 4 (Red) and band 8 (NIR) are used with the following equation:

N DV I = N IR -Red N IR + Red (3.1)
Band "QA60", which is a bit-mask band containing cloud mask information, is applied in order to remove areas covered by cloud.

• Sentinel-3

Sentinel-3 is an ocean and land mission based on a constellation of two satellites (Sentinel-3A and Sentinel-3B). Sentinel-3A was launched on 16 February 2016 with data available beginning in June 2016 and was followed by Sentinel-3B, which was launched on 25 April 2018. The Sentinel-3 Radar Altimeter (SRAL) instrument is the main topographic sensor used to provide water level measurements, and hence, it is used in this thesis. The detailed parameters of the Sentinel-3 SRAL are listed in Table 3.5. The SRAL tracks the surface in two different tracking modes, namely, closed loop and open loop tracking (Figure 3.3). For closed loop tracking, the altimeter range window is autonomously positioned based on an onboard near real-time (NRT) analysis of the previous SRAL waveform; in contrast, for open loop tracking, the altimeter range window is positioned using a priori knowledge of the surface height stored in the instrument in a one-dimensional along-track DEM. The tracking modes for the studied water bodies over the Ebro River basin includes both closed loop and open loop tracking modes.

Three levels of processed altimeter data are available: Level-0, Level-1 and Level-2 products. The Level-1 product is composed of data after instrumental corrections from the Level-0 product, which constitutes raw data. The Level-2 product is further processed by applying different retracking algorithms to the Level-1 waveforms to calculate the final altimeter range and the backscatter coefficient, including geophysical corrections.

Level-1 and Level-2 products are available to users with different levels of data latency related to the availability of the auxiliary or ancillary data (AVISO+, 2017):

Near-Real-Time (NRT): delivered less than 3 hours after data acquisition Slow Time Critical (STC): delivered within 48 hours after data acquisition Non-Time Critical (NTC): delivered within typically one month after data acquisition.

The NTC product is used with priority in order to derive a high precision.

In this PhD study, both Level-1 and Level-2 product are used, with Level-1 data used for the water level retrieval and Level-2 data for comparison (see Chapter 6).

MODIS data

MODIS (Moderate Resolution Imaging Spectroradiometer) is a key instrument NASA's Terra satellite (launched in 1999) and Aqua satellite (launched in 2002) as part of the Earth Observing System (EOS) Programme [START_REF] Barnes | Prelaunch characteristics of the moderate resolution imaging spectroradiometer (modis) on eos-am1[END_REF][START_REF] Justice | The moderate resolution imaging spectroradiometer (modis): land remote sensing for global change research[END_REF]. The MODIS instrument acquires data in 36 spectral bands ranging in wavelength from 0.4 µm to 14.4 µm. The spatial resolution of the products ranges from 250 m to 500 m and to 1 km, depending on the bands used, with daily time coverage.

In this study, the 250-m spatial resolution, 16-day composites of MODIS NDVI data (MOD13Q1, collection 5) were used to calculated NDVI with equation 3.1 for soil moisture study (see Chapter 4). This product is retrieved from atmosphere-corrected, daily bidirectional surface reflectance observations, using a compositing technique based on product quality assurance metrics to remove low-quality pixels. The 250-m spatial resolution is the finest available from the MODIS NDVI dataset, and the 16-day composite was selected in order to ensure a high probability of having the best quality pixel (reduced cloud effects) representing the NDVI within every 16 days.

DEM data

A DEM was used as an ancillary dataset for improving the accuracy of water level retrieval (see Chapter 6). The Shuttle Radar Topography Mission (SRTM) is an international research effort that obtains DEM data on a near-global scale from 56 • S to 60 • N (Farr et al., 2007). The SRTM provides global data at two resolutions: 1 arc-second ( 30 m) and 3 arc-seconds ( 90 m). In this thesis study, Version 2.1 of the 1 arc-second global elevation data, offering worldwide coverage of void-filled data at a resolution of approximately 30 m, is used. The vertical accuracy of SRTM DEM is about ±16 m (absolute) and ±6 m (relative) (Farr et al., 2007;Mukul et al., 2015Mukul et al., , 2017;;Rabus et al., 2003;Elkhrachy, 2017). However, the vertical accuracy of the data decreases with the increase in slope and elevation due to the presence of large outliers and voids (Mukul et al., 2017).

SMOS data

SMOS (Soil Moisture and Ocean Salinity) mission, launched in November 2009, is to observe the soil moisture over Earth's landmasses and the salinity over Earth's oceans. The SMOS satellite carries a new type of instrument called Microwave Imaging Radiometer using Aperture Synthesis (MIRAS), which measures Earth surface in protected L-band at a 1.4 GHz frequency [START_REF] Kerr | Soil moisture retrieval from space: the soil moisture and ocean salinity (smos) mission[END_REF][START_REF] Kerr | The SMOS mission : new tool for monitoring key elements of the global water cycle[END_REF]. SMOS has global coverage and a 3-day equatorial revisiting period. The top soil layer (5 cm) moisture is retrieved with an accuracy of 0.04 m 3 /m 3 (Kerr et al., 2012). The spatial resolution ranges from 35 to 55 km, depending on the incidence angle. SMOS low-resolution moisture products (SMOS Level 3 daily product), corresponding to the twoyear period of Sentinel-1 acquisitions, were used as ancillary data for soil moisture estimation (see Chapter 4).

In situ data 3.3.1 Ground soil moisture measurements

In-situ soil moisture measurements were acquired continuously over a period of several months, in two demonstration fields belonging to the new irrigation SG system: Foradada and Agramunt (Figure 3.4). Each measurement point was analysed at different depths. Table 3.6 lists the measured soil moisture and texture characteristics of the two test fields. The in situ soil moisture measurements serve as a validation dataset for the satellite retrieved soil moisture products (see Chapter 4). Since the electromagnetic (EM) waves can only measure the properties of the surface soil, the in situ moisture measurements with the depth nearest to the surface are considered. The soil moisture measurements obtained at a depth of 3 cm in the Foradada field, and at a depth of 5 cm in the Agramunt field, are used for comparison with the

Meteorological data

The meteorological data were collected from Ruralcat website [START_REF] Ruralcat | Dades agrometeorològiques[END_REF]. The precipitation along with other meteorological data including temperature and humidity were downloaded with hourly or daily resolution. The nearest meteorological station to the demonstration field is chosen as auxillary data. For the Foradada field, the Baldomar station, which is about 6 km away, is taken, while for Agramunt field, the Tornabous station, which is about 11 km, is taken. Meteorological data were used as supplementary data for the analysis of the satellite retrieved soil moisture (see Chapter 4).

SIGPAC

The Geographic Information System for Agricultural Parcels (SIGPAC) is a public registry of an administrative nature that contains information on the parcels, which is updated every year (Direcció General de Desenvolupament Rural , 2018). The supplied SIGPAC plotter contains both graphical and alphanumeric information. The graphics information includes the referenced geographic boundary of each plot of land. Alphanumeric information contains information about each one of the enclosures such as identification codes, surface area, land usage, irrigation coefficient, average slope, and others. The interface is shown in Figure 3.6, and further information can be found on the website: http://sig.gencat.cat/visors/Agricultura.html [START_REF] Darp | Mapa agricultura[END_REF]. In this study, the SIGPAC data were used as ground truth to validate irrigation mapping results with information regarding the land usage and irrigation coefficient (see Chapter 5).

SAIH Ebro

In situ data of water levels and river flows for the Ebro River basin are available in the Automatic Hydrological Information System (in Spanish Sistema Automático de Información Hidrológica -SAIH) Ebro data hub, with examples shown in Figure 3.7. SAIH Ebro is an online system providing hourly and daily hydrological information, including river gauge data, reservoir levels, rainfall amounts and temperatures, over the Ebro River basin (SAIH Ebro, 2018). The in situ data were collected from June 2016 to June 2018, used as validation for water level retrieval (see Chapter 6) and hydrological simulation involving the reservoirs (see Chapter 7).

Conclusion

This chapter presented the study area for water resources estimation, where the Ebro river basin is considered for water estimation in reservoirs and the analysis of dam influence (see Chapter 6), and the Urgell site within the Ebro river basin is considered for irrigation mapping

CONCLUSION

and soil moisture study(see Chapter 4 and 5). The Ebro River basin is within the main Hot-Spot of climate change, and irrigation is well developed there for agriculture, making it a typical area for the study of anthropogenic effects for the water cycle.

The dataset including the remote sensing data and the in situ data were presented as well. The Sentinel-1, 2, and 3 are the main satellite data used for methodology development, while DEM and SMOS data are used as supplementary data. The validation data for soil moisture retrieval comes from two demonstrative fields, which are located in Foradada and Agramunt, whilst for irrigation mapping, the SIGPAC system offers the irrigation coefficient indicating whether the field is irrigated or not. The SAIH Ebro system, which provides the water level, river flow and characteristics of plenty of reservoirs, serves as the validation data for water level retrieval and river flow simulation for the analysis of dams influences.

With all the valuable datasets, the methodology can be developed with confidence, the retrieved results can be validated, and the accuracy can be assessed.

Introduction

Soil moisture plays an important role in the water cycle. It is the direct water source for vegetation, and it indirectly influences the temperature and humidity close to the ground, via its effect on the surface water and energy balance [START_REF] Bramer | Chapter three -advances in monitoring and modelling climate at ecologically relevant scales[END_REF]. Soil moisture is also important for irrigation management. Before assessing the irrigation status, the soil moisture needs to be investigated in order to plan the irrigation events in the next step better. Soil moisture is a key variable in the rainfall-runoff process (Aubert et al., 2003). Regular evaluation of soil moisture can significantly improve flood and drought estimations (Wagner et al., 2007), since it affects the amount of water available for vegetation growth (Cook et al., 2006;Bezerra et al., 2013).

Before discussing the methodologies, the definition of the soil moisture needs to be presented. There are two ways of measurements for the amount of water or moisture in soil, either by (1) measuring gravimetric soil water content (GWC) or by (2) measuring volumetric soil water content (VWC) (Edaphic Scientific, 2018; [START_REF] Little | A comparison of three methods of soil water content determination[END_REF]. GWC is the mass of water per mass of dry soil in a given sample with equation (4.1):

GW C = (W et -Dry)/Dry (4.1)
where Wet is the weight of the soil sample from the field and Dry is the weight of the dry soil sample.

VWC is the volume of water per volume of soil defined as equation (4.2):

V W C = GW C * BD (4.2)
where BD is the bulk density of the soil sample. In this PhD study, the VWC is considered for the soil moisture assessment to be consistent with in situ probes which determine VWC by measuring the dielectric constant of the soil.

Microwave soil moisture estimations are based on the strong contrast between the dielectric properties of water ≈ 80) and dry soil (< 5) (Zhuo and Han, 2016). With the current Sentinel-1 mission, the active onboard C-band sensor offers regular temporal coverage (about five days for Europe when both A and B satellites are considered), together with a spatial resolution of 10 m.

In the last twenty years, different empirical, semi-empirical and physical models have been proposed for the retrieval of soil moisture from various sources of SAR data (such as ERS, RADARSAT, ENVISAT, TerraSAR-X, COSMO-SkyMed). At the field scale, inversion models often take into account the effects of roughness and vegetation, due to their significant influence on radar signals. The most widely used techniques for the retrieval of soil moisture from SAR data include techniques which make use of a statistical approach such as Neural Network (NN), the direct inversion of semi-empirical or physical models with the most commonly used Water Cloud Model (WCM), and the Change Detection method (CD).

The artificial NN technique involves nonlinear parameterised mapping from an input vector to an output vector (Pulvirenti et al., 2009) With the multi-temporal data available, the CD approach can be advantageously used, and several studies have applied it mainly over semi-arid areas (Thoma et al., 2006(Thoma et al., , 2004 Wagner et al., 2008). The advantage of this method is that it can retrieve soil moisture in the absence of prior knowledge of the study area. However, this approach still needs to be investigated over vegetation covered areas.

This chapter mainly aims at developing methodologies for soil moisture retrieval under dense vegetation cover at field scale, without using field calibrations. The soil moisture is retrieved from the synergetic interpretation of Sentinel-1 and Sentinel-2 data. The first proposed method is based on the interpretation of backscattering statistics from Sentinel-1 observations, using the minimum and maximum values of this parameter throughout the full period of observation, whereas the second method is based on the analysis of backscattering differences on two consecutive acquisition days. With both methodologies, the Sentinel-1 data is combined with the normalised difference vegetation index (NDVI) computed from Sentinel-2 data.

After a short description of the two CD methodologies in Section 4.2, the results and validation in two demonstrative fields at 1 km resolution are presented in Section 4.3. Better resolution results (100 m) are presented in Section 4.4, along with the article corresponding to this study integrated in Section 4.5, providing a more detailed description of the background review, the methodologies, results, as well as detailed discussions and analysis. A summary and general conclusions can be found in Section 4.6.

The author contributions can be found at the end of the published paper. As the first author of the article, Q. Gao performed the whole processing under the supervision of prof. M. Zribi and with the help of all co-authors.

Methodologies

Backscatter difference with the driest date

The radar signals backscattered by the surface can be modelled as the sum of the radar signals scattered by the bare soil and attenuated by vegetation effects, and the signals scattered by the vegetation cover. These two contributions can be expressed as:

σ 0 cover = σ 0 veg + γ 2 (θ)σ 0 soil (4.3)
where γ 2 (θ) = exp[-2τ /cos(θ)] is the two-way vegetation canopy transmissivity, θ is the incidence angle, and τ is the optical thickness parameter that depends on the type of geometrical structure and vegetation water content of the canopy (Gupta et al., 2013).

Temporal variations in soil moisture can be directly related to the dynamics of the radar signal. When radar signals are considered for the same 100 x 100m cell, and for approximately the same NDVI index, the roughness effect can be considerably reduced by computing the difference between two radar signals recorded at two dates.

For a given NDVI (retrieved from Sentinel-2 data), by taking all of the corresponding radar data into account, the minimum value of σ 0 corresponding to the driest signal can be determined for each cell. The radar signal difference for a given cell (i,j), between one radar signal at date d and the driest signal, can be written as follows:

∆σ N DV I (i,j) = σ 0 (i,j),N DV I (d) -σ 0 dry,(i,j),N DV I = H (i,j) (N DV I, M v) (4.4)
where σ 0 (i,j),N DV I (d) is the lowest backscattered signal, corresponding to the driest conditions, and computed using the Sentinel-1 time-series using the same NDVI as for the data recorded on date d (σ 0 (i,j),N DV I (d)), and H (i,j) (N DV I, M v) is a function of the NDVI and soil moisture M v in cell (i,j).

Various experimental studies have shown that a linear relationship exists between radar signal differences and changes in soil moisture (Zribi et al., 2007;Baghdadi et al., 2007a), in the case of bare soils and vegetation-covered surfaces. For a given NDVI, the radar signal difference, σ N DV I , can thus be written as:

∆σ N DV I = α(N DV I)δM v (4.5)
where ∆M v is the change in soil moisture between the date d and the date when the soil was at its driest, and where the parameter α depends on the NDVI.

When the NDVI increases, the moisture sensitivity of the signal can be expected to decrease (El Hajj et al., 2014;Baghdadi et al., 2007b), as shown in Figure 4.1. This means that the difference between surface backscattering at a given date d, and that observed on the driest date, decreases as a function of NDVI.

The maximum variation in backscattered signal (for a fixed value of NDVI), can be written as:

∆σ N DV I max = α(N DV I)δM v max = f (N DV I) (4.6)
where The soil moisture for pixel (i,j) can thus be retrieved using the following function:

∆M v max = M v max -M v min
M v(i, j, N DV I, d) = δσ N DV I i,j f (N DV I) (M v max -M v min ) + M v min (i, j, d) (4.7)

Backscatter difference between two consecutive dates

A second change detection approach is based on the difference in backscattered signals observed on two consecutive days of Sentinel-1 data (12 days). Under these conditions, the temporal change in vegetation cover is generally very small such that, for a nearly constant value of roughness and constant vegetation conditions, the difference between the backscattered signals depends mainly on the change in soil moisture (Baghdadi et al., 2011).

When the value of the NDVI increases, the radar signals' sensitivity to temporal variations in moisture decreases. This means that the absolute value of the Sentinel-1 radar signal difference decreases over two consecutive days, provided that the NDVI remains approximately stable on these two dates. In the present case, the latter parameter is taken to be the average of the NDVI values observed for the two consecutive dates. Figure 4.2 shows this change (for both negative and positive) in radar signal behaviour for successive Sentinel-1 acquisitions, as a function of NDVI. δσ N DV I max is the maximum change in a radar signal, corresponding to the maximum value of soil moisture change δM v max , for a given value of NDVI. This can be modelled by the empirical function g:

δσ N DV I max = g(N DV I) (4.8)
In the case of the maximum value of soil moisture change δM v N DV I max , the function g can be written as:

δσ N DV I max = g(N DV I) = bN DV I + δσ bare max (4.9)
where δσ bare max is the maximum radar signal difference between two consecutive measurements over bare soil, associated with the highest value of moisture change. When the NDVI is equal to zero, δσ N DV I max is equal to δσ bare max , where b is the slope of the empirical function g. This describes the decrease in radar signal sensitivity to soil moisture. We observe an approximately symmetrical result in the computed values for the upper and lower limits. For a given NDVI, the backscatter difference δσ(t 1 , t 2 ), with t 1 and t 2 being adjacent Sentinel-1 acquisition dates, is assumed to be linearly correlated with the soil moisture difference. The soil moisture difference δM v(t 1 , t 2 ) for each cell (i,j), between successive acquisition dates t 1 and t 2 , can be retrieved using the following function:

M v(i, j, t 2 ) = H(δσ(t 1 , t 2 )) + M v(i, j, t 1 ) (4.10)
where H is equal to :

H(δσ(t, t + 1)) = δσ N DV I g(N DV I) (δM v max ) (4.11)
From a starting date t 1 , which in the present case is a date corresponding to a ground measurement, an iterative calculation is used to determine the soil moisture for the following dates:

M v(i, j, t 2 ) = M v(i, j, t 1 ) + H(δσ(t 1 , t 2 )) (4.12) M v(i, j, t 3 ) = M v(i, j, t 2 ) + f (δσ(t 2 , t 3 )) (4.13)
......

Soil moisture at 1 km resolution

Before retrieving the soil moisture at 100 m resolution, which is applicable for agriculture study, 1 km resolution is tested first to evaluate the operation of the two methodologies, which are the backscatter difference with the driest date and the backscatter difference between two consecutive dates. The study was done within the 60 km by 60 km study area using three different tracks (110, 30, 132) of the Sentinel-1A. At the time the study was being conceived, the Sentinel-2 data was not yet released. In the alternative, the MODIS NDVI was used, with a resolution of 250m, which is suitable for 1 km resolution soil moisture study.

SOIL MOISTURE AT 100 M RESOLUTION

With the ground measurements available from May to August 2015 and February to October 2016 in Foradada field, and from May to October 2015 and July to November 2016 in Agramunt field, the validation is done to evaluate the performances of both methods. The Root Mean Square Error (RMSE) is 0.08 m 3 /m 3 with method 1 considering all three satellite tracks in two years and both demonstrative fields, while the RMSE is 0.07 m 3 /m 3 for method 2, as shown in Figure 4.3.

The retrieved soil moisture for both methods is mapped for our studied site. Blue colour represents wettest pixels and red colour relatively dry. All the retrieved soil moisture values are scaled from 0 to 0.3 m 3 /m 3 .

In Figure 4.4, the up two figures (a) and (b) are the results of method 1, and figures (c) and (d) are the results of method 2. A high similarity is observed between the two method products. We retrieve approximately the same moisture spatial variations for the two analysed dates. August 21st in 2015 is a very dry day, showing clearly the irrigated area. September 2nd of 2015 estimations are just after a rainfall more than 10 mm, resulting the soil moisture increased a lot in the whole area. The results at 1 km resolution show the potential of soil moisture retrieval with the two proposed methodologies. In the next step, both methods were applied at a better spatial resolution.

Soil moisture at 100 m resolution

The soil moisture retrieval at 100 m resolution was done with the same two methodologies. To reduce the calculation time, the study area was chosen at 20 km by 20 km area. Instead of MODIS NDVI, Sentinel-2 NDVI was used with a resolution of 10 m.

The retrieved soil moisture is also validated in Agramunt and Foradada fields. The RMSE for method 1 in volumetric soil moisture is 0.087 m 3 /m 3 , with a bias of approximately 0.026 m 3 /m 3 , while for method 2 the RMSE is 0.059 m 3 /m 3 , and the unbiased RMSE is 0.053 m 3 /m 3 . More details about the 100 m resolution study can be found in the published paper integrated in the next section.

Published paper

Introduction

Surface soil moisture plays an essential role in numerous environmental studies related to hydrology, meteorology and agriculture. For hydrological and agricultural applications, accurate soil moisture estimations are essential, since the hydric state of the soil is a key variable in the rainfall-runoff process [1]. Regular evaluation of this parameter can significantly improve flood and drought estimations [2], since it affects the amount of water available for vegetation growth [3,4]. In situ networks represent single point locations, and usually cover relatively short periods of observation [5], whereas the acquisitions of satellite data make it possible to continuously retrieve surface soil moisture, at regional and global scales. Various approaches have been developed for the retrieval of soil moisture, using optical, thermal infrared (TIR), and microwave (MW) sensors [6,7]. Optical sensors in the thermal spectrum are able to identify temperature differences, which can be related to surface soil moisture. Microwave soil moisture estimations are based on the strong contrast between the dielectric properties Sensors 2017, 17, 1966 2 of 21 of water (≈80) and dry soil (<5) [8]. The depth at which the moisture is sensed depends on the sensor frequency but usually does not exceed several centimeters, in order to access rootzone soil moisture a more or less complex model is needed and several approaches have been developed such as techniques based on the energy balance approach based on thermal infrared soil moisture [9] or simplified water balanced approaches [10,11].

Traditional passive remote sensing instruments can be used to determine the surface soil moisture with a temporal resolution of 2-3 days. However, these instruments, which include the European Space Agency (ESA) Soil Moisture and Ocean Salinity (SMOS) mission [12,13] and the National Aeronautics and Space Administration (NASA) Soil Moisture Active Passive (SMAP) mission [14], have a low spatial resolution (around 40 km [15]). With the current Sentinel-1 mission, the active onboard C-band sensor offers regular temporal coverage (about five days for Europe when both A and B satellites are considered), together with a spatial resolution of 10 m.

In recent decades, SAR imagery has been shown to be advantageous for the estimation of soil surface characteristics, in particular surface roughness and soil moisture [16][17][18][19][20][21][22][23]. SAR data in the L, C and X bands is widely used for soil moisture retrieval [15,16,19,[21][22][23][24][25][26][27][28][29][30][31][32][33][34], and the C-band sensor carried by Sentinel-1 has demonstrated its ability to retrieve soil characteristics over vegetation-covered surfaces [35][36][37][38][39][40][41][42]. The Sentinel-1 data can either be used to retrieve soil moisture or for downscaling SMOS or SMAP soil moisture. By using active and passive microwave data fusion method [43], it could be possible to retrieve soil moisture at a higher accuracy.

Radar remote sensing measurements of bare soil are very sensitive to the surface-layer water content, due to a pronounced increase in the dielectric constant of the soil with increasing water content [43]. In the last twenty years, different empirical, semi-empirical and physical models have been proposed for the retrieval of soil moisture from various sources of SAR data (ERS, RADARSAT, ENVISAT, TerraSAR-X, etc.). At the field scale, inversion models often take into account the effects of roughness and vegetation, due to their significant influence on radar signals. The most widely used techniques for the retrieval of soil moisture from SAR data include the Neural Network (NN) approach [17,[44][45][46][47][48][49][50][51][52][53][54], the Water Cloud Model (WCM) [39,41,[55][56][57][58][59], and the Change Detection method (CD) [11,[60][61][62][63][64][65][66][67][68][69].

The artificial NN technique involves nonlinear parameterized mapping from an input vector to an output vector [70,71]. Santi et al. [46] reported retrieved soil moisture measurements derived from ENVISAT/ASAR data, using an artificial neural network (ANN) technique. The neural network was trained using satellite backscattering coefficients and soil parameters measured during simultaneous ground-truth campaigns, characterized by an RMSE of 0.023 m 3 /m 3 . Baghdadi et al. [52] retrieved soil moisture values from C-band SAR data using the NN technique, with an RMSE close to 0.098 m 3 /m 3 in the absence of a priori information related to the soil parameters, and an RMSE of 0.065 m 3 /m 3 when a priori soil moisture data obtained over bare agricultural areas was included in the analysis.

The Water Cloud Model (WCM) approach can be used over densely vegetated areas, since it relates the backscattering coefficient to soil moisture content and the presence of vegetation. He et al. [57] estimated the soil moisture of an alpine grassland area, using the Integral Equation Method (IEM) and the modified Water Cloud Model (WCM), leading to R 2 = 0.71 and RMSE = 0.0332 m 3 /m 3 . Two-thirds of the data points derived from field surveys were used to parameterize the backscattering model, and the remainder were used for validation. Zribi et al. [41] estimated soil moisture values from C-band ASAR data using the WCM, leading to an RMSE of approximately 0.06 m 3 /m 3 over a semi-arid region. Typical input vegetation parameters include the albedo of the vegetation and the attenuation factor, both of which are difficult to define. Laboratory-based measurements of the vegetation water content can be used when a high level of accuracy is required [72].

When multi-temporal SAR data is available, the Change Detection (CD) approach can be advantageously used, in the absence of prior knowledge of the study area. Zribi et al. [68] mapped soil moisture in a semi-arid region using ASAR/Wide Swath satellite data, based on the CD approach, with an RMSE equal to 0.13 (approximately 0.035 m 3 /m 3 as volumetric moisture) over a semi-arid 4.5. PUBLISHED PAPER region. This approach makes the assumption that changes in vegetation and soil roughness have only a minor influence on variations in backscattering coefficient, which are dominated mainly by changes in the value of soil moisture [73].

In the present paper, two simple CD methodologies are applied, without field calibrations, and validated over an area characterized by a dense vegetation cover (irrigated crop fields). The soil moisture is retrieved from the synergetic interpretation of Sentinel-1 and Sentinel-2 data. The mean soil moisture is computed with a resolution of 100 m, which is compatible with agricultural applications. Despite S1 spatial resolution being around 10 m, the soil moisture is estimated a lower resolution (100 m) in order to decrease uncertainties caused by different types of heterogeneities in agricultural fields such as local changes in roughness, heterogeneities in vegetation cover etc. The first proposed method is based on the interpretation of backscattering statistics from Sentinel-1 observations, using the minimum and maximum values of this parameter throughout the full period of observation, whereas the second method is based on the analysis of backscattering differences on two consecutive acquisition days. With both methodologies, the Sentinel-1 data is combined with the normalized difference vegetation index (NDVI) computed from Sentinel-2 data. Site calibration is not mandatory for these two methods.

Our paper is organized as follows: in Section 2, the studied site and database are presented. Section 3 describes the two proposed methodologies. Section 4 presents our validation methodology, the ground measurements, and the resulting soil-moisture maps. Finally, our discussion and conclusions are presented in the last section.

Study Site and Database

Study Site

The study area covers a 20 km by 20 km area and is located in Urgell (Catalunya, Spain). The Urgell climate is typically Mediterranean, with a continental influence: mild in winter and warm in summer, with a very dry summer season and two rainier seasons in autumn and spring [74].

The average annual temperature lies between 13 and 14 • C for most of the region, where winter temperatures are low and the summer climate is rather warm, with frequent cool nights. The average yearly rainfall is less than 500 mm for most of the region, with the winters having low rainfalls, and the summers being very dry.

Sensors 2017, 17,1966 3 of 21 a minor influence on variations in backscattering coefficient, which are dominated mainly by changes in the value of soil moisture [73].

In the present paper, two simple CD methodologies are applied, without field calibrations, and validated over an area characterized by a dense vegetation cover (irrigated crop fields). The soil moisture is retrieved from the synergetic interpretation of Sentinel-1 and Sentinel-2 data. The mean soil moisture is computed with a resolution of 100 m, which is compatible with agricultural applications. Despite S1 spatial resolution being around 10 m, the soil moisture is estimated a lower resolution (100 m) in order to decrease uncertainties caused by different types of heterogeneities in agricultural fields such as local changes in roughness, heterogeneities in vegetation cover etc. The first proposed method is based on the interpretation of backscattering statistics from Sentinel-1 observations, using the minimum and maximum values of this parameter throughout the full period of observation, whereas the second method is based on the analysis of backscattering differences on two consecutive acquisition days. With both methodologies, the Sentinel-1 data is combined with the normalized difference vegetation index (NDVI) computed from Sentinel-2 data. Site calibration is not mandatory for these two methods.

Our paper is organized as follows: in Section 2, the studied site and database are presented. Section 3 describes the two proposed methodologies. Section 4 presents our validation methodology, the ground measurements, and the resulting soil-moisture maps. Finally, our discussion and conclusions are presented in the last section.

Study Site and Database

Study Site

The study area covers a 20 km by 20 km area and is located in Urgell (Catalunya, Spain). The Urgell climate is typically Mediterranean, with a continental influence: mild in winter and warm in summer, with a very dry summer season and two rainier seasons in autumn and spring [74].

The average annual temperature lies between 13 and 14 °C for most of the region, where winter temperatures are low and the summer climate is rather warm, with frequent cool nights. The average yearly rainfall is less than 500 mm for most of the region, with the winters having low rainfalls, and the summers being very dry.

CHAPTER 4. SOIL MOISTURE RETRIEVAL

More than 80% of the Urgell area is cultivated and there is very little natural vegetation, except in the mountainous areas towards the northern and southern extremities of the region. The most common crops are corn, fruit trees, wheat and alfalfa. In the mountains, the natural vegetation consists of oak forests.

An old irrigated district located in this region has an open channel (built in 1862) transporting water towards the agricultural fields, thus allowing the vegetation to flourish in this specific area. Although the land surrounding the irrigated area is much drier, and a new irrigation network is being developed to augment the coverage of the old irrigation system, its influence is not yet visible on satellite imagery. The locations of two demonstration fields (inside the new irrigation area, at Foradada and Agramunt) are shown in Figure 1.

Database

Ground Measurements

In-situ soil moisture measurements were acquired continuously (5 min sampling frequency) over a period of several months, in two demonstration fields belonging to the new irrigation district: Foradada and Agramunt (Figure 1). Each measurement point was analyzed at different depths. The precipitation data comes from the nearest meteorological station to the demonstration field. For Foradada field, the Baldomar station, which is about 6 km away, is taken, while for Agramunt field, the Tornabous station, which is about 11 km, is considered. Table 1 lists the measured soil moisture and texture characteristics of the two test fields. The Sentinel-1 satellites are equipped with C-band Synthetic Aperture Radar (SAR) instruments, providing data in dual or single polarizations. Sentinel-1 provides data with a spatial resolution of 10 m and a temporal resolution of 12 days, in both VV and VH polarizations. In the present study, signals recorded in the VV polarization were used to compute the soil moisture estimations. Only one ground track (110) was considered, for which the incidence angle was approximately 40.3 • . Previous studies have shown that VH data has only a limited potential for the estimation of soil moisture, in particular as a consequence of its high sensitivity to volume scattering, which depends strongly on the geometrical alignment and characteristics of the vegetation [75][76][77]. The Sentinel-1 satellite database corresponds to the period from July 2015 to November 2016 (Table 2). All of the Sentinel-1 data was pre-processed using the Sentinel-1 Toolbox, in three steps:

-Thermal noise removal -Radiometric calibration -Terrain correction using SRTM DEM at 30 m.

The last step is needed to average the data over 100 m pixels or cells. As discussed above, the methodologies proposed in this study were developed on the basis of this spatial resolution, which also has the advantage of eliminating speckle effects in the radar data. The Sentinel-2A satellite was launched in June 2015, and was followed by Sentinel-2B in March 2017. It is a wide-swath, high-resolution, multi-spectral imaging mission, and is designed to provide full and systematic coverage of the Earth's land surfaces [78]. The Sentinel-2 database corresponds to the period from July 2015 to November 2016 (Table 3). The Sentinel-2 data corresponds to images recorded in 13 spectral bands, with a spatial resolution of 10 m. In the present study, band 4 (Red) and band 8 (NIR) are used to calculate the NDVI:

NDVI = (NIR -Red) (NIR + Red) (1)
Band "QA60", which is a bit-mask band containing cloud mask information, is applied in order to remove areas covered by cloud. Figure 2 -

The first of these corresponds to dry (non-irrigated) land, revealing an NDVI cycle that occurs between April and July, with a low NDVI for the remaining periods of the year. This trend is confirmed for all of the pixels observed at this location. -

The second site corresponds to an irrigated area, which is characterized by a broad range of spatial and temporal variations in NDVI. (2) Sentinel-2 data

The Sentinel-2A satellite was launched in June 2015, and was followed by Sentinel-2B in March 2017. It is a wide-swath, high-resolution, multi-spectral imaging mission, and is designed to provide full and systematic coverage of the Earth's land surfaces [78]. The Sentinel-2 database corresponds to the period from July 2015 to November 2016 (Table 3). The Sentinel-2 data corresponds to images recorded in 13 spectral bands, with a spatial resolution of 10 m. In the present study, band 4 (Red) and band 8 (NIR) are used to calculate the NDVI:

(NIR -Red) NDVI = (NIR + Red) (1)
Band "QA60", which is a bit-mask band containing cloud mask information, is applied in order to remove areas covered by cloud. Figure 2 -

The first of these corresponds to dry (non-irrigated) land, revealing an NDVI cycle that occurs between April and July, with a low NDVI for the remaining periods of the year. This trend is confirmed for all of the pixels observed at this location. -

The second site corresponds to an irrigated area, which is characterized by a broad range of spatial and temporal variations in NDVI. In order to develop suitable soil moisture algorithms, a mask is used to remove high density vegetation areas with an NDVI > 0.8, corresponding to forests that are not encountered in the agricultural pixels, and low density vegetation areas with an NDVI < 0.1, corresponding to water surfaces.

Proposed Methodologies

Method 1 Description

The first method involves retrieval of soil moisture using the radar signal CD technique. This approach to soil water content estimations has already been applied to data recorded by the ERS Scatterometer over West Africa [15]. In the present study, this method was adapted to the characteristics of the Sentinel-1 observations, and the inversion algorithm was optimized to take advantage of the high repeat rate of this data. The radar signals backscattered by the surface can be modeled as the sum of the radar signals scattered by the bare soil and attenuated by vegetation effects, and the signals scattered by the vegetation cover. These two contributions can be expressed as:

σ 0 cover = σ 0 veg + γ 2 (θ)σ 0 soil ( 2
)
where γ 2 (θ) = exp[-2τ/ cos(θ)] is the two-way vegetation canopy transmissivity, θ is the incidence angle and τ is the optical thickness parameter that depends on the type of geometrical structure and vegetation water content of the canopy [79].

Temporal variations in soil moisture can be directly related to the dynamics of the radar signal. When radar signals are considered for the same 100 × 100 m cell, and for approximately the same NDVI index, the roughness effect can be considerably reduced by computing the difference between two radar signals recorded at two dates.

For a given NDVI (retrieved from S2 data), by taking all of the corresponding radar data into account, the minimum value of σ 0 , corresponding to the driest signal, can be determined for each cell. The radar signal difference for a given cell (i,j), between one radar signal at date d and the driest signal, can be written as follows:

∆σ NDVI (i,j) = σ 0 (i,j),NDVI (d) -σ 0 dry,(i,j),NDVI = H (i,j) (NDVI, Mv) (3)
where σ 0 (i,j),NDVI (d) is the backscattered signal from cell (i,j) at date d, with the corresponding NDVI computed from the (S2) optical images; σ 0 dry,(i,j),NDVI is the lowest backscattered signal, corresponding to the driest conditions, and computed using the S1 time-series using the same NDVI as for the data recorded on date d (σ 0 (i,j),NDVI (d)), and H (i,j) (NDVI, Mv) is a function of the NDVI and soil moisture Mv in cell (i,j).

As our radar database covers a period of only 1.5 years (due to the later launch date of the Sentinel-2 satellite, i.e., June 2015), it was not possible to retrieve this relationship for each value of NDVI. We thus consider NDVI classes for the computation of σ 0 dry,(i,j),NDVI , using intervals of 0.1 (0.1-0.2, 0.2-0.3, 0.3-0.4, etc.). In the present case, the NDVI over the studied agricultural site ranges between a minimum of 0.1 and a maximum of 0.8.

Various experimental studies have shown that a linear relationship exists between radar signal differences and changes in soil moisture [19,80], in the case of bare soils and vegetation-covered surfaces. For a given NDVI, the radar signal difference, ∆σ NDVI , can thus be written as:

∆σ NDVI = α(NDVI)∆Mv ( 4
)
where ∆Mv is the change in soil moisture between the date d and the date when the soil was at its driest. The parameter α depends on the NDVI. When the NDVI increases, the moisture sensitivity of the signal can be expected to decrease [22,81], as shown in Figure 3. This means that the difference between surface backscattering at a given date d, and that observed on the driest date, decreases as a function of NDVI.

The strongest variation in moisture, ∆Mv max , corresponding to the difference between the driest value (Mv min ) and the wettest conditions (Mv max ), can be written as:

∆Mv max = Mv max -Mv min (5)
Under the conditions for which ∆Mv max is defined, the maximum variation in backscattered signal (for a fixed value of NDVI), can be written as:

∆σ NDVI max = α(NDVI)∆Mv max = f(NDVI) (6)
The predicted values of backscattered signal difference, corresponding to S1 data, are shown as a function of NDVI in Figure 4. The backscattering difference calculations were carried out for all cells and all S1 acquisition dates (over a period of approximately two years).

∆σ NDVI max is modeled as [15]:

∆σ NDVI max = f(NDVI) = a NDVI + ∆σ bare max ( 7
)
When NDVI = 0, ∆σ NDVI max = ∆σ bare max , which corresponds to the maximum value of backscattering difference under the driest, bare-soil conditions.

In order to minimize the influence of noise when estimating f(NDVI), for each selected value of NDVI, we excluded the upper 1% of the corresponding values of radar signal difference, as well as all data points having a radar signal lower than -15 dB, since these are known to correspond to water [46,82] (Figure 4). As our radar database covers a period of only 1.5 years (due to the later launch date of the Sentinel-2 satellite, i.e., June 2015), it was not possible to retrieve this relationship for each value of NDVI. We thus consider NDVI classes for the computation of 0 dry,(i,j),NDVI σ

, using intervals of 0.1 (0.1-0.2, 0.2-0.3, 0.3-0.4, etc.). In the present case, the NDVI over the studied agricultural site ranges between a minimum of 0.1 and a maximum of 0.8.

Various experimental studies have shown that a linear relationship exists between radar signal differences and changes in soil moisture [19,80], in the case of bare soils and vegetation-covered surfaces. For a given NDVI, the radar signal difference, NDVI Δσ , can thus be written as:

NDVI Δσ = α(NDVI)ΔMv ( 4
)
where ∆Mv is the change in soil moisture between the date d and the date when the soil was at its driest. The parameter α depends on the NDVI. When the NDVI increases, the moisture sensitivity of the signal can be expected to decrease [22,81], as shown in Figure 3. This means that the difference between surface backscattering at a given date d, and that observed on the driest date, decreases as a function of NDVI.

The strongest variation in moisture, max ΔMv , corresponding to the difference between the driest value (Mvmin) and the wettest conditions (Mvmax), can be written as:

max max min ΔMv = Mv -Mv (5)
Under the conditions for which max ΔMv is defined, the maximum variation in backscattered signal (for a fixed value of NDVI), can be written as:

  NDVI max max Δσ = α(NDVI)ΔMv = f NDVI (6)
The predicted values of backscattered signal difference, corresponding to S1 data, are shown as a function of NDVI in Figure 4. The backscattering difference calculations were carried out for all cells and all S1 acquisition dates (over a period of approximately two years). In order to minimize the influence of noise when estimating   f NDVI , for each selected value of NDVI, we excluded the upper 1% of the corresponding values of radar signal difference, as well as all data points having a radar signal lower than -15 dB, since these are known to correspond to water [46,82]. (Figure 4). The soil moisture for each pixel can thus be retrieved using the following function: , since the ground measurements were recorded for a limited period of time. The mean S1 radar signal is estimated over a SMOS pixel (40 km × 40 km). Figure 5 plots the relationship between this mean radar signal and the SMOS moisture values, for dates that are common to both SMOS and S1 acquisitions. An approximately linear relationship is found between the values of volumetric soil moisture and the backscattered radar signal, up to Mvmax ≈ 0.32 m 3 /m 3 , following which it saturates with a constant radar signal strength of approximately -9.5 dB. This result confirms the findings of several scientific studies, which have revealed radar signal saturation for soil moisture levels in the range (0.3-0.35 m 3 /m 3 ) [21,81]. From this result, when using method 1 we consider for a cell (i,j).

          NDVI i,
The soil moisture for each pixel can thus be retrieved using the following function:

Mv(i, j, NDVI, d) = ∆σ NDVI (i,j) f(NDVI) (Mv max -Mv min ) + Mv min (i, j, d) (8)
SMOS low-resolution moisture products (SMOS Level 3 daily product), corresponding to the two-year period of S1 acquisitions, were used to estimate Mv max and Mv min , since the ground measurements were recorded for a limited period of time. The mean S1 radar signal is estimated over a SMOS pixel (40 km × 40 km). Figure 5 plots the relationship between this mean radar signal and the SMOS moisture values, for dates that are common to both SMOS and S1 acquisitions. An approximately linear relationship is found between the values of volumetric soil moisture and the backscattered radar signal, up to Mvmax ≈ 0.32 m 3 /m 3 , following which it saturates with a constant radar signal strength of approximately -9.5 dB. This result confirms the findings of several scientific studies, which have revealed radar signal saturation for soil moisture levels in the range (0.3-0.35 m 3 /m 3 ) [21,81]. From this result, when using method 1 we consider Mv max = 0.32 m 3 /m 3 . As shown shown in Figure 5, the value of Mv min is taken to be ≈ 0.05 m 3 /m 3 .

Method 2 Description

A second change detection approach is proposed in this paper. This is based on the difference in backscattered signals observed on two consecutive days of Sentinel-1 data (12 days). Under these conditions, the temporal change in vegetation cover is generally very small such that, for a nearly constant value of roughness and constant vegetation conditions, the difference between the backscattered signals depends mainly on the change in soil moisture [24].

When the value of the NDVI increases, the radar signals' sensitivity to temporal variations in moisture decreases. This means that the absolute value of the S1 radar signal difference decreases over two consecutive days, provided that the NDVI remains approximately stable on these two dates. In the present case, the latter parameter is taken to be the average of the NDVI values observed for the two consecutive dates. Figure 6 shows this change (either negative or positive) in radar signal behavior for successive S1 acquisitions, as a function of NDVI.

Method 2 Description

A second change detection approach is proposed in this paper. This is based on the difference in backscattered signals observed on two consecutive days of Sentinel-1 data (12 days). Under these conditions, the temporal change in vegetation cover is generally very small such that, for a nearly constant value of roughness and constant vegetation conditions, the difference between the backscattered signals depends mainly on the change in soil moisture [24].

When the value of the NDVI increases, the radar signals' sensitivity to temporal variations in moisture decreases. This means that the absolute value of the S1 radar signal difference decreases over two consecutive days, provided that the NDVI remains approximately stable on these two dates. In the present case, the latter parameter is taken to be the average of the NDVI values observed for the two consecutive dates. Figure 6 shows this change (either negative or positive) in radar signal behavior for successive S1 acquisitions, as a function of NDVI. δσ NDVI max is the maximum change in radar signal, corresponding to the maximum value of soil moisture change δMv max , for a given value of NDVI. This can be modeled by the empirical function g:

δσ NDVI max = g(NDVI)
Figure 7 shows the difference in radar signal between two consecutive dates, as a function of NDVI, for all cells (i,j) and all NDVI levels at the Urgell site. The radar signal difference (negative or positive) between two adjacent days decreases in absolute value, when the NDVI increases. The negative or positive radar signal differences, resulting from respectively increasing or decreasing values of soil moisture, can be seen to follow a symmetrical, linear pattern. Figure 7 shows the difference in radar signal between two consecutive dates, as a function of NDVI, for all cells (i,j) and all NDVI levels at the Urgell site. The radar signal difference (negative or positive) between two adjacent days decreases in absolute value, when the NDVI increases. The negative or positive radar signal differences, resulting from respectively increasing or decreasing values of soil moisture, can be seen to follow a symmetrical, linear pattern. Figure 7 shows the difference in radar signal between two consecutive dates, as a function of NDVI, for all cells (i,j) and all NDVI levels at the Urgell site. The radar signal difference (negative or positive) between two adjacent days decreases in absolute value, when the NDVI increases. The negative or positive radar signal differences, resulting from respectively increasing or decreasing values of soil moisture, can be seen to follow a symmetrical, linear pattern. In the case of the maximum value of soil moisture change δMv max , the function g can be written as:

δσ NDVI max = g(NDVI) = b NDVI + δσ bare max ( 9
)
where δσ bare max is the maximum radar signal difference between two consecutive measurements over bare soil, associated with the highest value of moisture change.

When the NDVI is equal to zero, δσ NDVI max is equal to δσ bare max , where b is the slope of the empirical function g. This describes the decrease in radar signal sensitivity to soil moisture. We observe an approximately symmetrical result in the computed values for the upper and lower limits. This is due to the fact that for a given value of mean soil moisture, a very similar behavior results from either a decrease or an increase in soil moisture, as these are linearly related to the radar signal.

In order to minimize the influence of noise arising from rare events, when estimating the function g(NDVI), for each selected value of NDVI we exclude the upper 1% of the corresponding values.

For a given NDVI, the backscatter difference δσ(t 1 , t 2 ), with t 1 and t 2 being adjacent S1 acquisition dates, is assumed to be linearly correlated with the soil moisture difference. The soil moisture difference δMv(t 1 , t 2 ) for each cell (i,j), between successive acquisition dates t 1 and t 2 , can be retrieved using the following function:

Mv(i, j, t 1 ) = H(δσ(t 1 , t 2 )) + Mv(i, j, t 1 )

where H is equal to:

H(δσ(t, t + 1)) = δσ NDVI g(NDVI) (δMv max )
From the ground measurement statistics, the maximum soil moisture difference between two adjacent dates of Sentinel-1 data, δMv max , is assumed to be equal to 0.15 m 3 /m 3 .

From a starting date t 1 , which in the present case is a date corresponding to a ground measurement, an iterative calculation is used to determine the soil moisture for the following dates t 1 , t 2 , t 3 , . . . :

Mv(i, j, t 2 ) = Mv(i, j, t 1 ) + H(δσ(t 1 , t 2 )) Mv(i, j, t 3 ) = Mv(i, j, t 2 ) + H(δσ(t 2 , t 3 )) ...... (11)

Results and Discussion

Results

Using ground measurements recorded in the Foradada field from May to August 2015, and from February to October 2016, and in the Agramunt field from May to October 2015, and from July to November 2016, the values of retrieved soil moisture were validated with Sentinel-1 data, using the two approaches described in the previous section. We compare the satellite estimations with surface moisture measurements obtained at a depth of 3 cm in the Foradada field, and at a depth of 5 cm in the Agramunt field.

Method 1 Validation with Ground Measurements

Figure 8 compares the ground measurements with the values of soil moisture modeled using method 1. The Root Mean Square (RMS) error in volumetric soil moisture is 0.087 m 3 /m 3 , with a bias of approximately 0.026 m 3 /m 3 . For Agramunt field, the RMSE is 0.074 m 3 /m 3 , with a bias of -0.019 m 3 /m 3 and for Foradada field, the RMSE is 0.095 m 3 /m 3 , with a bias of 0.057 m 3 /m 3 . The RMSE can be estimated more reliably by defining an unbiased RMSE [83]:

ubRMSE = E{[(Mv retrieved -E[Mv retrieved ]) -(Mv insitu -E[Mv insitu ])] 2 } ( 12
)
where E[•] is the expectation operator. The unbiased RMSE corresponding to the first method is 0.083 m 3 /m 3 , which is equal to 0.071 m 3 /m 3 for Agramunt field, and to 0.076 m 3 /m 3 for Foradada field.

It can be seen that the errors are particularly high in the case of high moisture levels. This is due to possible variations in saturation moisture levels, and/or to spatial variations in soil roughness at the studied site. The statistical analysis should be improved by using a larger number of data acquisitions from the S1 time series. This can be expected to improve calibration of the function f. Figure 9 compares the soil moisture estimations with the ground measurements, as a function of time. The soil moisture levels retrieved from the satellite data are well correlated with precipitation events: a strong increase in soil moisture is observed, following each significant rainfall event.

where E[•] is the expectation operator. The unbiased RMSE corresponding to the first method is 0.083 m 3 /m 3 , which is equal to 0.071 m 3 /m 3 for Agramunt field, and to 0.076 m 3 /m 3 for Foradada field.

It can be seen that the errors are particularly high in the case of high moisture levels. This is due to possible variations in saturation moisture levels, and/or to spatial variations in soil roughness at the studied site. The statistical analysis should be improved by using a larger number of data acquisitions from the S1 time series. This can be expected to improve calibration of the function f. Figure 9 compares the soil moisture estimations with the ground measurements, as a function of time. The soil moisture levels retrieved from the satellite data are well correlated with precipitation events: a strong increase in soil moisture is observed, following each significant rainfall event. where E[•] is the expectation operator. The unbiased RMSE corresponding to the first method is 0.083 m 3 /m 3 , which is equal to 0.071 m 3 /m 3 for Agramunt field, and to 0.076 m 3 /m 3 for Foradada field.

It can be seen that the errors are particularly high in the case of high moisture levels. This is due to possible variations in saturation moisture levels, and/or to spatial variations in soil roughness at the studied site. The statistical analysis should be improved by using a larger number of data acquisitions from the S1 time series. This can be expected to improve calibration of the function f. Figure 9 compares the soil moisture estimations with the ground measurements, as a function of time. The soil moisture levels retrieved from the satellite data are well correlated with precipitation events: a strong increase in soil moisture is observed, following each significant rainfall event.

Method 2 Validation with Ground Measurements

Figure 10 compares the ground measurements with the estimated values of soil moisture obtained with method 2. From this regression, the RMSE in volumetric soil moisture is 0.059 m 3 /m 3 , and the unbiased RMSE is 0.053 m 3 /m 3 . The RMSE is respectively equal to 0.048 m 3 /m 3 and 0.066 m 3 /m 3 for Agramunt and Foradada field, with a bias of 0.028 m 3 /m 3 and 0.026 m 3 /m 3 separately. The unbiased RMSE is 0.04 m 3 /m 3 for Agramunt field and 0.06 m 3 /m 3 for Foradada field. Figure 11 compares the moisture estimations with ground measurements, as a function of time. The soil moisture values retrieved from satellite data are also well correlated with the observed precipitation events, with the soil moisture increasing after each significant rainfall event. As both method1 and method 2 have a relatively high RMSE, the small number of ground measurements and the relatively small size of the radar signal database could explain this high error.

Method 2 Validation with Ground Measurements

Figure 10 compares the ground measurements with the estimated values of soil moisture obtained with method 2. From this regression, the RMSE in volumetric soil moisture is 0.059 m 3 /m 3 , and the unbiased RMSE is 0.053 m 3 /m 3 . The RMSE is respectively equal to 0.048 m 3 /m 3 and 0.066 m 3 /m 3 for Agramunt and Foradada field, with a bias of 0.028 m 3 /m 3 and 0.026 m 3 /m 3 separately. The unbiased RMSE is 0.04 m 3 /m 3 for Agramunt field and 0.06 m 3 /m 3 for Foradada field. Figure 11 compares the moisture estimations with ground measurements, as a function of time. The soil moisture values retrieved from satellite data are also well correlated with the observed precipitation events, with the soil moisture increasing after each significant rainfall event. As both method1 and method 2 have a relatively high RMSE, the small number of ground measurements and the relatively small size of the radar signal database could explain this high error.

Method 2 Validation with Ground Measurements

Figure 10 compares the ground measurements with the estimated values of soil moisture obtained with method 2. From this regression, the RMSE in volumetric soil moisture is 0.059 m 3 /m 3 , and the unbiased RMSE is 0.053 m 3 /m 3 . The RMSE is respectively equal to 0.048 m 3 /m 3 and 0.066 m 3 /m 3 for Agramunt and Foradada field, with a bias of 0.028 m 3 /m 3 and 0.026 m 3 /m 3 separately. The unbiased RMSE is 0.04 m 3 /m 3 for Agramunt field and 0.06 m 3 /m 3 for Foradada field. Figure 11 compares the moisture estimations with ground measurements, as a function of time. The soil moisture values retrieved from satellite data are also well correlated with the observed precipitation events, with the soil moisture increasing after each significant rainfall event. As both method1 and method 2 have a relatively high RMSE, the small number of ground measurements and the relatively small size of the radar signal database could explain this high error.

Discussion

The soil moisture at two study sites has been computed and mapped, using Equations ( 5)-( 11) and data produced by Sentinel-1 radar observations. Figure 12 provides two illustrations of moisture mapping, using methods 1 and 2, for two cases: a very dry day (21 August 2015) and a wet day (2 September 2015). All cells with NDVI < 0.1 or NDVI > 0.8, which are associated with water bodies and forests respectively, are masked out. A high similarity is observed between the products obtained with these two methods. We retrieve nearly the same spatial variations in moisture on the two analyzed dates. The first dry case clearly reveals the irrigated fields inside dry area. The second wet case shows high values of soil moisture over most of the observed area.

Figure 13 shows the difference between method 1 and method 2 for date 21 August 2015 and date 2 September 2015. A limited difference is illustrated for the two dates between the two methods. The highest differences correspond to high vegetation density covers. Figure 14 plots the variation in RMSE between methods 1 and 2, as a function of NDVI. This is estimated with a sliding NDVI window, with an NDVI width = 0.2, for NDVI values lying in the range between 0.1 and 0.8. The RMSE can be seen to increase with NDVI. In practice, a high vegetation density can significantly attenuate the signals, thus leading to correspondingly higher errors in soil moisture estimation. Areas with higher vegetation cover are with higher uncertainties. Rapid vegetation change and soil properties such as surface roughness will contribute to uncertainties as well. As change detection approaches, these two methods are very applied operationally since the ground measurements are not prerequisite and that they can be improved with the size of time series.

Discussion

The soil moisture at two study sites has been computed and mapped, using Equations ( 5)- (11) and data produced by Sentinel-1 radar observations. Figure 12 provides two illustrations of moisture mapping, using methods 1 and 2, for two cases: a very dry day (21 August 2015) and a wet day (2 September 2015). All cells with NDVI < 0.1 or NDVI > 0.8, which are associated with water bodies and forests respectively, are masked out. A high similarity is observed between the products obtained with these two methods. We retrieve nearly the same spatial variations in moisture on the two analyzed dates. The first dry case clearly reveals the irrigated fields inside dry area. The second wet case shows high values of soil moisture over most of the observed area.

Figure 13 shows the difference between method 1 and method 2 for date 21 August 2015 and date 2 September 2015. A limited difference is illustrated for the two dates between the two methods. The highest differences correspond to high vegetation density covers. Figure 14 plots the variation in RMSE between methods 1 and 2, as a function of NDVI. This is estimated with a sliding NDVI window, with an NDVI width = 0.2, for NDVI values lying in the range between 0.1 and 0.8. The RMSE can be seen to increase with NDVI. In practice, a high vegetation density can significantly attenuate the signals, thus leading to correspondingly higher errors in soil moisture estimation. Areas with higher vegetation cover are with higher uncertainties. Rapid vegetation change and soil properties such as surface roughness will contribute to uncertainties as well. As change detection approaches, these two methods are very applied operationally since the ground measurements are not prerequisite and that they can be improved with the size of time series.

Conclusions

In this study, two inversion approaches are developed for the interpretation of high repeat frequency Sentinel-1 radar data in synergy with Sentinel-2 optical data. Change detection techniques in proposed methodologies, are validated with ground measurements carried out in two demonstration fields. The estimated (volumetric) RMS soil moisture errors are approximately 0.087 m 3 /m 3 for method 1 and 0.059 m 3 /m 3 for method 2. Both methods are found to predict soil moisture variations that are well correlated with rainfall events. Method 1 models the backscattering difference with the driest value, whereas method 2 is based on the difference between radar signals observed on two consecutive dates, meaning that the radar signals are influenced by much smaller changes in vegetation. Method 2 is found to be more robust than method 1, since it does not require searching for the minimum value in each pixel, which can introduce larger errors under extreme local conditions. The backscattered radar contributions produced by the vegetation are small in the case of method 2. However, as the retrieved value of soil moisture depends on the soil moisture determined at an earlier date with this method, the iterative process can lead to the accumulation of errors. Comparing to other types of inversion algorithms (e.g., NN or calibrated model inversion), both of these methods allow soil moisture to be estimated, with no need for calibrations based on ground measurements, and have led to the production of similar, 100 m resolution soil moisture maps of the study area. SMOS data is used for limiting the maximum soil moisture retrieved by satellite. However, the interest of SMOS (or SMAP or other low-resolution soil moisture sensor) is that it is available globally and needs no local maintenance, which makes our method applicable globally in contrast with methods that require in-situ data such as NN or calibrated model inversion.

These results demonstrate the potential of Sentinel-1 data for the retrieval of 100 m (or even better) resolution soil moisture. Both methods can be applied to any vegetation-covered area for which time-series of SAR and optical data have been recorded. In future, the statistical analysis should be improved by using a larger number of data acquisitions from the S1 time series.

In the present study, data derived from the VV polarization was analyzed, since it is more sensitive to soil conditions. However, Sentinel-1 provides data in both VV and VH polarization modes, and it is planned to include VH polarization analyses in future studies, since this operational mode is highly sensitive to the influence of vegetation, and can be used to discriminate between the effects of vegetation.

Conclusions

In this study, two inversion approaches are developed for the interpretation of high repeat frequency Sentinel-1 radar data in synergy with Sentinel-2 optical data. Change detection techniques in proposed methodologies, are validated with ground measurements carried out in two demonstration fields. The estimated (volumetric) RMS soil moisture errors are approximately 0.087 m 3 /m 3 for method 1 and 0.059 m 3 /m 3 for method 2. Both methods are found to predict soil moisture variations that are well correlated with rainfall events. Method 1 models the backscattering difference with the driest value, whereas method 2 is based on the difference between radar signals observed on two consecutive dates, meaning that the radar signals are influenced by much smaller changes in vegetation. Method 2 is found to be more robust than method 1, since it does not require searching for the minimum value in each pixel, which can introduce larger errors under extreme local conditions. The backscattered radar contributions produced by the vegetation are small in the case of method 2. However, as the retrieved value of soil moisture depends on the soil moisture determined at an earlier date with this method, the iterative process can lead to the accumulation of errors. Comparing to other types of inversion algorithms (e.g., NN or calibrated model inversion), both of these methods allow soil moisture to be estimated, with no need for calibrations based on ground measurements, and have led to the production of similar, 100 m resolution soil moisture maps of the study area. SMOS data is used for limiting the maximum soil moisture retrieved by satellite. However, the interest of SMOS (or SMAP or other low-resolution soil moisture sensor) is that it is available globally and needs no local maintenance, which makes our method applicable globally in contrast with methods that require in-situ data such as NN or calibrated model inversion.

These results demonstrate the potential of Sentinel-1 data for the retrieval of 100 m (or even better) resolution soil moisture. Both methods can be applied to any vegetation-covered area for which time-series of SAR and optical data have been recorded. In future, the statistical analysis should be improved by using a larger number of data acquisitions from the S1 time series.

In the present study, data derived from the VV polarization was analyzed, since it is more sensitive to soil conditions. However, Sentinel-1 provides data in both VV and VH polarization modes, and it is planned to include VH polarization analyses in future studies, since this operational mode is highly sensitive to the influence of vegetation, and can be used to discriminate between the effects of vegetation.

projects. The authors wish to thank the technical teams for their support and contributions to the successful ground campaigns.

SUMMARY AND CONCLUSIONS

71

Summary and Conclusions

In this chapter, two inversion approaches are developed for soil moisture retrieval using Sentinel-1 radar data in synergy with Sentinel-2 optical data. Change detection techniques in proposed methodologies are validated with ground measurements carried out in two demonstration fields. Both methods are found to predict soil moisture variations that are well correlated with rainfall events. Method 1 models the backscattering difference with the driest value, whereas method 2 is based on the difference between radar signals observed on two consecutive dates, meaning that the radar signals are influenced by much smaller changes in vegetation. Method 2 is found to be more robust than method 1, as it does not require searching for the minimum value in each pixel, which can introduce larger errors under extreme local conditions. The backscattered radar contributions produced by the vegetation are small in the case of method 2. However, as the retrieved value of soil moisture depends on the soil moisture determined at an earlier date with this method, the iterative process can lead to the accumulation of errors. Both of these methods allow soil moisture to be estimated, with no need for calibrations based on ground measurements and have led to the production of similar, 100 m resolution soil moisture maps of the study area.

This study only used VV polarization with better sensitivity to the soil moisture compared with VH polarization. However, VH polarization will be included in future research to better characterize the role of vegetation in soil moisture retrieval using SAR data. Also, since vegetation effects on radar signal are not only due to the amount of biomass but also the vegetation structure, a better vegetation index which can reveal both characters need to be investigated instead of NDVI. A land cover map could be considered along with the vegetation index to better characterize the relation between backscattered radar signals and vegetation characters. In the present study, the ground measurements are taken in different depth (3 cm in Foradada field and 5 cm in Agramunt field), the impacts of sensor depth need to be evaluated in the future.

Besides the possibilities of improving and enriching the study, the results in the current research already demonstrated the potential of Sentinel-1 data for the retrieval of 100 m (or even better) resolution soil moisture. Both methods can be applied to any vegetation-covered area for which time-series of SAR and optical data have been recorded.

With the soil moisture map available at a good spatial and temporal resolution, irrigation events can be managed in a better way.

Introduction

Irrigation, as one of the most important human interventions of the water cycle, is of great interest in evaluating the food yield and water usage. The irrigated fields are mostly stable, with fields having access to water always irrigated, but fields that do not have access to water are not. Irrigation events in irrigated fields directly increase the soil moisture, which is the direct water source of vegetation. Therefore, the spatially explicit irrigated area information needs to be investigated to provide a clear image for the agricultural water usage and irrigation management.

This chapter is devoted to developing the methodology of irrigation mapping using remote sensing data. The most commonly used methodologies include the optical multi-bands classification method, the vegetation index method, and the radar-based approach. The first two methodologies, which are based on optical remote sensing, are limited by weather conditions and may not be applicable over frequently cloudy regions. SAR remote sensing technology provides a new potential tool for monitoring irrigation fields under any weather conditions. However, few studies have used SAR data alone to do irrigation mapping, except for paddy rice fields (Ribbes and Toan, 1999; Shao et al., 2001). Since irrigation increases the soil moisture in the fields and helps cultures flourish, the SAR data that is affected both by soil moisture and vegetation have the ability to discriminate the irrigated and non-irrigated fields.

Sentinel-1 SAR mission data, with its high spatial and temporal resolution, allows to distinguish irrigated and non-irrigated areas at the field scale. In this study, both the VV and VH polarizations of multi-temporal SAR data are used to map irrigated fields over an agricultural site in Urgell, Catalunya. The database consists of Sentinel-1A from June 2015 until September 2017, and Sentinel-1B from September 2016 until September 2017. By analysing the mathematical metrics derived from the SAR backscatter time series, including the mean value, the signal variance, the correlation length, irrigated and non-irrigated fields can be separated. Furthermore, the SAR characteristics of backscatter show a difference for different irrigation types; based on this aspect, the irrigated trees and irrigated crops can be separated as well. The classification CHAPTER 5. IRRIGATION MAPPING is done both with the Support Vector Machine (SVM) and Random Forest (RF) machine. The validation using the ground truth from SIGPAC over the whole study area shows good overall accuracy with 81.08% for SVM and 82.2% for RF classification.

In this chapter, a short description of the methodology for classification is presented in Section 5.2. Then in Section 5.3, the published article corresponding to this study can be found, providing a more detailed literature review, the methodologies, results, as well as detailed discussions and analysis. Finally, the summary and conclusions are presented in Section 5.4.

The author contributions can be found at the end of the published paper. As the first author of the article, Q. Gao performed the whole processing under the supervision of Prof. M. Zribi, and with the help of all co-authors.

Methodology

This study used multi-temporal SAR data with field segmentation information from SIGPAC (the Geographic Information System for Agricultural Parcels). The backscatter intensities are averaged within each field to reduce the speckle effects and at the same time maintain the field scale. From the Sentinel-1 time series for each field, the statistics and metrics, including the mean value, the variance of the signal, the correlation length, and the fractal dimension, are analysed. Then, based on the SVM or RF, the classification is done over the whole study area.

Four metrics of backscatter time series both in VV and VH polarization were analysed over three groups of selected fields including irrigated croplands, irrigated trees, and non-irrigated crops:

• Mean value: shows the average scale of the backscatter intensity. In principle, the value should be higher for irrigated than non-irrigated fields in VV polarization, with higher soil moisture in irrigated fields makes more contribution to the backscatter intensity. For VH polarization, which is more sensitive to vegetation, the mean value should be higher over irrigated fields than non-irrigated fields, since irrigated helps vegetation to flourish.

• Variance: measures how far a dataset is spread out. Its value over the forest and urban areas are expected to be lower than the agricultural areas, with limited temporal variations.

• Correlation length: the temporal interval in which the autocorrelation function decays to half of the power the fastest, providing the information regarding the temporal variance of the backscatter signal, with limited values for fast-changing fields (e.g., irrigated fields) in principle.

• Fractal dimension: measures the self-similarity of the backscatter time series. Its value is expected to be higher over the irrigated field in principle, since irrigation should bring faster change to SAR signals. However, the results show no big differences between irrigated and non-irrigated fields. The fractal dimension was analysed but not used for classification.

With the metrics being analysed, SVM was used for classification. First, the irrigated and nonirrigated fields were separated using mean value and variance; then, the irrigated fields were further separated into irrigated trees and irrigated crops using variance and correlation length.

With the classification functions derived in this step, a tree classification method was applied over the whole study area. RF is also tested with half of the fields for learning and the other half for validating. A detailed methodology description can be found in the published article integrated with Section 5.3.

Introduction

Irrigated agriculture is essential for the global food yield. In the past 40 years, global agricultural production has more than doubled, while the cropland has only increased by 12%, which reveals that irrigation has made a great contribution [1][2][3]. With the constant increase of the global population, the shortage of water resources, and the mismatch between irrigation needs and the actual amount of water used for irrigation [4,5], irrigation needs to be better planned in order to fulfill the high demand for food and water. Spatially explicit irrigated area information is needed for cropland irrigation management [6]. The amount of irrigation surfaces, which is of great importance to water resource management, is still unclear. Remote sensing technology leads to a new direction for mapping irrigated areas to better support water resources and agricultural development [7]. However, studies using remote sensing to map irrigated fields remain relatively rare [3]. With the diverse range of irrigated field size and scattered distribution, it is difficult to map irrigation fields by satellite remote sensing because of its relatively coarse spatial resolution in comparison to field scale.

The most commonly used methodologies are the optical multi-bands classification method, the vegetation index method, and the radar-based approach. Several studies [8][9][10][11][12][13][14][15] have shown the capability of optical multi-bands remote sensing for irrigation mapping. Thiruvengadachari [14] demonstrated that areas irrigated by surface water could be distinguished from those irrigated with groundwater by using single-date Landsat imagery, based on supplementary information regarding major surface irrigation projects, canal network maps, drainage patterns, and recorded groundwater utilization. Even though the single-date image can be used for distinguishing irrigation fields, it is not always reliable, since single-date analysis in visible cropping intensity often does not take into account planting dates that vary from year to year [3]. Therefore, multi-temporal analysis provides a better potential to define irrigated fields [16]. Thenkabail et al. [8] used Moderate Resolution Imaging Spectroradiometer (MODIS) time-series data to generate Land Use/Land Cover (LULC) and a map of irrigated areas for the Ganges and Indus river basins.

Various vegetation indices, such as the NDVI (Normalized Difference Vegetation Index), the NDWI (Normalized Difference Wetness Index), and the GVI (Green Vegetation Index), derived from multi-bands satellite data, are proven to be able to map irrigated areas [10,12,[17][18][19][20]. Boken et al. [18] demonstrated the potential of NOAA-AVHRR (Advanced Very High Resolution Radiometer) for estimating irrigated areas using NDVI and the Vegetation Health Index (VHI) with coefficients of determination (R2) of approximately 0.49 and 0.80, respectively. Xiao et al. [19] developed a paddy rice mapping algorithm that uses the time series of three vegetation indices, namely, the Land Surface Water Index (LSWI), the Enhanced Vegetation Index (EVI), and the NDVI, derived from MODIS images.

Studies combining the two methods above increased the spatial and temporal resolution for irrigation mapping. Gumma et al. [21] developed a decision tree approach using Landsat 30-m one-time data fusion with MODIS 250-m time series data. Fuzzy classification accuracy assessment for the irrigated classes varied between 67-93%.

Methods based on optical data rely heavily on weather conditions. For areas with frequent cloud cover, these methods may not be adaptable. The availability of Synthetic Aperture Radar (SAR) data offers a new potential for irrigation monitoring by providing the ability to observe under any weather conditions. The radar remote sensing measurements of soil are very sensitive to the water content in the surface layer, due to the pronounced increase in the soil dielectric constant with increasing water content [22][23][24][25]. Although primarily affected by soil moisture, active radar backscatter is also influenced by vegetation and surface roughness [26][27][28][29][30][31][32], and the proportions depend on the different polarization modes. The Sentinel-1 mission also proved that it can be used to retrieve soil moisture under vegetation cover [33,34]. Studies also demonstrated that radar data can provide unique characteristics of irrigated croplands such as rice fields [35,36]. Ribbes and Toan [35] demonstrated that the radar backscatter coefficient of rice fields had a significant temporal variation, and that this variation can be used to identify paddy rice fields. However, rice fields are special since they are inundated, which can be distinguished more easily by SAR than general irrigated croplands. Few studies have used SAR data alone to do irrigation mapping, except for paddy rice fields. Since irrigation events change the soil moisture in the fields and help cultures flourish, the SAR data that is affected both by soil moisture and vegetation should have the ability to discriminate the irrigation fields. Irrigation mapping using SAR data needs to be studied to enrich remote sensing applications in agricultural and hydrological fields. Since irrigation is a time dynamic activity, more multi-temporal datasets are needed.

In this study, multi-temporal SAR data is used to map irrigated crops, irrigated trees, and non-irrigated fields. Sentinel-1 SAR mission data, with its high spatial and temporal resolution, allows for more possibilities for distinguishing irrigated areas at the field scale. Unlike optical data, which is restricted by cloud coverage, SAR data can work under any weather conditions. By analyzing the metrics derived from the SAR backscatter time series, irrigated and non-irrigated fields can be separated. Furthermore, the SAR characteristics of backscatter show a difference for different irrigation types; based on this aspect, the irrigated trees and irrigated crops can be separated as well.

In this paper, we introduce four metrics, including the mean value, the signal variance, the correlation length, and the fractal dimension derived from Sentinel-1 multi-temporal data. By combining different metrics, the irrigated trees, irrigated crops, and non-irrigated fields are separated. In Section 2, the studied site, the database, and the land use are presented. Section 3 describes the methodologies and metrics that were studied. Section 4 shows the results of the classification mapping and comparison with ground truth. Section 5 contains a discussion about the results. Finally, the conclusions are presented in the last section.

Database and Study Area

Database

Sentinel-1 Data

The Sentinel-1 mission provides data from a dual-polarization C-band Synthetic Aperture Radar (SAR) instrument. In this study, we use both Sentinel-1A and Sentinel-1B backscatter data, with a resolution of 10 m and a temporal resolution of 12 days, for each satellite. Previous studies have shown that VV (vertical-vertical) polarization data, in comparison to VH (vertical-horizontal) polarization, show high sensitivity to soil moisture [37][38][39][40][41][42]. VH data, in its turn, has a higher sensitivity to volume scattering, which depends strongly on the geometrical alignment and characteristics of the vegetation. Thus, VH data has a limited potential for the estimation of soil moisture compared to VV data, but higher sensitivity to vegetation [38,43,44]. Since irrigation makes a difference for vegetation density, and in order to make full use of the satellite data, we use both VV and VH polarizations for our analysis. The period of Sentinel-1A is from June 2015 until September 2017, while Sentinel-1B is available from September 2016 until September 2017. The database is shown in Figure 1. All of the Sentinel-1 data was pre-processed using the Sentinel-1 Toolbox, in three steps:

Ground Truth (SIGPAC)

The Geographic Information System for Agricultural Parcels (hereinafter SIGPAC) is a public registry of an administrative nature that contains information on the parcels [45]. The supplied SIGPAC plotter contains both graphical and alphanumeric information. The graphic information

78

CHAPTER 5. IRRIGATION MAPPING includes the referenced geographic boundary of each plot of land. Alphanumeric information contains information about each one of the enclosures such as identification codes, surface area, land usage (updated annually), irrigation coefficient, average slope, and others. Information can be found on the website: http://sig.gencat.cat/visors/Agricultura.html [46]. In this study, the SIGPAC data is used as ground truth to validate our classification results with information regarding the land usage and irrigation coefficient.

Study Area

The study area covers a 20 km by 20 km area, located in Urgell, Catalonia, Spain (Figure 2). The Urgell climate is typically Mediterranean, with continental influence; it is mild in winter and warm in summer, with a very dry season in summer and two rainier seasons in autumn and spring [47]. The average annual precipitation is around 376 mm (347 mm in 2015, 385 mm in 2016, and 397 mm in 2017). Precipitation events at the Baldomar meteo station, which is the nearest one to our study area, are illustrated in Figure 3. The study area is mainly irrigated using two irrigation methods: inundation, which is the old irrigation method, and sprinkler or drip, which are the new irrigation systems. The status of the irrigation is mostly static from year to year. Fields that have access to water are always irrigated, and fields that do not have access to water are not. In the old irrigated district, the open channel leads water into the agricultural fields, which causes vegetation flourishing in this area (lower left of the study area). Surrounding the irrigated area, the land is much drier without irrigation. Currently, a new irrigated system is being developed, surrounding the old irrigated system, using sprinkler or dripping irrigation; scatter irrigation fields are located in the middle and north of the study area. The irrigation period mainly occurs in summer, from May to September, and the frequency depends on the irrigation district. In the old district, the irrigation frequency is mainly every two weeks, and in the new area it may be daily, but it depends on the farmers.

In the fields with the new irrigation method, the drippers are always used to irrigate fruit trees, and the sprinkler is always used to irrigate crops. In our study, the fields are separated into irrigated trees, irrigated crops, and non-irrigated fields in the end.

Figure 2 shows the land-use map from SIGPAC in our study area. The croplands dominate most parts of the area, including wheat, corn, and alfalfa. The croplands include fields irrigated using inundation, sprinkler, or subsurface drippers, and rainfed fields. The productive trees are mainly fruit trees, olive trees, or vineyards. Forest and wild trees are spreading throughout the whole area, mostly in the mountainous places in the northwest part. In total, there are 22,775 fields of croplands and 3659 fields of productive trees that are used for the mapping and validation.

Methodology

The first step is the data pre-processing. Second, we use metrics of a backscatter signal time series over the selected fields to perform modeling to build the classification function. Third, irrigation mapping is performed using the metrics maps, the classification function, and the fields boundary information. The accuracy assessment will be explained in the next section. The overall workflow using the Support Vector Machine (SVM) is shown in Figure 4.

Data Pre-Processing

To build the classification method, we first selected 10 areas for irrigated crops, five areas for irrigated trees, and 10 areas for non-irrigated crop lands, as shown in Figure 5, to perform learning and modeling. The surface size of the selected areas varies from 3 ha to 25 ha. All of the pixels within the selected reference areas are pre-processed by averaging every 10 by 10 pixels, starting from the left-top corner into 100 m by 100 m cells in order to reduce the speckle effect.

Analyzed Metrics

The mathematical statistics of the backscatter time series for the selected areas, including irrigated trees, irrigated croplands, and non-irrigated crops, will be analyzed. In principle, the backscatter values are different for irrigated and non-irrigated areas and for different irrigation methods, since different water content in the soil creates differences in soil dielectric constant, which relates to the backscatter characters. Furthermore, the irrigation makes a difference for vegetation density, which affects the backscatter coefficient, especially for VH polarization. Additionally, the changing pattern and other mathematical statistics will be analyzed for irrigated trees, irrigated crops, and non-irrigated crops.

Analysis of the temporal-related metrics is the first step in the analysis of the possibility that the three types of classes will be separated. The analyzed metrics include the mean value of the backscatter time series, the temporal variance of the signal, the signal correlation length, and the fractal dimension. First, we analyze these metrics solely within the selected areas over irrigated croplands, irrigated trees, and non-irrigated crops. Cells with hundred-meter resolution fell into the selected fields are considered for the analysis.

Mean Value of σ •

The mean value, which is the basic statistical parameter of a time series, shows the average scale of the data. In the first analysis, the mean value of the backscatter time series is taken for each 100-m cell within the selected irrigated and non-irrigated areas. Analysis is based on the whole period of our data acquisition, as shown in Figure 1. The mean value reflects the general level of signal intensity for different land use.

For VV polarization, the backscatter intensity should be higher for irrigated than non-irrigated fields, with more water content in the fields because the increase of water increases reflectivity and thus, the backscatter value. For VH polarization, the sensitivity to water content in the soil is less pronounced, while the impact of vegetation increases. Figure 6 illustrates the case of VV and VH polarization signal dynamics through radar signal histograms; the approximate mean value of the selected irrigated cells for VV polarization is -11.724 dB, as shown in Figure 6a, and -14.137 dB for non-irrigated cells, as shown in Figure 6b. The mean values are smaller for VH polarization than VV polarization, and smaller for non-irrigated fields than irrigated fields.

Signal Variance

Variance is the expectation of the squared deviation of a random variable from its mean. Informally, it measures how far a set of numbers is spread out from its average value. The variance value over the forest and urban areas are expected to be lower over the agricultural areas, with limited temporal variations. Figure 6 shows that the variance values are smaller for irrigated areas compared

82

CHAPTER 5. IRRIGATION MAPPING with non-irrigated areas, while VH polarization shows in average higher values than VV polarization, but similarly to VV, it is also smaller for the irrigated area.

Signal Correlation Length

Autocorrelation is a measure of the internal correlation within a time series. Autocorrelation is the correlation of a signal (µ) with a delayed copy of itself as a function of delay (τ). It is a way of measuring and explaining the internal association between observations in a time series.

ρ(τ) = E [ƒ(µ) ƒ (µ -τ)], (1)
Autocorrelation length represents the temporal interval in which the autocorrelation function decays to half of the power the fastest. It provides the information regarding the temporal variance of the backscatter signal, with limited values for fast-changing fields (e.g., irrigated fields) in principle. In our analysis, Edelson and Krolik's Discrete Correlation Function (DCF) is used [48]. Figure 7 gives two examples for irrigated cells and non-irrigated cells for VV and VH polarization; an exponential function [y = A exp(x/B)] is used to fit the temporal discrete correlation function. The correlation length is larger for irrigated cells in this case, since the existence of temporal variations is linked principally to soil moisture dynamic, and for non-irrigated areas, local temporal variations due to roughness and vegetation could generate reduction on correlation length; however, the overall results show that a large range of autocorrelation length varies from less than 10 days to more than 60 days for both irrigated and non-irrigated cells.

Fractal Dimension

To quantitatively characterize backscatter dynamics, fractal geometry can be used to extract robust features hidden in the fluctuations. Fractals have the property of self-similarity [49][50][51]. Our fractal characterization is based on the power spectrum dependence of a fractional Brownian motion. We determine the power spectrum of the backscatter time series by computing its Fourier transform. The scaling behavior of the data is revealed by a power-law dependence of the spectrum as a function of frequency [52,53]:

S(f) ≈ 1/f β
The fractal dimension is derived from the slope β of a least-squares regression linear fit to the data points in the log-log plot of the power spectrum, leading to the dimension [50]: D = 7/2β/2.

PUBLISHED PAPER

Over irrigated fields, the fractal dimensions are expected to be higher in principle, since irrigation should bring faster change in the irrigation period to SAR signals. However, the results show that the values are almost the same for irrigated and non-irrigated fields, since the backscatter is influenced by many factors, including soil roughness and precipitation. Figure 8 provides two examples of power spectral density plots of irrigated cells and non-irrigated cells for VV and VH polarization. The difference of slopes is not evident for the separation of irrigated and non-irrigated fields. In the implementation of our method, we discard fractal dimensions, which are not robust enough for classification.

Modeling

Support Vector Machine

A Support Vector Machine (SVM) is a typical learning machine for two-group classification problems [54]. SVM is a discriminative classifier that is formally defined by a separating hyperplane. In other words, given labeled training data (supervised learning), the algorithm builds an optimal hyperplane, which can be used to categorize a new dataset. In a two-dimensional space, the hyperplane becomes a line dividing a plane into two classes, which lay in each part of the plane. Our method of separating non-irrigated fields, irrigated crops, and irrigated trees is based on the linear Support Vector Machine.

In our process, first, we separate irrigated and non-irrigated fields using mean value and variance, then, the irrigated fields are further separated into irrigated trees and irrigated crops using variance and correlation length. The reason to separate two classes each time instead of three is that the mean values of irrigated trees and irrigated crops are a bit mixed, since the backscatter difference between the flouring crops and trees are not evident, and the correlation length for non-irrigated crops varies a lot depending on the soil roughness and vegetation status, which makes it difficult to be separated from irrigated fields.

Classification Function

The classification functions are derived using SVM as shown in Figures 9 and10. First, a function for the irrigated fields and non-irrigated fields classification is built by the mean value and variance. Second, a function for irrigated crops and irrigated trees classification is built by the variance and correlation length.

84

CHAPTER 5. IRRIGATION MAPPING Considering the selected irrigation and non-irrigated masks, we analyzed the mean value and the variance of the backscatter time series for each 100 m by 100 m cell. The results are shown in Figure 9. The mean values can be easily distinguished for irrigated and non-irrigated cells, with large values for irrigated cells. Over the selected fields, the mean variance for VV polarization is 3.96 dB for irrigated areas and 5.96 dB for non-irrigated areas. The function of the VV polarization increases, while for the VH polarization, it decreases, because VV is more sensitive to soil moisture and VH is more sensitive to vegetation.

Based on the two metrics, the irrigated and non-irrigated fields can be separated by the function of the mean value and the variance, using both VV and VH polarization. Then, for irrigated areas, the SVM derived linear function of correlation length and variance is used for the separation of the irrigated trees and irrigated crops (as shown in Figure 10). The characteristics of its self-correlation vary substantially, as shown in Figure 10. The effect of vegetation coverage limits the variance for the case of trees. Apart from this, we see a limited increase of variance with correlation length.

By using the function of the correlation length and variance, we can separate irrigated crops and irrigated trees.

Tree Classification

Our study uses a tree classification, considering metrics from both VV and VH polarization, including mean value, variance, and autocorrelation length. Over the whole study area, the two linear SVM functions in Figure 9 are first used to separate irrigated fields and non-irrigated fields. Second, the two linear SVM functions in Figure 10 are used to separate irrigated crops and irrigated trees.

Random Forest (RF) Classification

Random Forest (RF) is an ensemble learning method for classification, regression, and other tasks, that operates by constructing a multitude of decision trees at training time and outputting the class that is the type of the classes (classification) or mean prediction (regression) of the individual trees [55,56].

Our study targets the classification process using half of the fields for learning, and the other half for validating, with the mean value, the variance, and autocorrelation length of both VV and VH polarization as input.

Results and Validation

With the geoinformation from SIGPAC, we know the boundary of all of the fields in our study area, which changed very slightly over the study period. The latest version of the fields' boundary is used. The surface area of the fields varies from less than 0.5 ha to more than 50 ha.

To apply our method over the whole study area, we average the backscatter into each field, so the speckle effect can be reduced. Instead of validating a 100-m resolution map, we validated the backscatter at each field scale to be consistent with the ground truth. Some of our fields are smaller than the 100 m by 100 m size, so it makes more sense to consider a field scale for which we have ground information.

Metrics Mapping

For each segmented field, we calculate these four metrics using the averaged backscatter signals for both VV and VH polarization. The maps of mean value, variance, correlation length, and fractal dimension for VV polarization over the study area are shown in Figure 11. Areas including forests, urban areas, and water bodies are masked out, so all of the productive fields are left to be illustrated. In the mean value map, areas with low values relate to non-irrigated fields, which can be easily distinguished from the irrigated fields. The same case exists in the variance map, where areas with high values relate to non-irrigated fields. In the correlation length map, the non-irrigated fields have relatively lower values, but it is becoming more difficult to separate them from the irrigated fields. Finally, the fractal dimension map shows a mixed situation for these two classes, which we therefore decided not to use.

Classification Map

Using the tree classification method based on the SVM with the mean value, variance, and correlation length for both VV and VH polarization as input, three classes, including non-irrigated fields, irrigated crops, and irrigated trees are classified after masking out forest, urban areas, and water bodies, as shown in Figure 12. The blue color represents the irrigated croplands, while the green color represents the irrigated trees. Yellow-colored areas are classified as non-irrigated fields.

Validation

Validation is done using SIGPAC information over all of the productive fields in our study area (26,434 fields in total). Each field is checked if it is correctly classified through being compared with the ground truth. The confusion matrix is listed in Table 1. The overall accuracy, which is the ratio of correctly classified fields and the total number of fields, is 81.08%, with the accuracy for non-irrigated crops slightly higher, with a value of 83.27%. The accuracy for irrigated crops is 77.53%, while for irrigated trees, it is 73.49%.

The most wrongly classified fields are irrigated crops into non-irrigated fields, with a percentage of 20.6%. The crops are mainly wheat. The different types of irrigation techniques may be a reason for the error, since fields using a new irrigated system are irrigated daily, while fields in the old irrigation area are irrigated about every two weeks, which makes the contribution of the irrigation events to the signal less important. The non-irrigated fields that were wrongly classified into irrigated crops were also high, with a percentage of 14.5%; most of them are located in the mountainous areas, where they may be influenced by topography. For irrigated trees, the most wrongly classified was into irrigated crops, which are the green-colored fields, with a percentage of 22.3%. The difference in the behaviors of backscatter signals for irrigated trees and very flouring crops may not be very evident.

The RF classification is also tested using half of the fields for learning and the other half for validating. The overall accuracy is around 82.2%, with the tree depth at three. The result shows the best accuracy for the irrigated crop, followed by irrigated trees and non-irrigated crops, which is the worst. The accuracy may vary from 80% to 83% when using different ground truth fields for learning, since RF is a very sensitive method. When using mean value alone as the input feature, the accuracy is about 80.9%, and using variance alone, the accuracy is about 66.4%, while when using correlation length alone, the accuracy is about 70.3%. Table 2 shows an example when using different combinations of metrics as input for RF classification. Combined the mean value with the variance, the accuracy increased by 1% above using the mean value only. After adding the correlation length, the accuracy increased by 0.3%.

Discussions

This study uses Sentinel-1 SAR multi-temporal data from June 2015 to September 2017 and focused on the statistics and metrics including the mean value, the variance, the correlation length, and the fractal dimension. In the classification, only the first three metrics are used as input. The SVM and RF are tested, and both methods showed a good accuracy.

The mean value shows a clear difference for irrigated and non-irrigated fields in VV polarization, as shown in Figure 11a. The non-irrigated fields have smaller mean values of the backscatter time series than irrigated fields, because of the backscatter signal intensity increasing with the soil dielectric constant, which is proportional to the water content in the soil until it is flooded.

The variance values of non-irrigated fields are relatively higher than those for irrigated fields, as shown in Figures 6,9 and 11b. This behavior is because non-irrigated areas can generally reach extreme values, with high levels due to precipitation events, and extremely low values in the summer period due to a long absence of precipitation. This is not the case for irrigated fields that receive water frequently, and they rarely reach those extreme low moisture levels. The variance values are somehow mixed, but we can still find the pattern, as the signals of non-irrigated areas tend to be more spread out. There is a class with dense blue points for irrigated fields with low variance values (less than 1 for both VV and VH) and high mean values (greater than -11 dB for VV and -18 dB for VH) in Figure 9, which consists of the irrigated trees, whose relatively high volumes yield higher reflectivity with small changes.

The function of the VV polarization increases, while for the VH polarization, it decreases as shown in Figure 9, because VV is more sensitive to soil moisture, and VH is more sensitive to vegetation. For VV polarization in a relatively dry region such as in our study area, small changes in soil moisture with a low variance value will have a relatively low average backscatter mean value, while the increase of soil moisture on some days will lead to a high variance value and a relatively high average mean value. However, for the irrigated trees, the vegetation volumes contribute more to the backscatter signal, which leads to higher mean values with low variances. The case for VH polarization, which is more sensitive to vegetation, is more complex. For irrigated fields, the trees always have a high backscatter mean value with small changes, and for crops, the increase of variance may due to the flourishing of the vegetation as a consequence of irrigation. This could also be due to the decrease of canopy caused by harvest events, which in the end leads to relatively stable backscatter mean values. For non-irrigated fields, the vegetation volume will not be able to reach to a very dense situation because of water shortage, so the increase of variance can only be the consequence of harvest events, which of course will lead to a low backscatter mean value.

The correlation length shows a relatively higher value for irrigated fields than non-irrigated fields, as shown in Figures 7 and11c. As the irrigation event happens in the summertime, the backscatter time series show a more seasonal pattern, which brings a higher correlation length value. Since the irrigation status is quite different for different croplands (different irrigation method, different irrigation period, different amount of applied water, etc.), the characteristics of its self-correlation may vary substantially, as shown in the x-axis of Figure 10.

The fractal dimension didn't show a clear difference between irrigated fields and non-irrigated fields (Figures 8 and11d). The temporal backscatter signal is influenced by many factors such as the irrigation frequency, precipitation, the vegetation type, etc.; the self-similarity of the time series is hardly observed. We discarded it as input in the end.

For the classification step, we tested both the SVM and RF method, and both show good accuracy. In SVM classification, we use the linear SVM functions to separate irrigated fields and non-irrigated fields. Then, the irrigated fields are further separated into irrigated crops and irrigated trees. RF classification used half of the fields for learning and the other half for testing; the fields are also separated into irrigated crops, irrigated trees, and non-irrigated fields. RF can give a better accuracy with a larger number of ground truth fields for learning and a bigger tree depth. However, the accuracy of RF classification with the tree depth at three can change from 80% to 83% when using different fields for learning. Through increasing the input metrics, the RF method shows an increase of accuracy. The SVM method also gives a good accuracy, and is more robust. Besides, SVM does not need as many ground truth fields as RF, but of course, more ground truth fields may bring a better accuracy.

Conclusions

In this paper, a methodology for irrigation mapping using the Sentinel-1 SAR data time series is introduced. The backscatter mean value, the signal variance, and the correlation length are derived from the backscatter signal time series and are used for classification. The classification result of irrigated and non-irrigated fields is compared with the ground truth from SIGPAC database over 26,434 fields in the whole study area. The result shows a good overall accuracy, which is 81.08% using the SVM. The mean value of the backscatter time series is the key to separating irrigated and non-irrigated fields, while the correlation length and variance are used for the separation of irrigated trees and irrigated crops.

The results of the SVM show good accuracy for irrigated croplands and non-irrigated fields, with soil moisture change more dominant to multi-temporal radar signals. Classification results for irrigated trees are slightly poor, and the effect of vegetation cover over backscatter in the case of trees is more obvious than the effect of soil moisture change.

The RF classification gives a similar accuracy as SVM, but it depends a lot on the number and location of the ground truth fields for machine learning, and is less robust compared with SVM.

In areas where there is no ground field information available, instead of using field boundary information, field segmentation was performed using Sentinel-2 NDVI data [57][58][59], but only for areas without frequent cloud cover. Another alternative method is using 100-m resolution cells directly to calculate the metrics only from Sentinel-1 SAR data, which is more consistent with our original intention to use the methodology under any weather conditions. This approach can be used in any agricultural areas when SAR data is available, and is adapted to regions with a limited use of S2 because of climate conditions. For areas that are more humid than the Catalunya study area during irrigation season, soil moisture contribution will be more limited, which may make the method less robust. Additionally, this method does not need to develop operational algorithms to estimate soil moisture before application; it is based directly on radar signal analysis. It is unrestricted by weather conditions and the location of the fields. The results demonstrated the potential of using Sentinel-1 data for irrigation mapping at the field scale.

94

CHAPTER 5. IRRIGATION MAPPING

Summary and Conclusions

In this chapter, an approach is developed for irrigation based on Sentinel-1 radar data considering both VV and VH polarization. Three statistics and metrics are found to be operational for classification of irrigated crops, irrigated trees and non-irrigated fields. The validations with ground measurements show good performances of both SVM and RF classifications. However, RF depends too greatly on the number and location of the ground truth fields for machine learning and is less robust compared with SVM.

The results of classification demonstrate the potential of Sentinel-1 data for irrigation mapping at the field scale. This approach can be used in any agricultural areas when SAR data is available, and it is adapted to regions where the use of optical data is limited because of climate conditions. Additionally, this method does not need to develop operational algorithms to estimate soil moisture before application; it is based directly on radar signal analysis.

In the present study, individual fields were identified thanks to the SIGPAC system. In areas where there was no ground field information available, an alternative method was to calculate the metrics from Sentinel-1 SAR data at the product spatial resolution. With the irrigation status monitored, the inclusion of it within an LSM needs to be investigated in a further study.

CHAPTER 6. WATER LEVEL RETRIEVAL

In this chapter, a short description of the methodology for water level retrieval is presented in Section 6.2, including the retrackers and waveform portion selection using DEM. The published article corresponding to this study can be found in Section 6.3, providing a more detailed literature review of inland water study and different retrackers, the methodologies, results, as well as detailed discussions and analysis. Finally, the summary and conclusions are presented in Section 6.4.

The author contributions can be found at the end of the published paper. As the first author of the article, Q. Gao performed the processing and analysis under the supervision of M.J. Escorihuela and based on the work of E. Makhoul (Makhoul et al., 2018).

Methodologies

Retracking is a process that is applied to the waveforms to get a resulting estimated geophysical parameters. In this study, three retrackers are modelled to pick out the tracking point, which is corresponding to the satellite range, and furthermore to the water level:

• Threshold retracker, which is a simple retracker based on an estimation of the epoch (the leading edge position) as a percentage of the maximum peak (Davis, 1997(Davis, , 1995)).

• Offset Centre of Gravity (OCOG) retracker, which is an empirical retracker that implements a combination of the centre of gravity (COG) and a conventional offset to estimate the related epoch (the leading edge position). Three main parameters are estimated for the OCOG retracker: the COG, the window size (W) and the amplitude (A), as shown in Figure 6.1. It is integrated into a two-step fitting procedure to meet the diversity of the backscattered echo waveforms over inland waters. In a first step, there is the fitting of epoch, amplitude, and SWH whilst keeping the roughness parameter fixed, then a second run is the fitting of the epoch, amplitude and roughness whilst keeping the SWH fixed. The estimate of the surface height was extracted from the fit with the best correlation coefficient between the measured waveform and the modelled waveform.

Retracking waveforms over inland water bodies are complex since the heterogeneity of the land covers potentially contaminates the waveform from the inland water bodies of interest. An example of the variation of waveforms in continental waters is shown in Figure 6.1. The presence of multiple peaks represents a challenge for the retracker when trying to obtain an accurate estimation of the continental water level. Taking into account the complexity of the multi-peaks in the backscattered waveforms, a methodology of the waveform portion selection using DEM information is proposed to select the portion of the waveform from nadir before the retrackers being applied, which can be seen as a pre-processing in the retracking chain. An illustration of the waveform portion selection process is shown in Figure 6.2. With the DEM height (SRTM DEM at a resolution of 30 m was used in this study) and the window delay as input, the method selects the waveform portion from nadir return based on a peak-valley approach. It selects the peak whose window delay is close to the nadir return (calculated from DEM), and selects the samples defined by the valleys surrounding this peak in addition to some guard samples to the left and right. The water level for each footprint is retrieved based on the selected portion of the waveform with the application of the retracker. The time series of the water levels are calculated using a water mask polygon. The results show good agreement with the in situ measurements with ubRMSE over the Ribarroja Reservoir at about 16 cm, and that over the Ebro Reservoir at about 28 cm. A detailed methodology description can be found in the published article integrated with Section 6.3.

Published paper

Introduction

Inland water systems constitute crucial resources of fresh water necessary for human survival [1]. Water within rivers and reservoirs represent the primary supply of drinking water, agricultural irrigation and industrial water usage globally [2], and also used to produce renewable hydrological energy. In addition, floods that periodically occur in every region of the world represent threats to crops, settlements, infrastructure, and most importantly life. Therefore, it is important to monitor the water levels of inland water bodies to provide early-warning alerts for water shortages or flood predictions. Unfortunately, in situ gauging stations are not always available in many parts of the world, or not publicly available, being maintained by local authorities [3][4][5]. An alternative to this is the use of satellite radar altimeters. Satellite radar altimeters are essential tools for monitoring the oceans, a task they have performed for over 25 years [6]. Satellite radar altimetry has also proven to be valuable tools for monitoring the water levels within inland water systems, including lakes and rivers [3,7,8]. Satellite radar altimetry can help overcome the lack of data in many parts of the world, and contribute to monitoring the water levels both within inland water systems where no in situ data are available, and also within transboundary river basins.

A list of satellite radar altimetry missions and their resolutions is shown in Table 1. Until CryoSat-2 was launched, the satellite radar altimetry was pulse limited to low-resolution mode (LRM). CryoSat-2 is the first altimetry mission with the synthetic aperture radar (SAR) mode available, followed by Sentinel-3. Unlike classical pulse-limited altimetry, SAR altimetry exploits coherent processing of groups of transmitted pulses to make the most efficient use of the power reflected from the surface [9], therefore significantly improving along-track resolution [9][10][11], as shown in Table 1. Sentinel-3, the newest in-orbit satellite, has a temporal resolution of approximately 27 days and an inclination angle of approximately 98.6 degrees, meaning it can cover almost all Earth latitudes except the poles. Operating in high-resolution SAR mode [10], the ground track separation of Sentinel-3 at the equator is approximately 104 km. In contrast, the CryoSat-2 ground track separation at the equator reaches 7.7 km at the expense of its temporal resolution, at 369 days. Sentinel-3 exhibits the best along-track resolution at 300 m, close to the resolution of CryoSat-2 in SAR and interferometric SAR (SARIn) modes. However, CryoSat-2 has a focus on the cryosphere, and as a result covers most continental surfaces in LRM. Therefore, Sentinel-3 constitutes the first altimetry mission that covers the globe completely in SAR mode, and is thus the most ideal tool for inland water level monitoring thanks to its good global coverage and sufficient temporal and spatial resolutions. 1 Geodetic Satellite. 2 European remote sensing satellite. 3 Topography Experiment/Poseidon. 4 GEOSAT Follow-On satellite. 5 Environmental Satellite. 6 Earth Explorer Opportunity Mission. 7 Haiyang-2/Ocean-2. 8 Satellite with Argos (Data Collection System) and ALtiKa (Altimeter in Ka-band). 9 Surface Water Ocean Topography Mission. 10 Jason Continuity of Service Mission.

The majority of previous studies [12][13][14][15][16][17][18][19][20][21][22][23][24][25][26][27] focus on relatively large water bodies with a scale of several kilometers. One of the first studies employing satellite altimetry for inland water level extraction was performed by Koblinsky et al. [12], who processed Geodetic Satellite (Geosat) waveforms to estimate the water levels at four sites in the Amazon basin. A root-mean-square 100 CHAPTER 6. WATER LEVEL RETRIEVAL error (RMSE) of 70 cm was estimated between the satellite and in situ measurements. Since then, with new generations of satellites, the spatial resolution and the accuracy of the orbit determination have both improved, thereby increasing the accuracies of the estimation results. Accordingly, several studies have since focused on lakes and large rivers with scales, ranging from 1000 km 2 to 80,000 km 2 , and demonstrated good results [12][13][14][15][16][17][18][19][20][21][22][23][24][25][26][27]. For example, Nielsen et al. [17] used Environmental Satellite (Envisat) and CryoSat-2 SAR data to retrieve the water levels over Vaenern (Sweden, 5650 km 2 ) and Lake Okeechobee (Florida, 1900 km 2 ), and obtained good accuracies with RMSEs, varying from 4 cm to 9 cm between the satellite-retrieved water levels and in situ measurements. Currently, several online data hubs, such as the Database for Hydrological Time Series of Inland Waters (DAHITI) (https://dahiti.dgfi.tum.de/en/) [28], Hydroweb (http://hydroweb.theia-land. fr/) [29], Global Reservoir and Lake Monitoring (https://ipad.fas.usda.gov/cropexplorer/global_ reservoir/) [30], and Thematic Exploitation Platform (TEP) for Hydrology (https://hydrology-tep. eo.esa.int/) [31], also provide time series of water levels over large lakes and rivers. The web platform named SAR Versatile Altimetric Toolkit for Ocean Research & Exploitation (SARvatore) (http://gpod.eo.esa.int/services/CRYOSAT_SAR/) can also provide Level-2 water level product of Sentinel-3 by exploiting 80 Hz data on line and on demand using the G-POD (Grid-Processing On Demand) service [32].

Numerous studies [33][34][35][36][37][38][39][40][41][42][43] have also focused on middle-sized water bodies. For example, Birkinshaw et al. [33] measured the water level of the Mekong River, the width of which varies from 400 m to 1700 m; their results show RMSEs of approximately 44-65 cm for Envisat and 46-76 cm for European Remote-Sensing Satellite-2 (ERS-2). Da Silva et al. [34] studied water level time series over the Amazon River basin, which exhibits widths ranging from several kilometers to less than a hundred meters, using both Envisat and ERS-2; they showed a RMSE varying from 12 cm in the best cases, 40 cm in most cases, to several meters in the worst cases. Furthermore, Michailovsky et al. [35] employed Envisat to monitor the Zambezi River basin and reported river widths reaching 80 m with RMSEs of 32-72 cm relative to in situ measurements at different locations. In addition, Maillard et al. [36] measured the water levels of medium-sized rivers with widths between 100 and 1000 m over the São Francisco River, Brazil, with RMSEs lower than 60 cm, and better than 30 cm in some cases using Envisat and Satellite with ARgos and ALtiKa (SARAL) data. Most studies of middle-sized water bodies used mainly Envisat and ERS. In addition, previous results referred to LRM altimeters, which provide accuracies on the order of tens of centimeters. RMSEs over different water targets vary with size and morphology, with smaller size and complexity of the shape bring more uncertainties. With the launch of Sentinel-3, equipped with a SAR Altimeter (SRAL), there is a need to evaluate its performance. Because of its higher spatial along-track resolution, SRAL is naturally better suited to eliminate land contaminations found within LRM altimetry footprint [44]. Our study objective was to retrieve water levels in challenging environments exploiting SAR altimeter data capabilities. We focused on middle-sized and small-sized water bodies whose width varies from hundreds of meters (close to along-track resolution) to a few kilometers, and/or surrounded by rapid changing topography.

As the backscattered waveform depends on the surface reflecting the signal, altimeter data are usually applied by different retrackers adapted to different surfaces to better locate the height of the reflective surface. A lot of research has been carried out regarding this, and is still ongoing in the retracking modeling framework, specifically considering the evolution from the LRM to SAR altimetric operation [45]. There are two types of retrackers, based on either empirical observations and practical experience, or theoretical knowledge of microwave scattering at nadir, namely empirical retracker and physical-based retracker, respectively. Most notable empirical retrackers include the threshold [46], beta-5 [47], and the Offset Center of Gravity (OCOG) [46,48], with the OCOG being the core of ice-1 retracker [49]. The physical-based retrackers include the LRM ocean retrackers based on the well-known Brown ocean waveform model [50] (e.g., ice-2 retracker [51,52]), and the altimetry SAR mode ocean retrackers based on the original work of Ray et al. [53] (e.g., the SAMOSA retracker [53,54]).

PUBLISHED PAPER

Whilst well-developed retrackers have been defined for operation over the oceans and ice surfaces, called ocean retrackers and ice retrackers, respectively, no retracker has been adapted for inland water bodies as of yet. Jarihani et al. [15] compared the results from different satellite missions for inland water bodies, where most previous studies used ice retrackers [15][16][17][18]24,34,36], and some others included ocean retrackers [15,16,24]. Nevertheless, different retrackers need to be developed and compared for inland water bodies. Due to the relatively recent advance towards SAR altimetry operation, started by the CryoSat-2 mission [55], and continued with the Sentinel-3 [56] mission, most of the studies over inland water are quite limited to the classical LRM operation. When trying to retrack waveforms over continental water bodies, the radar altimetric community has widely applied different retrackers including threshold, ice-1 and ice-2 retrackers. The robust OCOG retracker shows a consistent behavior over the different types of waveforms found in [57] with SARAL/Altika mission over different rivers and reservoirs in India. The SAMOSA retracker and its further adaptations have proven well fitting with CryoSat-2 data not only over the open-ocean and coastal zones [58], but also over ice sheets and inland waters where waveforms show specular characteristics [25,59].

Studies also show the possibilities of retrieving water levels over very small water bodies using techniques that allow precisely monitoring water surfaces. Investigation of individual echoes of Envisat LRM altimetry [60][61][62][63] has shown the possibility of retrieving water levels over water bodies with widths down to 50 m [63]. The SARvatore service provides 80 Hz data dedicated to inland water bodies. Another retracker, "fully focused SAR", has shown the potential of monitoring 40 m wide water surfaces by reducing the along-track resolution down to the theoretical limit, equaling half the antenna length [64]. The fully focused processing is a very promising processing technique that can provide very high resolution Level-1B products, and has the ability to open a new paradigm in SAR altimetry when it becomes fully operational.

The objective of our study was to develop an operational methodology for water level retrieval from altimetry SAR mode Level-1B waveforms. We focused on small water bodies where continuous clean waveforms rarely exist, and efforts were made to better retrack the contaminated waveforms to improve the accuracy of the water level. First, we modeled three different retrackers including the most widely used threshold retracker, the OCOG retracker, and our in-house SAR physical-based retracker [65] (isardSAT SAR retracker in [65] is an implementation based on the original physical-based model developed by Ray et. al [53], which is also the origin of the sometimes-called SAMOSA model). The latter was integrated in a two-step fitting procedure, exploiting both its normal operation over ocean-like waveforms, and the adaptation to peaky-like waveforms, expected from inland waters. For the second typology, the mean squared slope (MSS) describing the sea surface roughness [66] was used as a fitting parameter instead of the significant wave height (SWH), as proposed and initially evaluated in [59] over Tibetan lakes when using SAMOSA model. (In isardSAT retracker, a rough surface of Gaussian distributions is assumed to define the scattering model, such that the normalized radar cross section (NRCS) as a function of incidence angle is modeled as a Gaussian [67], where the MSS is a parameter defining the directivity of surface radiation pattern, i.e. small values of MSS lead to very narrow and directive radiation patterns related to very specular surfaces. Therefore, considering this parameter (MSS) as fitting term rather than a constant one would allow better fitting specular waveforms.) Additionally, to limit land contamination within the received waveforms, we included Digital Elevation Model (DEM) information to select the waveform portion from nadir, which was performed as a pre-processing in the whole retracking chain. The remainder of this paper is organized as follows. In Section 2, the studied area and the database are presented. In Section 3, the methodologies are introduced, including the selection of the waveform portion, different retrackers and waveform filtering approach. In Section 4, the performances of the different retrackers and the results of the water level time series are shown in comparison with the Sentinel-3 Level-2 ocean retracked results from the European Space Agency (ESA) and in situ measurements. In Section 5, the method and results are discussed. In Section 6, the conclusions are presented.

Study Area and Database

Study Area

Our study area encompassed the reservoirs and rivers in the Ebro River basin (Figure 1). The Ebro River, which has a length of approximately 910 km and a drainage basin with an area of approximately 85,362 km 2 , is one of the most important rivers on the Iberian Peninsula [68]. The river flow is irregular throughout the year, with low levels at the end of summer and high levels during the spring due to melt runoff in the Pyrenees, leading to a danger of flooding. The annual variation of the river flow is more than 800 m 3 s -1 at Ribarroja dam, which is the last dam controlling the water flow entering the Mediterranean Sea [69]. The Ebro River is of great importance for agriculture in summer, during which drought often occurs due to the continental Mediterranean climate. Nevertheless, the mean annual flow has decreased by approximately 29% during the 20th century due to many causes: the construction of dams, the increasing demands for irrigation, and evaporation, the latter being higher than the precipitation due to low rainfall amounts, high sunshine intensities and strong, dry winds, from reservoirs within the river basin [68].

from the European Space Agency (ESA) and in situ measurements. In Section 5, the method and results are discussed. In Section 6, the conclusions are presented.

Study Area and Database

Study Area

Our study area encompassed the reservoirs and rivers in the Ebro River basin (Figure 1). The Ebro River, which has a length of approximately 910 km and a drainage basin with an area of approximately 85,362 km 2 , is one of the most important rivers on the Iberian Peninsula [68]. The river flow is irregular throughout the year, with low levels at the end of summer and high levels during the spring due to melt runoff in the Pyrenees, leading to a danger of flooding. The annual variation of the river flow is more than 800 m 3 s -1 at Ribarroja dam, which is the last dam controlling the water flow entering the Mediterranean Sea [69]. The Ebro River is of great importance for agriculture in summer, during which drought often occurs due to the continental Mediterranean climate. Nevertheless, the mean annual flow has decreased by approximately 29% during the 20th century due to many causes: the construction of dams, the increasing demands for irrigation, and evaporation, the latter being higher than the precipitation due to low rainfall amounts, high sunshine intensities and strong, dry winds, from reservoirs within the river basin [68]. The Pyrenees mountain range lies to the north of the Ebro River basin. Consequently, it is challenging to retrieve the water surface heights of the studied area as the water bodies are relatively small and are greatly influenced by the mountainous terrain. The selected water bodies, for which we have gauging stations for validation, are shown in Figure 1. The coordinates, the widths of the water bodies covered with satellite tracks, the approximate distances between the gauging station and satellite tracks, and the average slope within 5 km along the satellite track are listed in Table 2. The water bodies that can be retracked properly, along with their DEMs, are shown in Figure 2 along with the water masks, Sentinel-3 tracks and the location of gauging stations. Only tracks overpassing the upper stream of the river (before the dam) were considered, ensuring the results could be validated against the in situ measurements. The Pyrenees mountain range lies to the north of the Ebro River basin. Consequently, it is challenging to retrieve the water surface heights of the studied area as the water bodies are relatively small and are greatly influenced by the mountainous terrain. The selected water bodies, for which we have gauging stations for validation, are shown in Figure 1. The coordinates, the widths of the water bodies covered with satellite tracks, the approximate distances between the gauging station and satellite tracks, and the average slope within 5 km along the satellite track are listed in Table 2. The water bodies that can be retracked properly, along with their DEMs, are shown in Figure 2 along with the water masks, Sentinel-3 tracks and the location of gauging stations. Only tracks overpassing the upper stream of the river (before the dam) were considered, ensuring the results could be validated against the in situ measurements.

(e) (f) (g) (h) (i) (j)

Data Base

Sentinel-3

Sentinel-3 is an ocean and land mission based on a constellation of two satellites (Sentinel-3A and Sentinel-3B). Sentinel-3A was launched on 16 February 2016, with data available beginning in June 2016. It was followed by Sentinel-3B, which was launched on 25 April 2018, with data available beginning in December 2018. In this study, two years of Sentinel-3A data were used from June 2016 to June 2018. The SRAL instrument is the main topographic sensor used to provide water level

Data Base

Sentinel-3

Sentinel-3 is an ocean and land mission based on a constellation of two satellites (Sentinel-3A and Sentinel-3B). Sentinel-3A was launched on 16 February 2016, with data available beginning in June 2016. It was followed by Sentinel-3B, which was launched on 25 April 2018, with data available beginning in December 2018. In this study, two years of Sentinel-3A data were used from June 2016 to June 2018. The SRAL instrument is the main topographic sensor used to provide water level measurements, and, hence, it was used in our research. The detailed parameters of the Sentinel-3 SRAL can be found in [70]. To acquire altimeter measurements, the Sentinel-3 SRAL transmits pulses at a Ku-band frequency, which is complemented by a C-band frequency to correct range delay errors due to the varying density of electrons in the ionosphere [70]. Sentinel-3 has two operational modes: SAR mode and LRM. As the SAR mode is available globally, we could retrieve inland water levels over any area tracked by Sentinel-3. The SRAL tracks the surface in two different tracking modes, namely, closed loop and open loop tracking. For closed loop tracking, the altimeter range window is autonomously positioned based on an on-board near real-time (NRT) analysis of the previous SRAL waveform; in contrast, for open loop tracking, the altimeter range window is positioned using a priori knowledge of the surface height stored in the instrument in a one-dimensional along-track DEM. The tracking modes for the studied water bodies over the Ebro River basin are listed in Table 2.

Three levels of processed altimeter data are available: Level-0, Level-1 and Level-2 products. In our study, Level-1 non-time critical (NTC) 20 Hz data were used for the water level retrieval by the three retrackers (i.e., the threshold, OCOG and two-step physical-based retrackers), and Sentinel-3 Level-2 ocean retracked data from ESA were used for comparison.

The main objectives of the Level-2 processing of SAR mode data are to provide elementary retracked altimeter estimates of the oceans, coastal zones, ice sheets and sea ice elevations [71]. Different retracking algorithms are more suited to specific surfaces, that is, for ocean retracking, OCOG retracking, ice sheet retracking, ice retracking and sea ice retracking. Unfortunately, ice-related retracker results are not available for inland water bodies; hence, we used Sentinel-3 Level-2 ocean retracker results for comparison. Ocean retrackers try to fit theoretically modeled multi-look Level-1B (L1B) waveforms to real L1B SAR waveforms, thereby providing estimates of the epoch, composite sigma, amplitude and the mispointing angle [71].

In Situ Data

In situ data for the Ebro River basin are available in the Automatic Hydrological Information System, in Spanish known as the "Sistema Automático de Información Hidrológica", or the SAIH Ebro data hub [72]. SAIH Ebro is an online system providing hourly and daily hydrological information, including river gauge data, reservoir levels, rainfall amounts and temperatures, over the Ebro River basin. Together with the Sentinel-3 passes, we studied eight water bodies in our research, as listed in Table 2. The in situ data were collected from June 2016 to June 2018 on an hourly basis.

Digital Elevation Model

A DEM was used as an ancillary dataset for the selection of the waveform portion. The Shuttle Radar Topography Mission (SRTM) is an international research effort that obtains DEM data on a near-global scale from 56 • S to 60 • N [73]. The SRTM provides global data at two resolutions: 1 arc-second (~30 m) and 3 arc-seconds (~90 m). In our study, 1 arc-second global elevation data, which offer a worldwide coverage of void-filled data at a resolution of approximately 30 m, were used. The DEMs of water bodies are shown in Figure 2. The vertical accuracy of SRTM DEM is about ±16 m (absolute) and ±6 m (relative) [73][74][75][76][77]. However, vertical accuracy of the data decreases with the increase in slope and elevation due to presence of large outliers and voids [75].

Methodology

Geophysical Corrections

The space-borne radar altimeter is an essential tool for monitoring the oceans, but it can also be used for inland water surfaces, including lakes and rivers. The principle of altimetry can be found in [78] with application of geophysical corrections. In this study, corrections for the wet troposphere, dry troposphere, ionosphere, solid earth tide, geocentric pole tide and ocean loading 106 CHAPTER 6. WATER LEVEL RETRIEVAL tide were considered, as discussed by Fernandes et al. [79], and the geoid correction was applied.

As the objective of the study was to directly exploit official Sentinel-3 L1B data and compare against Level-2 official products, the use of the geophysical corrections available in the Level-2 product were considered. The detailed corrections and their ranges are listed in Table 3 [80]. The Geoid model used is the Earth Gravitational Model 2008 (EGM2008) [81]. To improve the altimeter range accuracy, which is related to the water level measurement accuracy, the waveform needs to be retracked precisely to determine the accurate tracking point located on the leading edge [START_REF] Deng | A coastal retracking system for satellite radar altimeter waveforms: Application to ERS-2 around Australia[END_REF]. We tested three different retrackers, which are listed in Section 3.2, to find the tracking points precisely from the land-contaminated waveforms.

Retrackers

Threshold Retracker

The threshold retracker is a simple retracker based on an estimation of the epoch (the leading edge position) as a percentage of the maximum peak [46,[START_REF] Davis | Growth of the Greenland ice sheet: A performance assessment of altimeter retracking algorithms[END_REF]. In principle, it works well when the maximum peak originates from a nadir water body. The epoch is calculated as the first sample or range bin that is above a percentage of the waveform peak. To provide a higher precision on the estimated epoch, a linear interpolation between the range bins adjacent to the threshold crossing is used as suggested in [46]:

epoch = n th -1 + A peak • η th -y(n th -1) /(y(n th ) -y(n th -1)) (1)
where n th is the first range bin of the waveform whose power, y(n th ) is right above the threshold (A peak • η th ), A peak is the amplitude of the peak of the waveform, and η th is defined as a threshold percentage (50% in our case). The term A peak • η th -y(n th -1) /(y(n th ) -y(n th -1)) in Equation ( 1) corresponds to the linear proportionality coefficient.

Offset Center of Gravity (OCOG) Retracker

The OCOG retracker [46,48] was designed to calculate the center of gravity of the reflected waveform based on the power levels in the bins. It is an empirical retracker that implements a combination of the center of gravity (COG) and a conventional offset to estimate the related epoch (the leading edge position). Three main parameters are estimated for the OCOG retracker: the COG, the window size (W) and the amplitude (A).

Frappart's definition [24] was considered in our analysis using squares of the power samples:

COG = ∑ n 2 n=n 1 n•y 2 (n) ∑ n 2 n=n 1 y 2 (n) W = (∑ n 2 n=n 1 y 2 (n)) 2 ∑ n 2 n=n 1 y 4 (n) A = ∑ n 2 n=n 1 y 4 (n) ∑ n 2 n=n 1 y 2 (n) (2)

PUBLISHED PAPER

where y(n) is the n th power sample of the input waveform, and n 1 and n 2 are the initial and end range bins or positions, respectively, in the waveform used for the OCOG parameter estimation.

The estimated epoch can then be implemented using the conventional OCOG definition of the epoch:

epoch = COG -W/2 (3) 3.2.3. Two-Step Physical-Based Retracker
The in-house isardSAT SAR ocean retracker in [65] is integrated in a two-step fitting procedure to operate over inland waters, providing robust surface height estimation with a minimal modification to the SAR ocean retracker model. The waveforms reflected from an inland water body might exhibit ocean-like shapes (winds blowing over lakes) or, most commonly, a peakier shape (smoother surfaces) for small water bodies. To adapt to the latter ones, the MSS, describing the sea surface roughness [66], was used instead of the SWH as a fitting parameter; this approach is initially considered in [59], when exploiting the SAMOSA model over Tibetan plateau lakes.

The SAR ocean retracker in [65] considers a scattering model of rough surface with Gaussian statistics, and so the normalized backscattering coefficient (surface radiation pattern) as function of the incidence angle θ inc can be modeled as [67]:

σ 0 (θ inc ) = σ 0 (0) • exp - (tan θ inc ) 2 MSS ( 4
)
where the MSS represents the variance of the surface slopes. From the relationship in Equation ( 4), it can determined that, when the variance of the surface slopes (MSS) is large, a broad surface radiation pattern is obtained (rough surfaces), while, when MSS tends to small values, the pattern is more directed as expected for specular returns (smooth surfaces). Therefore, the MSS can be accordingly adjusted as a fitting parameter to properly fit specular returns. With the objective of exploiting the SAR ocean retracker [65] with minimal changes over inland waters, the same retracker model can be used, but the SWH is fixed (to a value close to zero), while the MSS or surface roughness is fitted. A similar approach with the SAMOSA model is exploited in [START_REF] Jain | Evaluation of SAMOSA3 adapted retracker using Cryosat-2 SAR altimetry data over the Arctic ocean[END_REF] to fit waveform returns from leads.

To ensure that the different typology of waveforms (ocean-and specular-like) can be fitted, a two-step retracking is implemented. In a first run, the normal operation of the physical-based retracker in [65] was considered, i.e., the SWH was one of the fitting parameters (keeping MSS constant). Then, a second run of the fitting was implemented but this time the SWH was fixed (to a very small value around 1e-5) and the MSS was fitted so that the radiation pattern of the rough surface was adjusted. Then, the estimate of the surface height was extracted from the fit with the best correlation coefficient between the measured waveform and the modeled waveform. In this way, it was ensured that peaky-like waveforms could be properly fitted using the MSS as a fitting parameter, since the normal SAR ocean retracker would fail to adjust the model. An example of the operation of the two-step retracker is shown in Figure 3 over a peaky waveform.

The in-house isardSAT SAR ocean retracker integrated in the two-step procedure [65] is an implementation based on the SAR ocean model proposed by Chris Ray et al. [53]. It is a physically-based model of the altimetric power waveform over the ocean received by the altimetric sensor after delay-Doppler (or SAR) processing. Details regarding the definition and implementation of this retracker can be found in [65].

In SAR or delay-Doppler altimetry, different Doppler beams, also known as looks (l), are synthesized to focus on a given surface on the ground [65]. Thus, averaging of all these Doppler power beams would lead to the final multilook power waveform to be retracked. Therefore, for this physical-based retracker, each single-look power waveform was modeled as follows:

108

CHAPTER 6. WATER LEVEL RETRIEVAL

P (k,l) (P u , k 0 , SW H/MSS) = P u • B (k,l) (MSS) • (g l (SW H))[ f 0 (g l (SW H) • (k -k 0 )) + T (k,l) (MSS) • g l (SW H) • σ 2 z L 2 z • f 1 (g l (SW H) • (k -k 0 ))] (5)
where k refers to the range index vector within the received waveform window, l is the Doppler index or beam pointed to the specific surface, P u is the amplitude, k 0 is the epoch (related to the surface height), SWH is the significant wave height (referred also as H s ), and MSS is the mean square slope related to the surface roughness. In the first iteration of the two-step retracker, P u , k 0 , SW H are the fitting parameters, while in the second iteration the parameters being adjusted are P u , k 0 , MSS. The in-house isardSAT SAR ocean retracker integrated in the two-step procedure [65] is an implementation based on the SAR ocean model proposed by Chris Ray et al. [53]. It is a physically-based model of the altimetric power waveform over the ocean received by the altimetric sensor after delay-Doppler (or SAR) processing. Details regarding the definition and implementation of this retracker can be found in [65].

In SAR or delay-Doppler altimetry, different Doppler beams, also known as looks (𝑙), are synthesized to focus on a given surface on the ground [65]. Thus, averaging of all these Doppler power beams would lead to the final multilook power waveform to be retracked. Therefore, for this physical-based retracker, each single-look power waveform was modeled as follows:

where 𝑘 refers to the range index vector within the received waveform window, 𝑙 is the Doppler index or beam pointed to the specific surface, 𝑃 is the amplitude, 𝑘 is the epoch (related to the surface height), SWH is the significant wave height (referred also as 𝐻 ), and MSS is the mean square slope related to the surface roughness. In the first iteration of the two-step retracker, 𝑃 , 𝑘 , 𝑆𝑊𝐻 are the fitting parameters, while in the second iteration the parameters being adjusted are 𝑃 , 𝑘 , 𝑀𝑆𝑆.

The terms 𝐵 , and 𝑇 , encompass the information related to the antenna pattern, antenna mis-pointing, and surface radiation patterns. They are obtained, respectively, as the constant and first order terms of a Taylor approximation of the antenna and surface radiation pattern's product. The closed expression for these two terms can be found in [65]. A Gaussian-like antenna and surface radiation (as indicated in Equation ( 5)) patterns were assumed.

In Equation ( 5), 𝑔 refers to the Doppler-dependent (𝑙) dilation term, and takes into account The terms B k,l and T k,l encompass the information related to the antenna pattern, antenna mis-pointing, and surface radiation patterns. They are obtained, respectively, as the constant and first order terms of a Taylor approximation of the antenna and surface radiation pattern's product. The closed expression for these two terms can be found in [65]. A Gaussian-like antenna and surface radiation (as indicated in Equation ( 5)) patterns were assumed.

In Equation ( 5), g l refers to the Doppler-dependent (l) dilation term, and takes into account the instrument configuration and the significant wave height:

g l = 1 σ 2 ac + 2 • σ al • L 2 x L 2 y • l 2 + σ 2 z /L 2 z ( 6
)
where σ ac and σ al refer to the widths of the Gaussian functions that approximate the point target response in the across-track (or range) and along-track (azimuth) dimensions, respectively, and σ z corresponds to the standard deviation of the sea surface height probability density function (PDF) and is related to the SWH as σ z = H s 4 . The definition of the additional terms in Equation ( 6) can be found in Table 4, where the main parameters of the SAR ocean model are summarized.

Finally, f 0 and f 1 represent the range-dependent functions modulated by the Doppler-dependent dilation term and can be obtained from the general expression (for n = 0 and 1): Standard deviation of the height PDF 1 For a chirp pulse the bandwidth can be expressed in terms of the chirp rate (k r ) and pulse duration (T p ) as follows: BW = K r • T p .

f n (ψ) = ∞ 0 v 2 -Ψ n • e -(v 2 -Ψ) 2 2 • dv (7)
Hs 4 α R Orbital factor α R = R+H orb R L x Along-track resolution c0•H orb •PRF 2•vsat• fc•Np L y Across-track resolution c0•H orb αR•BW L z Vertical resolution c0 2•BW

Waveform Exclusion

For very small water bodies, the waveform may be highly contaminated by the land. Thus, we added additional conditions to exclude the waveforms with a poor quality.

First, if the epoch of the reference SRTM DEM is outside the waveform window, which has a measuring range of approximately 60 m, it is likely that the on-board tracker window was unable to focus on the reflecting surface. Second, within the vicinity of land, reflections from both water and land will contribute to the received waveform, resulting in more than one peak. If the number of outstanding peaks is larger than a threshold, the waveform may be contaminated substantially by land or contains no water information at all, making it not suitable for water level retrieval. In our case, the threshold was set to 5 based on the analysis of waveforms. The number of outstanding peaks is normally less five, even with land contaminations. Finally, sigma_0, which is the backscatter coefficient for water from nadir, should be much larger than the backscatter coefficient for land, where the value is normally less than 50 dB based on the analysis in our study area. With these conditions (summarized below), we could filter out highly contaminated waveforms that were unsuitable for detecting the levels of water bodies:

•

The epoch of the reference SRTM DEM must be within the waveform window

• Number of outstanding peaks < 5 • Sigma_0 (backscatter coefficient) > 50 dB
After waveform exclusion, only the remaining waveforms (including both slightly contaminated and clean waveforms) were considered in our study to maximize the usage of the available measurements in the water body of interest. For a small water body, the quality of the waveform largely depends on the satellite tracks which drift daily, and the surface area changes with seasons, resulting in no clean waveform from time to time.

110

CHAPTER 6. WATER LEVEL RETRIEVAL

Selection of the Waveform Portion

Over inland water bodies, the on-board tracker may not properly locate the retrieval window to acquire the signal from the nadir corresponding to the water body. In addition, specific across-track surface contamination will be present in the waveforms, making it difficult to properly track the desired portion of the waveform. In our study, we used DEM information (SRTM DEM at a resolution of 30 m) to locate the waveform portion that comes from nadir within the receiving window, and then retracked the waveform using the selected portion only. This pre-processing approach, which prepares the portion of waveforms, can be applied and integrated with any kind of retracker.

The waveform portion selection method isolates the nadir return within the receiving window as follows:

•

From the geo-located surface of interest within the water body and beneath the satellite track, the associated height was interpolated from the DEM information and referred to the geodetic ellipsoid, H DEM .

• This height was subtracted from the satellite height (H orb ) at the geo-located surface to obtain a rough estimation of the range (or equivalently window delay) from where the nadir returns were expected.

• This range was linked to a specific bin within the received waveform, by comparing it with the vector of ranges associated to each bin in the receiving window, which could be obtained from the measured range by the radar and the range sampling.

•

The peak location closest to the previous range bin position was taken, and the portion of the waveform was selected around this peak considering the valley positions to the left and right of the peak plus some guard samples on top:

A built-in Matlab function (findpeaks) was used to compute the prominent or outstanding peaks within the waveform. Prominent peaks are those peaks that drop more than a given threshold value on either side of the peak before the signal attains a higher value. The associated valley locations can be extracted using this built-in function, but taking as input the maximum of the whole waveform minus the waveform itself.

The water level for each footprint was retrieved based on the selected portion of the waveform applying the three retrackers. The time series of the water levels were calculated using a strict water mask polygon as shown in Figure 2, within which the water levels of the altimeter footprints were considered, selected and averaged for each track. The water mask polygon was created manually according to the shoreline displayed in the Google Earth image.

The workflow of the water level retrieval using the Sentinel-3 Level-1 product included four main steps: waveform exclusion, selecting the waveform portion using a DEM determining the retracking range using different retrackers, geophysical corrections and water level averaging, as shown in Figure 4. Geophysical corrections were applied to the retracked range using the Sentinel-3 Level-2 product, and the accuracy was assessed for each water body using in situ measurements from the gauging stations in the SAIH system.

Results

Selecting the Waveform Portion Using a DEM

When a waveform is influenced by land-based contamination, the waveform will contain small peaks due to off-nadir land or land coverage that follow the large peak, which corresponds to the water reflection within the altimeter footprint. In this case, as shown in Figure 5 (two outstanding peaks in Figure 5 waveform), the portion of the selected waveform provided a more precise tracking point for the two-step physical-based retracker and OCOG retracker, especially for the OCOG retracker, whose retracking point was the front edge of the COG window. The footprint of the signal was very close to the shore, thereby introducing a large peak at the trailing edge. The DEM information was used to select the portion of the waveform closer to nadir. The closest peak and valley after nadir were selected for retracking. Figure 5a,c shows the location of nadir in a black vertical line. The portion being selected is under the SWH fitting or MSS fitting line in Figure 5a with the two-step physical-based retracker, and is marked out in red crosses as the L1B-Waveform portion in Figure 5c with OCOG retracker. For the two-step physical-based retracker, the fitting using the selected waveform portion (Figure 5a) resulted in a 99% accuracy, while the fitting using the whole waveform (Figure 5b) resulted in an accuracy of 83%, making the location of the tracking point on the leading edge slightly different. In this study, the objective of the two-step physical-based retracker was to ensure proper retrieval of the surface height over inland waters, when exploiting the SAR ocean retracker in [65] with minimal changes on its operation; the

Results

Selecting the Waveform Portion Using a DEM

When a waveform is influenced by land-based contamination, the waveform will contain small peaks due to off-nadir land or land coverage that follow the large peak, which corresponds to the water reflection within the altimeter footprint. In this case, as shown in Figure 5 (two outstanding peaks in Figure 5 waveform), the portion of the selected waveform provided a more precise tracking point for the two-step physical-based retracker and OCOG retracker, especially for the OCOG retracker, whose retracking point was the front edge of the COG window. The footprint of the signal was very close to the shore, thereby introducing a large peak at the trailing edge. The DEM information was used to select the portion of the waveform closer to nadir. The closest peak and valley after nadir were selected for retracking. Figure 5a,c shows the location of nadir in a black vertical line. The portion being selected is under the SWH fitting or MSS fitting line in Figure 5a with the two-step physical-based retracker, and is marked out in red crosses as the L1B-Waveform portion in Figure 5c with OCOG retracker. For the two-step physical-based retracker, the fitting using the selected waveform portion (Figure 5a) resulted in a 99% accuracy, while the fitting using the whole waveform (Figure 5b) resulted in an accuracy of 83%, making the location of the tracking point on the leading edge slightly different. In this study, the objective of the two-step physical-based retracker was to ensure proper retrieval of the surface height over inland waters, when exploiting the SAR ocean retracker in [65] with minimal changes on its operation; the significance and interpretation of other geophysical parameters that can be inferred (SWH and MSS) needs to be further investigated and it is out of the scope of the current study.

For the OCOG retracker, the COG when using the selected waveform portion was located inside the first peak. However, when using the whole waveform, it was located in between peaks, resulting in different tracking points that were the front edge of the COG window. For the threshold retracker, the tracking point did not change in this case. The threshold retracker depended only on the largest peak, which, as shown in Figure 5, was the same with or without waveform portion selection. However, in the case of the largest peak not being from the nadir water body but from the surroundings or other objectives, the threshold retracker with and without waveform portion selection led to different results.

Remote Sens. 2019, 11, x FOR PEER REVIEW 15 of 26 significance and interpretation of other geophysical parameters that can be inferred (SWH and MSS) needs to be further investigated and it is out of the scope of the current study.

For the OCOG retracker, the COG when using the selected waveform portion was located inside the first peak. However, when using the whole waveform, it was located in between peaks, resulting in different tracking points that were the front edge of the COG window. For the threshold retracker, the tracking point did not change in this case. The threshold retracker depended only on the largest peak, which, as shown in Figure 5, was the same with or without waveform portion selection. However, in the case of the largest peak not being from the nadir water body but from the surroundings or other objectives, the threshold retracker with and without waveform portion selection led to different results. The differences between the performance using the waveform portion selection method, and the method using the whole waveform, are shown in the time series results over the Ebro Reservoir in Figure 6. The Level-2 data were obtained from the Sentinel-3 Level-2 product using the ocean retracker, which does not include any waveform portion selection. The waveform portion selection method eliminated most outliers compared with the method of using the whole waveform. The differences between the performance using the waveform portion selection method, and the method using the whole waveform, are shown in the time series results over the Ebro Reservoir in Figure 6. The Level-2 data were obtained from the Sentinel-3 Level-2 product using the ocean retracker, which does not include any waveform portion selection. The waveform portion selection method eliminated most outliers compared with the method of using the whole waveform.

Time Series Validation

The time series of the water levels were calculated using a strict water mask polygon. The water levels of the altimeter footprints within the mask were considered, selected and averaged for each date. Figures 6a and7 show the water levels for different water bodies derived using the three retrackers: the two-step physical-based retracker, OCOG retracker and threshold retracker, combined with the waveform portion selection method. The retracked results were compared with the in situ measurements and the Level-2 results using the ocean retracker.

Time Series Validation

The time series of the water levels were calculated using a strict water mask polygon. The water levels of the altimeter footprints within the mask were considered, selected and averaged for each date. Figures 6a and7 show the water levels for different water bodies derived using the three retrackers: the two-step physical-based retracker, OCOG retracker and threshold retracker, combined with the waveform portion selection method. The retracked results were compared with the in situ measurements and the Level-2 results using the ocean retracker.

Among the eight water bodies monitored in this study, the water levels over the Mequinenza Reservoir, Cavallers Reservoir and San Salvador Reservoir, all of which were tracked in an open loop tracking mode (on-board DEM dependent), could not be retracked. The difference between the on-board tracking heights derived from the received waveforms and the heights from the SRTM DEM with reference to the geodetic ellipsoid was almost 50 m over the Mequinenza Reservoir, more than 1000 m over the Cavallers Reservoir, and approximately 80 m over the San Salvador Reservoir. These findings indicate that the SRAL sensor lost its track and was unable to acquire the waveforms reflected from these water body surfaces. Among the eight water bodies monitored in this study, the water levels over the Mequinenza Reservoir, Cavallers Reservoir and San Salvador Reservoir, all of which were tracked in an open loop tracking mode (on-board DEM dependent), could not be retracked. The difference between the on-board tracking heights derived from the received waveforms and the heights from the SRTM DEM with reference to the geodetic ellipsoid was almost 50 m over the Mequinenza Reservoir, more than 1000 m over the Cavallers Reservoir, and approximately 80 m over the San Salvador Reservoir. These findings indicate that the SRAL sensor lost its track and was unable to acquire the waveforms reflected from these water body surfaces.

The RMSE and unbiased RMSE (ubRMSE) (Equation ( 8)) were calculated, as listed in Table 5. The unbiased RMSE is more reliable as it removes the bias.

𝑢𝑏𝑅𝑀𝑆𝐸 = 𝐸{[(𝑊𝐿 -𝐸[𝑊𝐿]) -(𝑖𝑛𝑠𝑖𝑡𝑢 -𝐸[𝑖𝑛𝑠𝑖𝑡𝑢])] }

where E[•] is the expectation operator and WL is the water level derived from satellite data. The RMSEs for different water bodies varied from 17 cm to more than 1 m. A comparison of the performances for the different retrackers over different water bodies is shown in The RMSE and unbiased RMSE (ubRMSE) (Equation ( 8)) were calculated, as listed in Table 5. The unbiased RMSE is more reliable as it removes the bias.

ubRMSE = E [(WL -E[W L]) -(insitu -E[insitu])] 2 (8)
where E[•] is the expectation operator and WL is the water level derived from satellite data. The RMSEs for different water bodies varied from 17 cm to more than 1 m. A comparison of the performances for the different retrackers over different water bodies is shown in

Discussion

In this study, the size of the target water bodies vary from 130 m to 4.5 km and change with seasons. The geometries of the water bodies are complex with irregular lake shores and various sizes, 116 CHAPTER 6. WATER LEVEL RETRIEVAL making the along-track width (the width as seen by the satellite flying along the water body) very different from one satellite overpass to another due to the drifting satellite ground tracks. As a result, continuous clean waveforms rarely exist in our study area. The shore of the Ebro Reservoir is irregular, and lands exist within water, resulting in no clean waveforms available over many tracks. Figure 6b shows the water level time series derived using the whole waveform, while Figure 6a shows the water level time series derived using the waveform portion selection. Retracking of the waveform portion accordingly reduced land contamination effects on the waveform. Figure 6a shows higher consistency with in situ time series, and thus shows that most of the outliers were caused by the land contamination. The land contaminated waveforms can be originated by surface area changes or satellite track drifting, making the altimeter footprint closer to the side bank, introducing outliers in the time series.

Over the water bodies where Sentinel-3 is operated in closed loop, the water levels of all water bodies could be tracked. Over the water bodies where Sentinel-3 is operated in open loop, three water bodies were missed and only one was tracked accurately. The results over the Ebro Reservoir, Sotonera and Ribarroja Reservoirs are satisfactory. The Ebro Reservoir is located in the north near the Pyrenees mountain range with elevations exceeding 820 m, and an average slope within 5 km along track of about 4%. The RMSE was 30 cm compared with the in situ measurements. The Sotonera Reservoir is the biggest reservoir studied, and is located in a relatively flat area with a slope of about 3% and an elevation of about 410 m. The result is good as there is no influence from the mountainous terrain, and it is big enough to be well monitored by satellite. The Ribarroja Reservoir is only 400 m wide and with a slope of more than 20%, but an elevation of only around 70 m. Its small width results in only one altimeter footprint per track available within the water body. The ubRMSE of Ribarroja Reservoir shows very good performance with a value of about 16 cm. However, the results for a few tracks were excluded due to poor waveform quality, which may be influenced by the rapid changes of surrounding terrain. The RMSEs over the Itoiz Reservoir and Irabia Reservoir exceeds 1 m. One reason for this could be the rapid change in the elevations along the satellite track as both water bodies are located in the Pyrenees, making it easy for a satellite to lose its track. The width of the Itoiz Reservoir with satellite data varies from 400 m to more than 2 km depending on the different satellite tracks, which drift daily. In addition, the Irabia Reservoir is a very small water body with a width of approximately 130 m, making the retrieval of its surface height more challenging. Besides the waveform quality, the size of the water body and the surrounding terrain, another reason for the differences with in-situ data from the gauging station is the distance between the gauging station and the Sentinel-3 tracks, varying from a few meters to several kilometers. Satellite tracks are always drifting (up to 1.5 km), making the footprints always changing in their locations for each water body. With the altimetry mission originally designed for ocean, the accuracy of geophysical corrections over inland water ranges from a few centimeters to tens of centimeters, depending on the size of the water body and wind conditions [79,[START_REF] Crétaux | Lakes Studies from Satellite Altimetry[END_REF][START_REF] Zhang | Application of retracked satellite altimetry for inland hydrologic studies[END_REF]. These effects will most definitely influence the final results, and need to be further considered quantitatively.

The results of the three retrackers, in combination with the waveform portion selection method, always give better accuracies than the results of the Sentinel-3 Level-2 ocean retracker (the results of both sets of retrackers are shown in Figure 6 and Table 5). Therefore, given the lack of an ice retracker over inland water bodies, retrieving the water levels by retracking the Level-1 waveforms and using the waveform portion selection method constitutes an alternative approach to obtain more accurate water levels.

The waveform portion selection method using DEM information greatly improved the results shown in both Figures 5 and6. As mentioned above, the vertical accuracy of SRTM DEM is approximately ±16 m (absolute) and ±6 m (relative). In the case of a 6 m error, the difference of choosing the right peak corresponding to nadir would be about 3 m. Therefore, if two large peaks were very close to each other, the portion of the waveform chosen might not be accurate. However, for most cases, the DEM error is small compared to the waveform tracking window, which is more than 60 m, and thus optimistic results can be obtained. Besides, waveform portion selection using DEM information improved the water level results over water bodies where only a few or no clean waveforms are available, and retained more results to have a better complete time series. In the case of the water level changing enormously and the waveform being heavily contaminated at the same time, the waveform portion selection approach may lose its capability of precise retracking.

The different retrackers obtained a good accuracy with in situ measurements. The ubRMSEs showed smaller values compared with RMSEs for all retrackers over different reservoirs, as listed in Table 5. The bias between the retrieved water levels and the in situ measurements existed in most cases, which may be related to the drift of the river course (the measurements are not co-located with the Sentinel-3 track). Other error sources include the waveform quality, which depends on the surrounding terrain. These include the geometry of the lake shore and the characteristics of the terrain, all of which contribute to the shape of the waveform. Nevertheless, Sentinel-3 has been proven to work for small water bodies.

We detected some issues related to the Sentinel-3 SRAL tracking mode. The open loop tracking mode, which depends on the on-board DEM to locate the tracking window, loses its track over the Cavallers Reservoir, which is surrounded by mountains, but also in Mequinenza and San Salvador Reservoirs, which are in relatively flat areas. We provided the correct water level height information to make Sentinel-3 able to capture the correct water level in the open loop tracking mode using the service altimeter open loop tracking command for hydrology monitoring (https://www.altimetry-hydro.eu/). We optimistically expect the water bodies that are under the Sentinel-3 coverage, but with loss of track, will be monitored in the near future. Another potential solution could be changing the tracking mode to a closed loop mode.

Conclusions

The objective of this study was to develop a methodology that can estimate water levels in small water bodies (similar size to the along-track resolution) and/or complex topography environments where many waveforms are expected to be contaminated by land. The study compared the performances of the threshold, OCOG and two-step physical-based retrackers, but no significant differences were found in the results of the three retrackers. However, the main improvements found show that selecting the portion of the waveform to be fitted significantly affected the retrievals.

The novel DEM-oriented waveform portion selection method could isolate the nadir peak from land-based contamination relatively well, and greatly improved the results with a reduction of RMSE by more than 0.5 m over Ebro Reservoir. Retracking the water levels from Level-1 waveforms using the combination of a retracker and the waveform portion selection method was more robust to land contamination than using Sentinel-3 Level-2 data directly. It also resulted in a better accuracy.

The water levels were compared with in situ measurements and the Sentinel-3 Level-2 ocean retracker from ESA. The results show good agreement with the in situ measurements. The ubRMSE over the Ribarroja Reservoir with a width of about 400 m could reach 16 cm, while that over the Ebro Reservoir whose width is about 1.8 km, it could reach 28 cm. Waveform portion selection using DEM information showed its ability to retrieve water levels over small-to medium-sized water bodies. In contrast, the accuracy over water bodies located in mountainous areas still needs to be improved, but the results demonstrate the possibility of retrieving the water levels over a very small water body with a width of approximately 130 m.

Overall, the Sentinel-3 SRAL has been proven to work over inland water bodies, even those with widths as small as 130 m. Retracking using the waveform portion selection method with an SRTM DEM as a pre-processing to the two-step physical-based retracker, the OCOG retracker and threshold retracker improved the accuracy with an optimal ubRMSE of 16 cm. Further steps need to be taken to explore the possibilities for both tracking water bodies that are not currently being tracked, and for improving the accuracies for very small water bodies with Sentinel-3 altimeter data.

In this chapter, a methodology has been developed to estimate water levels in small water bodies (similar size to the along-track resolution) and/or complex topography environments where many waveforms are expected to be contaminated by land. The study compared the performances of the threshold, OCOG and two-step physical-based retrackers, but no significant difference is found in the results of the three retrackers. However, waveform portion selection using DEM information makes great improvements on the results. The Sentinel-3 SRAL has been proven to work over inland water bodies, even those with widths as small as 130 m.

With the altimetry mission originally designed for the ocean, the accuracy of geophysical corrections over inland water ranges from a few centimetres to tens of centimetres, depending on the size of the water body and wind conditions (Fernandes et With the water levels retrieved, the integration of them with an LSM needs to be done to include the anthropogenic effects within models. In the next chapter (Chapter 7), a simple dam simulation example is given by including the buffering effect of dams in river flow simulations.

Introduction

With water level products available (Chapter 6), there is a need to include them in models to better characterise water processes with the anthropogenic intervention included. Humans play an important role in the hydrological process. However, most of the models, especially those applied at large scales (national, continental or global), do not include human activities like irrigation and dam construction [START_REF] Haddeland | Global water resources affected by human interventions and climate change[END_REF]. This chapter, as an expansion of Chapter 6, presents simple examples to evaluate the buffering effect of dams in river flow simulations.

Dams modify the water cycle in different ways: they are used to remove water from the river; they are used for irrigation; they change the seasonality of the streamflow; they remove floods and maintain low flows, which is what we call the buffering effect of the dam on the flow. In order to be able to include these effects, there are different options; one is to simulate the dam operating rules, which is always difficult; another option is to use observations to force the dam behaviour in the simulation. In the second case, remotely sensed dam levels can be converted to volumes, and can be used to force a simulation. This approach is beneficial for global simulations or for simulations in countries where the dam level data cannot be obtained from the dam operators.

A Land surface model (LSM) coupled to a river routing scheme is considered in this study to generate the simulated river flow as input. The idea of including the buffering effect of dams 125 126 CHAPTER 7. DAM SIMULATION is schematically illustrated in Figure 7.1. If the dam is not simulated in a river flow simulation model, the outflow of a dam is simply the simulated inflow; when considering the effect of the dam, the water volume in the dam is forced by the observed volume converted from the satellite water level observation, and the outflow is what is called modelled outflow. Then, both the simulated outflow and the modelled outflow are evaluated with the observed outflow. In this chapter, first, the models used to generate the simulated inflows are introduced in Section 7.2, including SURFEX LSM forced by SAFRAN meteorological model, who provide meteorological dataset as input, and RAPID river routing model, which is used under the framework of Eaudyssée. Second, the methodology of forcing volumes by observations is presented in Section 7.3. Then, the results are presented and discussed in Section 7.4. With the difficulties of having satellite observations and gauging stations' historical observations simultaneously, the results of the simulation using historical observation data are shown first. Finally, the conclusion of this chapter is presented in Section 7.5.

Models

SURFEX land surface model

LSMs physically simulate the exchange of surface water and energy fluxes at the soil-atmosphere interface, taking into account numerous physical processes at the interface between soil, vegetation and atmosphere. When coupled to a river routing scheme, LSMs are also able to simulate river flow. This converts them into physically based and distributed hydrological models, which may be used to quantify and understand hydrological processes at large scales.

In this study, the SURFEX (SURFace EXternalisée) land-surface modelling platform developed by Météo France [START_REF] Masson | The surfexv7.2 land and ocean surface platform for coupled or offline simulation of earth surface variables and fluxes[END_REF] is used. SURFEX uses the ISBA (Interaction Sol-Biosphère-Atmosphère) land surface scheme (Noilhan and Planton, 1989; Noilhan and Mahfouf , 1996) to describe the vertical processes in the soil column and the vegetation, and generates the outflows that will allow us to simulate the river flow. Using an offline method, the SURFEX LSM is forced by a gridded dataset of meteorological observations.

SAFRAN meteorological analysis system

The SAFRAN (Système d'analyse fournissant des renseignements atmosphériques à la neige) meteorological analysis was developed in order to provide an analysis of the atmospheric forcing in mountainous areas for the avalanche forecasting [START_REF] Durand | A meteorological estimation of relevant parameters for snow models[END_REF][START_REF] Durand | A computer-based system simulating snowpack structures as a tool for regional avalanche forecasting[END_REF]. SAFRAN analyses eight parameters including air temperature, relative air humidity, wind speed, incoming solar radiation, incoming atmospheric radiation, cloudiness, snowfall and rainfall.

SASER hydrological model

SURFEX does not simulate river flow directly. Thus, a modelling platform called Eau-dyssée [START_REF] Saleh | Modeling the impact of in-stream water level fluctuations on stream-aquifer interactions at the regional scale[END_REF]) is used to transport SURFEX's runoff and drainage to the river and then compute the river flow. Eau-dyssée couples existing specialised models to simulate surface and sub-surface water processes. With the lack of routing scheme within SURFEX, the river routing model called RAPID [START_REF] David | River network routing on the nhdplus dataset[END_REF], Routing Application for Parallel Computation of Discharge, is coupled to SURFEX by Eau-dyssée. The whole chain is called SASER (SAfran-Surfex-Eaudyssée-Rapid), which works as a hydrological model to generate river flows.

Methodology

Forcing model

To evaluate the buffering effects of dams, a simple model to force the volumes within dams to the target is developed with the inclusion of a simple dam operation rule considering the dam scheme included within the Soil and Water Assessment Tool (SWAT) model [START_REF] Douglas-Mankin | Soil and water assessment tool (swat) model: Current developments and applications[END_REF]. At each time step, the model modifies the outflows in order to guarantee that we attain the target volume.

The simple dam model is defined as follows: where vol target (d) is the target volume at date (d), qin sim is the simulation inflow, qout is the modelled outflow, and τ is a constant that determines the speed in which the volume converges to the target volume, which allows filtering peak flows. If vol(d) < vol target (d), qout is set to be zero. The water volume is controlled by the maximum volume of the reservoir to make sure the daily volume never exceeds the maximum, to avoid the risk of flooding.

This model has been coded in Python and run separately to the SASER modelling chain. This allows us to easily test the model individually for each dam, without dealing with the complexities of the SASER Fortran code.

Conversion from water level to volume

From satellite altimetry mission, only water levels can be observed. To convert water level products into volumes, the historical gauging measurements including water levels and volumes are used to build a polynomial function. Therefore, the target volumes to be the input of the forcing model can be derived from water levels observed by satellite. One example of the polynomial function is given over the Ebro Reservoir in Figure 7.2. However, this approach would not work in dams where no historical data is available, and thus, in those cases, other approaches should be used to relate volume and level. However, for simplicity's sake, it will be kept like this.

Results

With limited satellite data available from 2016, first, the test was done with historical observed in situ volumes to expand period in Section 7.4.1's simulation. Second, the test was performed with satellite altimetry data as the forcing target, and results are shown in Section 7.4.2. The modelled outflow is evaluated by the Root Mean Square Error (RMSE). In order to validate the model, the RMSE must be smaller than the RMSE of not having a dam.

Simulation

In this part, the simulated inflows were derived using the SAFRAN forced SURFEX, and the target volumes are simulated by historical database from 2003 to 2005. The volumes are extracted every 27 days from the historical observations to be consistent with the time resolution of the Sentinel-3 mission. Six reservoirs are chosen: the Itoiz, Yesa, Mequinenza, Barasona, Rialb and Oliana Reservoirs. The results are shown in Figure 7.2, and the RMSE values are listed in Table 7.1. Obtaining observed flows is not evident at times as sometimes the dams have different tributaries coming, and data is often lacking.

In Figure 7.2, the left column shows the comparison of the observed volumes (vol obs) and the modelled volumes (vol linear). The target volumes (target) are the smooth of the observed volumes extraction by every 27 days. The maximum volume (max) is for controlling the daily volume not exceeding the maximum value, which may bring a risk of flooding. The overfill is the absolute volume of the reservoir, just for showing the capability of maximum water storage.

RESULTS

129

The right column shows the comparison of stream flows including the simulated inflow (qin sim), the observed outflow (q out obs) and the modelled outflow (q out mod linear). The modelled outflow filtered out significant peaks in the simulated inflow and showed a closer trend to the observed outflow. This can also conclude with the RMSE value in Table 7.1. The first column of the RMSE values is between the modelled outflow and observed outflow, whilst RMSE-reference values are between the simulated outflow (which is equal to the simulated inflow) and observed outflow. For all tested reservoirs, the RMSEs of modelled outflow are smaller than the referenced RMSEs, which means the model improves the results of the outflow. One of the limitations of using monthly data is that the sub-monthly changes are sometimes massive and have a significant impact, which can be observed in Barasona Reservoir simulation in Figure 7.2.

In the simulation test, it must be noted that the simulated inflows are sometimes considered not good enough and, in this case, the dam model can do little to improve the simulation. Also, in some cases, the peaks are not well filtered (for example in Barasona the first peaks should not have been filtered). This is difficult with such a simple model, but in general what we see is an improvement in low flows and filtering of peaks, which are the main buffering effects of dams. , where vol obs is the observed volumes, target is the target volume, vol linear is the modelled volumes, max is the maximum volume, overfill is the absolute volume of the reservoir, qin sim is the simulated inflow, q out obs is the observed outflow, q out mod linear is the modelled outflow.

Satellite data forcing

In Ebro River basin, the water levels at dams can be observed by Sentinel-3 satellite including for the Ebro, Itoiz, Irabia, Sotonera and Ribbarroja Reservoirs. However, most of the reservoirs have problems in different aspects. For example, the Ebro Reservoir is very complex with more than one upper stream, but only one has observations. The LSM could not give near realistically simulated inflows. Over Irabia Reservoir, observations do not exist after 2016. Sotonera Reservoir is also very complex, which would require important modifications in our simple model. Ribbarroja Reservoir is not considered interesting as the dam is operated in a very simple way; the water level is almost constant, making the dam transparent (the inflow is equal to the outflow). With the only reservoir left, which is Itoiz, the satellite forcing data is involved. However, the simulated inflows are not satisfying for Itoiz Reservoir. Therefore, we decide to use the observed inflows instead of simulated inflows.

Figure 7.3 shows the results with target volumes converted from satellite water levels over Itoiz Reservoir. The results show that still there are some big peaks in the modelled outflow. In order to avoid big peaks and test the sensitivity of the constant release pace for the outflow, and instead of releasing the volumes within one week, different possibilities are tested. Table 7.2 shows the RMSEs of releases of the volumes at a constant speed within two weeks, three weeks and four weeks. When setting the time span with larger values from 1 week to 3 weeks, the RMSE decreases, whilst peaks in the modelled outflows are better filtered. However, the decrease no longer continues after that. To evaluate the improvement of forcing volumes using satellite data, the comparison between using the observed smoothing volumes as the target and using the satellite converted volumes as the target has been done. The results are shown in Figure 7.4. The RMSE of using satellite data is 1.654, and the RSME for using observed volumes as forcing data is 1.659. The performance of using volumes converted from satellite water level as target forcing is satisfied compared with results using observations as the target.

Conclusion

In this chapter, we have shown that a very simple model, that only tries to keep the volume close to the observed one (at the monthly scale), and that filters peak flows by releasing the water slowly (in a few days), is able to improve the outflows greatly, producing more realistic downstream river flows. It must be taken into account that our approach is still only a prototype, but it has shown the value of the remote sensing of dam levels for hydrological modellers.

Forcing the volumes using the observed volumes as the target and using volumes converted from satellite water levels as the target is tested. A straightforward linear model is considered as the dam operation rule, and results show that using target volume to force the dam volumes improved the simulation outflow. Due to differing circumstances (e.g. lack of observation data etc.) dam forcing using satellite data is only performed over Itoiz reservoir. Good performance is shown compared with using in situ gauging as target volumes. However, the simulated inflows are sometimes not good enough, and some of the peaks are difficult to fully filtered with such a simple model. In general though, the improvement in low flows and filtering of peaks, which are the main buffering effects of dams, are observed using target volume forcing compared to simulated inflows. The very simple model is still in development, and more tests need to be done to prove these findings in other river basins and to improve the performance of the dam model.

The following shows how this work should be expanded: 1. increasing the number of dams by including dams outside the Ebro basin. 2. integrating the extractions in those dams where they are important.

3. calibrate the time constant (τ ) using long time series of observations. 4. improve the model by making more realistic assumptions concerning its behaviour during dry, normal and wet periods (flood control).

5. rewrite the simple model in Fortran and integrate it into the SASER chain.

Chapter 8

Conclusions and Perspectives

Main conclusions

During the last few years, the Sentinel missions are to replace the current older Earth observation missions, which have reached retirement or are currently nearing the end of their operational life span. With the most up-to-date sensors on board, Sentinel missions bring more possibilities to observe the Earth surface at a higher resolution. In this framework, spaceborne synthetic aperture radar (SAR) has proven to be a suitable tool for inland water resources monitoring, regardless of daylight and weather conditions, and thus is considered primarily in this PhD study.

Human intervention in the water cycle is important both with respect to irrigation and dam construction. The water resources need to be assessed, and the impact of human activities need to be evaluated and included in the modelling. With the work described in this thesis, the feasibility to monitor the water resources in an environment heavily influenced by humans has been demonstrated, with the capability of soil moisture retrieval and water level retrieval from satellite data. This thesis intends to contribute to the development of methodologies to monitor soil moisture and water levels at a better resolution with the use of SAR and altimetry, both of which can work under any weather condition. Furthermore, to include the human interventions in the water cycle, this thesis has also contributed to evaluating irrigation by developing the methodology of irrigation mapping and involving dams using the Land Surface Model (LSM). The main and original contributions and conclusions of this thesis are outlined below:

• Soil moisture retrieval. Two methodologies to retrieve the soil moisture by using the relationship between the SAR backscatter and NDVI has been developed. Method 1 models the backscattering difference with the driest value, whereas method 2 is based on the difference between radar signals observed on two consecutive dates, meaning that the radar signals are influenced by much smaller changes in vegetation. Both methods are found to predict soil moisture variations that are well correlated with rainfall events. Method 2 is found to be more robust than method 1, since it does not require searching for the minimum value in each pixel, which can introduce larger errors under extreme local conditions. These results demonstrate the potential of Sentinel-1 data for the retrieval of 100 m (or even better) 136 CHAPTER 8. CONCLUSIONS AND PERSPECTIVES resolution soil moisture. Both methods can be applied to any vegetation-covered area for which time-series of SAR and optical data have been recorded.

• Irrigation mapping. A new methodology for irrigation mapping using the Sentinel-1 SAR data time series is developed. The backscatter mean value, the signal variance, and the correlation length are derived from the backscatter signal time series and are used for classification. The result shows a good overall accuracy both using the Support Vector Machine (SVM) and Random Forest (RF). The mean value of the backscatter time series is the key to separate irrigated and non-irrigated fields, while the correlation length and variance are used for the separation of irrigated trees and irrigated crops. With the advantage of SAR being able to work under any weather conditions, this method can be applied to areas with frequent cloud cover. Additionally, this method does not need to develop operational algorithms to estimate soil moisture before application; it is based directly on radar signal analysis. The results demonstrated the potential of using Sentinel-1 data for irrigation mapping at the field scale.

• Water level retrieval. A methodology of estimating water levels in small water bodies has been developed with the novel DEM-oriented waveform portion selection method. The performances of the threshold, OCOG and 2-step physical-based retrackers are compared, but no significant differences were found in the results of the three retrackers. However, the main improvements that have been found show that selecting the portion of the waveform to be fitted affects the retrievals in a significant way. The results show good agreement with the in situ measurements. Retracking the water levels from Level-1 waveforms using the combination of a retracker and the waveform portion selection method is more robust to land contamination than using Sentinel-3 Level-2 data directly and with better accuracies. The Sentinel-3 SAR altimeter (SRAL) has been proven to work over inland water bodies, even those with widths as small as 130 m.

• Dam simulation. The effects of dams have been tested in the river flow simulations by forcing the dam volumes to the target. We have shown that even with a very simple model that tries to keep the monthly dam volume very close to the observed one, and while filtering peak flows, this is enough to significantly improve streamflows in rivers affected by dams. However, due to different kinds of reasons (e.g. lack of observation data), dam forcing using satellite data is only performed over Itoiz reservoir, but we were able to show the value of the satellite-derived monthly dam level data in the context of hydrological modelling.

Future research lines

Various research perspectives in connection with the work carried in this doctoral research can be considered for the future:

• Soil moisture retrieval using both VV and VH polarization. In the present study, data derived from the VV polarization was analysed, since it is more sensitive to soil conditions. However, Sentinel-1 provides data in both VV and VH polarization modes, and it is planned to include VH polarization analyses in future studies, since this operational mode is highly sensitive to the influence of vegetation, and can be used to discriminate between the effects of vegetation. In future, the statistical analysis should be improved by using a larger number of data acquisitions from the Sentinel-1 time series. Besides, with the truth that vegetation effects on radar signal come both from the amount of biomass and vegetation 8.2. FUTURE RESEARCH LINES 137 structure, a better vegetation index needs to be investigated instead of NDVI which is only sensitive to the amount of photosynthetically active vegetation.

• Irrigation mapping using 100-m resolution cells. In the present study, the ground field information is used to average the backscatters within agriculture fields. However, the methodology could not directly be applied in areas where there is no ground field information available. Instead of using field boundary information, field segmentation could be performed using Sentinel-2 NDVI data (Comaniciu andMeer , 1999, 2002;Fukunaga and Hostetler , 1975), but only for areas without frequent cloud cover. Another alternative method is using 100-m resolution cells directly to calculate the metrics only from Sentinel-1 SAR data, which is more consistent with our original intention to use the methodology under any weather conditions. This need to be tested in the future study.

• Improving water level retrieval over very small water bodies and assess the influence of geophysical corrections. In the present water level study, waveform portion selection using DEM information shows its capability to retrieve water levels over small to medium-sized water bodies. In contrast, the accuracy over water bodies located in mountainous areas still needs to be improved. Further steps need to be taken to explore the possibilities for both tracking water bodies that are not currently being tracked, and for improving the accuracies for very small water bodies with Sentinel-3 altimeter data. The drifting effects of satellite tracks on a daily basis makes the footprints always changing in their locations for each water body, resulting in uncertainties in the validation. Also, with the altimetry mission originally designed for the ocean, the accuracy of geophysical corrections over inland water ranges from a few centimetres to tens of centimetres, depending on the size of the water body and wind conditions (Fernandes et al., 2014;[START_REF] Vignudelli | Coastal Altimetry[END_REF][START_REF] Zhang | Application of retracked satellite altimetry for inland hydrologic studies[END_REF]. These effects will most definitely influence the final results and need to be further considered quantitatively.

• Evaluate the buffering effects of dams with more datasets. In the present study of the inclusion of the buffering effects of dams in river flow simulations, there are only five reservoirs that can be retrieved from Sentinel-1 altimetry data. Moreover, because of differing reasons (e.g. lacking observations, the complexity of reservoirs with canals of irrigation, etc.), target volumes converted from satellite-derived water levels are only forcing in the Itoiz Reservoir. To better evaluate the performance of satellite data forcing, more tests need to be done in other river basins where both satellite data and gauging stations are available.

• Involving both irrigation and the buffering effects of dams in LSMs. Most of the recent LSMs do not include anthropogenic interventions in the water process. To improve the models to more realistic ones, both anthropogenic interventions of irrigation and dam construction should be considered. However, it is a very complex work since neither of them is easy to be modelled. This doctoral research has contributed to part of this goal, yet more works need to be done in the future study to integrate the results of this PhD study in the modelling.

In summary, the work developed in the present thesis has demonstrated the capability of monitoring water resources from space. This doctoral activity compiles multidisciplinary studies covering soil moisture study, altimetry study, data processing, signal processing, as well as hydrological modelling, showing the complexity and diversity of the hydrological related studies. The diversification allowed the identification of different issues to be tackled in future research, and also showed the potential of using satellites to understand our planet in a better way.

5 1. 1 . 4 6 1. 2 7 1. 2 . 1 8 1. 2 . 2 9 1. 2 . 3 9 1. 3 10 1. 3 . 1 10 1. 3 . 2

 514627218229239310311032 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.1.1 Climate change and interaction with water resources . . . . . . . . . . . 2 1.1.2 Water resources in global scale . . . . . . . . . . . . . . . . . . . . . . . 4 1.1.3 Water resources in regional scale . . . . . . . . . . . . . . . . . . . . . . Water cycle and anthropogenic interventions . . . . . . . . . . . . . . . Water resources estimation . . . . . . . . . . . . . . . . . . . . . . . . Hydrological models and land surface models . . . . . . . . . . . . . . . Measurements and remote sensing techniques . . . . . . . . . . . . . . . Challenges . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . About the thesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Objectives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Structure of the thesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11

2 CHAPTER

 2

  Figure 1.1 (b) (Pachauri and Meyer , 2014).

Figure 1 . 1 :

 11 Figure 1.1: Change in average surface temperature (a) and change in average precipitation (b) based on multi-model mean projections for 2081-2100 relative to 1986-2005 under the RCP2.6 (left) and RCP8.5 (right) scenarios. RCP is short for Representative Concentration Pathway, which is a greenhouse gas concentration (not emissions) trajectory adopted by the IPCC for its fifth Assessment Report in 2014 (Pachauri and Meyer , 2014).

Figure 1 . 3 :

 13 Figure 1.3: The global map of water risk index considering both water stress and flood risk (Gassert et al., 2015).

Figure 1 . 4 :

 14 Figure 1.4: Processes and pathways of the water cycle (The Drainage Basin Hydrological Cycle, 2018).

  . The passive satellite Soil Moisture and Ocean Salinity (SMOS) launched November 2nd, 2009 (Kerr et al., 2010) and the Soil Moisture Active Passive (SMAP) launched January, 31st 2015 (Entekhabi et al., 2010) provide a global mapping of surface soil moisture based on radiometric measurements at L-band (21 cm, 1.4 GHz). Alternative satellite-based soil moisture products, including products generated by the Advanced Microwave Scanning Radiometer for the Earth observing system (AMSR-E) onboard CHAPTER 1. INTRODUCTION

  al., 2015; Siebert et al., 2005; Thenkabail et al., 2009; Meier et al., 2018). The global irrigated area mapping (GIAM) undertaken by Thenkabail et al. (2009) uses the combination of meteorological data, land use classification information (forest) and remote sensing data from multiple satellite sensors, and gives spatial resolution about 10 km. Salmon et al. (Salmon et al., 2015) combined national surveys, climate and remote sensing data. However, all these studies are not able to depict individual farm fields. With the diverse range of irrigated field size and scattered distribution, it is difficult to map irrigation fields by satellite remote sensing because of its relatively coarse spatial resolution in comparison to field scale.

Figure 2 . 1 :

 21 Figure 2.1: The electromagnetic spectrum. Modified from GSP (2018).

Figure 2 . 2 :

 22 Figure 2.2: The EM wave illustration. Modified from Jensen (2005).

Figure 2 . 3 :

 23 Figure 2.3: One-way transmission rate (%) of microwave through vapor clouds, ice clouds, and rain as a function of frequency (and wavelength) (Ouchi , 2013).

Figure 2 . 4 :

 24 Figure 2.4: Geometry of radar equation (Ulaby et al., 1986b).

Figure 2 . 5 :

 25 Figure 2.5: SAR geometry (Moreira et al., 2013).

Figure 2 .

 2 Figure 2.5 illustrates the typical SAR geometry, where the platform moves in the azimuth or along-track direction, and slant range is the direction perpendicular to the radar's flight path. The swath width gives the ground-range extent of the radar scene, which depends on the measurement taken duration, i.e., how long the radar is turned on (Moreira et al., 2013). Compared to Real Aperture Radar (RAR), SAR synthetically increases the antenna's size to

Figure 2 . 6 :

 26 Figure 2.6: Summary of SAR processing steps. In the first step, the range of compressed data results from a convolution of the raw data with the range reference function. In a second step, the azimuth compression is performed through convolution with the azimuth reference function, which changes from near to far range. The "*" represents the convolution operation (Moreira et al., 2013).

Figure 2 . 7 :

 27 Figure 2.7: Scattering Mechanisms (ESA, 2007).

Figure 2 . 8 :

 28 Figure 2.8: Radar topographic effects (Braun and Hochschild , 2017).

  et al., 2007; Calman and Seyler , 2006; Cretaux and Birkett, 2006).

Figure 2 . 9 :

 29 Figure 2.9: Illustration of altimetry principle (Fricker , 2009).

Figure 2 . 10 :

 210 Figure 2.10: Waveform retracking. Created based on Deng and Featherstone (2006).

CHAPTER 2 .

 2 Figure 2.11: Beam-limited, pulse-limited and Doppler-limited altimetry (Sørensen, 2016).

Figure 3 . 1 :

 31 Figure 3.1: The Ebro river basin. Modified from Ebro Basin (2007).

Figure 3 . 2 :

 32 Figure 3.2: Study area located in Urgell, Catalonia. The bigger red square in the left figure shows the 60 km by 60 km square for soil moisture study at 1 km resolution; the inside smaller red square shows the 20 km by 20 km square for soil moisture study at 100m resolution and irrigation study at field scale. The blue polygon shows the Segarra-Garrigues (SG) system being developed.

CHAPTER 3 . 2 Figure 3 . 3 :

 3233 Figure 3.3: Sentinel-3 tracking modes (open loop and closed loop) per observed surface (ESA, 2018b).

Figure 3 . 4 :

 34 Figure 3.4: Demonstration fields (yellow stars) within the 20 km by 20 km study area: Foradada and Agramunt.

3. 3 .

 3 IN SITU DATA 39 satellite estimations.

Figure 3 . 5 :

 35 Figure 3.5: The meteorological stations of Catalonia in the Ruralcat database. Modified from Ruralcat (2018).

CHAPTER 3 .Figure 3 . 6 :

 336 Figure 3.6: The SIGPAC interface (DARP , 2018).

Figure 3 . 7 :

 37 Figure 3.7: Examples of the SAIH Ebro reservoirs (SAIH Ebro, 2018).

  , and several studies have obtained quite accurate estimations of the soil moisture(Baghdadi et al., 2002; Ulaby et al., 1986a; Sahebi et al., 2004; Santi et al., 2013; Pham-Duc et al., 2017; Notarnicola et al., 2008; Pasquariello et al., 1997; Baghdadi et al., 2016; El Hajj et al., 2016; Baghdadi et al., 2012; Santi et al., 2016; Paloscia et al., 2013). However, this methodology is highly sensitive to the input parameters (Ahmed et al., 2011; Ghedira et al., 2004), and over vegetated areas the results can become highly variable due to high rates of evaporation (Bousbih et al., 2018; Hassan-Esfahani et al., 2015; Mishra et al., 2018). The WCM approach (Attema and Ulaby, 1978) relates the backscattering coefficient to soil moisture content and the presence of vegetation, making it possible to be used over densely vegetated areas. Many studies have shown the possibilities of retrieving soil moisture at high accuracies (Zribi et al., 2011; Attema and Ulaby, 1978; Chai et al., 2015; He et al., 2014; Wang et al., 2011; Kumar et al., 2015). However, typical input vegetation parameters include the albedo of the vegetation, and the attenuation factor is difficult to define.

Figure 4 . 1 :

 41 Figure 4.1: Illustration of the relationship between NDVI and Mv used in method 1.

Figure 4 . 2 :

 42 Figure 4.2: Illustration of method 2.

Figure 4 . 3 :

 43 Figure 4.3: Intercomparison between ground measurements in the two demonstration fields of Agramunt and Foradada and Sentinel-1 moisture estimations based on method 1 (a) and method 2 (b).

Figure 4 . 4 :

 44 Figure 4.4: Retrieved soil moisture [m 3 /m 3 ] maps of method 1 and method 2.

Figure 1 .

 1 Figure 1. Study area in Urgell, Catalunya (20 km × 20 km) with ground measurements (yellow stars).

Figure 1 .

 1 Figure 1. Study area in Urgell, Catalunya (20 km × 20 km) with ground measurements (yellow stars).

  is a NDVI map of our study area, computed from Sentinel-2 data recorded on 25 August 2015, and characterized by a dynamic range between 0.1 and 0.8. The temporal variations in NDVI during 2016 are shown for two different locations:

  is a NDVI map of our study area, computed from Sentinel-2 data recorded on 25 August 2015, and characterized by a dynamic range between 0.1 and 0.8. The temporal variations in NDVI during 2016 are shown for two different locations:

Figure 2 .

 2 Figure 2. Sentinel-2 NDVI map in the study area (left), NDVI time series for the dry land site (upper right), and NDVI time series for the site in the old irrigated area (lower right).

Figure 2 .

 2 Figure 2. Sentinel-2 NDVI map in the study area (left), NDVI time series for the dry land site (upper right), and NDVI time series for the site in the old irrigated area (lower right).

  , which corresponds to the maximum value of backscattering difference under the driest, bare-soil conditions.

Figure 3 .Figure 3 .

 33 Figure 3. Illustration of the relationship between NDVI and Mv used in method 1.

Figure 4 .

 4 Figure 4. Illustration of the processed radar signal differences (dB) for all dates, with the driest radar signals shown as a function of NDVI for all (100 m × 100 m) cells in the Urgell area. Each point corresponds to a single radar signal difference (i, j) NDVI Δσ

  low-resolution moisture products (SMOS Level 3 daily product), corresponding to the two-year period of S1 acquisitions, were used to estimate max Mv and min Mv

  32 m 3 /m 3 . As shown shown in Figure 5, the value of min Mv is taken to be ≈ 0.05 m 3 /m 3 .

Figure 4 .

 4 Figure 4. Illustration of the processed radar signal differences (dB) for all dates, with the driest radar signals shown as a function of NDVI for all (100 m × 100 m) cells in the Urgell area. Each point corresponds to a single radar signal difference ∆σ NDVI(i,j)

Figure 5 .

 5 Figure 5. Mean S1 radar signal as a function of the SMOS soil moisture computed over a single SMOS pixel (40 km × 40 km). The radar signal saturates beyond soil moisture levels of 0.32 m 3 /m 3 .

Figure 5 .

 5 Figure 5. Mean S1 radar signal as a function of the SMOS soil moisture computed over a single SMOS pixel (40 km × 40 km). The radar signal saturates beyond soil moisture levels of 0.32 m 3 /m 3 .

Figure 6 .

 6 Figure 6. Illustration of method 2.

Figure 7 .

 7 Figure 7. Illustration of the radar signal difference (dB) computed for two consecutive dates, as a function of NDVI over the Urgell site. Each point corresponds to a single cell (i,j). For each value of NDVI, the green points indicate the upper decile of the corresponding differences in radar signal.In the case of the maximum value of soil moisture change

Figure 6 .

 6 Figure 6. Illustration of method 2.

Figure 6 .

 6 Figure 6. Illustration of method 2.

Figure 7 .

 7 Figure 7. Illustration of the radar signal difference (dB) computed for two consecutive dates, as a function of NDVI over the Urgell site. Each point corresponds to a single cell (i,j). For each value of NDVI, the green points indicate the upper decile of the corresponding differences in radar signal.In the case of the maximum value of soil moisture change

Figure 7 .

 7 Figure 7. Illustration of the radar signal difference (dB) computed for two consecutive dates, as a function of NDVI over the Urgell site. Each point corresponds to a single cell (i,j). For each value of NDVI, the green points indicate the upper decile of the corresponding differences in radar signal.

Figure 8 .

 8 Figure 8. Intercomparison between ground measurements and S1 moisture estimations based on method 1, for the case of two demonstration fields, at Agramunt and Foradada.

Figure 8 .

 8 Figure 8. Intercomparison between ground measurements and S1 moisture estimations based on method 1, for the case of two demonstration fields, at Agramunt and Foradada.

Figure 8 .Figure 9 .

 89 Figure 8. Intercomparison between ground measurements and S1 moisture estimations based on method 1, for the case of two demonstration fields, at Agramunt and Foradada.

4. 5 Figure 9 .

 59 Figure 9. Temporal variations in ground measurements and S1 estimations of soil moisture at the Agramunt site (a) and Foradada site (b).

Figure 10 .

 10 Figure 10. Intercomparison between ground measurements in the two demonstration fields of Agramunt and Foradada and S1 moisture estimations based on method 2.

Figure 9 .

 9 Figure 9. Temporal variations in ground measurements and S1 estimations of soil moisture at the Agramunt site (a) and Foradada site (b).

Sensors 2017 ,Figure 9 .

 20179 Figure 9. Temporal variations in ground measurements and S1 estimations of soil moisture at the Agramunt site (a) and Foradada site (b).

Figure 10 .

 10 Figure 10. Intercomparison between ground measurements in the two demonstration fields of Agramunt and Foradada and S1 moisture estimations based on method 2.

Figure 10 .Figure 11 .

 1011 Figure 10. Intercomparison between ground measurements in the two demonstration fields of Agramunt and Foradada and S1 moisture estimations based on method 2.

Figure 11 .

 11 Figure 11. Temporal variations in ground measurements and S1 estimations of soil moisture at the Agramunt site (a) and Foradada site (b), determined using method 2.

Figure 12 .Figure 13 .

 1213 Figure 12. Retrieved soil moisture maps obtained using methods 1 and 2.

Figure 12 .

 12 Figure 12. Retrieved soil moisture maps obtained using methods 1 and 2.

Figure 12 .

 12 Figure 12. Retrieved soil moisture maps obtained using methods 1 and 2.

Figure 13 .

 13 Figure 13. Difference of soil moisture retrieved by method 2 and 1 for date 1 (21 August 2015) and date 2 (2 September 2015).

Figure 13 .

 13 Figure 13. Difference of soil moisture retrieved by method 2 and 1 for date 1 (21 August 2015) and date 2 (2 September 2015).

Figure 14 .

 14 Figure 14. RMS between methods 1 and 2, as a function of NDVI for date 21 August 2015 (a) and date 2 September 2015 (b).

Acknowledgments:

  This study is part of the REC project partially funded by the European Commission Horizon 2020 Programme for Research and Innovation (H2020) in the context of the Marie Skłodowska-Curie Research and Innovation Staff Exchange (RISE) action under grant agreement No. 645642 and TOSCA/CNES ASCAS projects. The authors wish to thank the technical teams for their support and contributions to the successful ground campaigns.

Figure 14 .

 14 Figure 14. RMS between methods 1 and 2, as a function of NDVI for date 21 August 2015 (a) and date 2 September 2015 (b).

  using SRTM (Shuttle Radar Topography Mission) DEM (Digital Elevation Model) at 30 m.

Figure 1 .

 1 Figure 1. Temporal availability of Sentinel-1 time series data used in this study. Blue color indicates Sentinel-1A data and red color indicates Sentinel-1B data.

Figure 2 .

 2 Figure 2. Study area located in Urgell, Catalunya, Spain, and the land-use map.

Figure 3 .

 3 Figure 3. Precipitation event in Baldomar meteo station (1.5 km away from the study area).

Figure 4 .

 4 Figure 4. Workflow overview using a Support Vector Machine (SVM).

Figure 5 .

 5 Figure 5. Selected reference areas for irrigated crops, irrigated trees, and non-irrigated crops.

Figure 6 .

 6 Figure 6. Signal backscatter distribution over selected irrigated areas for vertical-vertical (VV) polarization (a) and vertical-horizontal (VH) polarization (b); signal backscatter distribution over selected non-irrigated areas for VV polarization (c) and VH polarization (d).

Figure 7 .

 7 Figure 7. Autocorrelation function of irrigated and non-irrigated cells (100 m × 100 m) for VV polarization (a) and VH polarization (b).

Figure 8 .

 8 Figure 8. Power spectral density plots of irrigated and non-irrigated cells (100 m × 100 m) for VV polarization (a) and VH polarization (b).

Figure 9 .

 9 Figure 9. The relation of the mean value and the variance for irrigated and non-irrigated areas for VV polarization (a) and VH polarization (b).

Figure 10 .

 10 Figure 10. The relation of the variance and the correlation length for VV polarization (a) and VH polarization (b).

Figure 11 .

 11 Figure 11. Cont.

Figure 11 .

 11 Figure 11. Examples of metrics maps for each field in VV polarization: mean value (a), variance value (b), correlation length (c), and fractal dimension (d).

Figure 12 .

 12 Figure 12. Irrigation mapping after masking out forest, urban areas, and water bodies.

Figure 6 . 1 :

 61 Figure 6.1: Schematic description of OCOG retracker with main parameters, where W is the window size, A is the amplitude, COG is the center of gravity of the window, and LEP is the leading edge position. (Vignudelli et al., 2011).

Figure 6 . 2 :

 62 Figure 6.2: Schematic representation of waveform variability along the satellite track in radar altimetry (Maillard et al., 2015).

CHAPTER 6 .Figure 6 . 3 :

 663 Figure 6.3: Illustration of waveform portion selection. Modified from Makhoul et al. (2018).

Figure 1 .

 1 Figure 1. The Ebro River basin: (a) the Ebro River network in the Iberian Peninsula; and (b) the water bodies covered by Sentinel-3 satellite tracks with all available gauging stations.

Figure 1 .

 1 Figure 1. The Ebro River basin: (a) the Ebro River network in the Iberian Peninsula; and (b) the water bodies covered by Sentinel-3 satellite tracks with all available gauging stations.

Figure 2 .

 2 Figure 2. Water bodies and DEMs with water masks (green polygons), Sentinel-3 tracks (red lines) and gauging stations (yellow location pins): (a) Ebro Reservoir, (b) Ebro Reservoir DEM, (c) Itoiz Reservoir, (d) Itoiz Reservoir DEM, (e) Irabia Reservoir, (f) Irabia Reservoir DEM, (g) Sotonera Reservoir, (h) Sotonera Reservoir DEM, (i) Ribarroja Reservoir, and (j) Ribarroja Reservoir DEM.

Figure 2 .

 2 Figure 2. Water bodies and DEMs with water masks (green polygons), Sentinel-3 tracks (red lines) and gauging stations (yellow location pins): (a) Ebro Reservoir, (b) Ebro Reservoir DEM, (c) Itoiz Reservoir, (d) Itoiz Reservoir DEM, (e) Irabia Reservoir, (f) Irabia Reservoir DEM, (g) Sotonera Reservoir, (h) Sotonera Reservoir DEM, (i) Ribarroja Reservoir, and (j) Ribarroja Reservoir DEM.

Remote 26 Figure 3 .

 263 Figure 3. Operation of the two-step physical-based retracker over peaky-like waveform. r (1st) is the correlation coefficient for the SWH fit with the original waveform, while r (2nd) is the correlation coefficient for the MSS fit.

  𝑃 ( , ) (𝑃 , 𝑘 , 𝑆𝑊𝐻/𝑀𝑆𝑆) = 𝑃 ⋅ 𝐵 ( , ) (𝑀𝑆𝑆) ⋅ 𝑔 (𝑆𝑊𝐻) [𝑓 (𝑔 (𝑆𝑊𝐻) ⋅ (𝑘 -𝑘 )) + 𝑇 ( , ) (𝑀𝑆𝑆) ⋅ 𝑔 (𝑆𝑊𝐻) ⋅ 𝜎 𝐿 ⋅ 𝑓 (𝑔 (𝑆𝑊𝐻) ⋅ (𝑘 -𝑘 ))]

Figure 3 .

 3 Figure 3. Operation of the two-step physical-based retracker over peaky-like waveform. r (1st) is the correlation coefficient for the SWH fit with the original waveform, while r (2nd) is the correlation coefficient for the MSS fit.

Figure 4 .

 4 Figure 4. Overview of the workflow.

Figure 4 .

 4 Figure 4. Overview of the workflow.

Figure 5 .

 5 Figure 5. Comparison of the waveform portion selection method with the method using the whole waveform for the two-step physical-based retracker and OCOG retracker: (a) two-step physical-based retracker with the selected waveform portion; (b) two-step physical-based retracker with the whole waveform; (c) OCOG retracker with the selected waveform portion; and (d) OCOG retracker with the whole waveform. Red stars mark the selected portion of the waveform. r (1st) is the correlation coefficient for the SWH fit with the original waveform, and r (2nd) is the correlation coefficient for the MSS fit.

Figure 5 .

 5 Figure 5. Comparison of the waveform portion selection method with the method using the whole waveform for the two-step physical-based retracker and OCOG retracker: (a) two-step physical-based retracker with the selected waveform portion; (b) two-step physical-based retracker with the whole waveform; (c) OCOG retracker with the selected waveform portion; and (d) OCOG retracker with the whole waveform. Red stars mark the selected portion of the waveform. r (1st) is the correlation coefficient for the SWH fit with the original waveform, and r (2nd) is the correlation coefficient for the MSS fit.

Figure 6 .

 6 Figure 6. Time series water levels over the Ebro Reservoir: using the waveform portion selection method (a); and using the whole waveform (b).

Figure 6 .

 6 Figure 6. Time series water levels over the Ebro Reservoir: using the waveform portion selection method (a); and using the whole waveform (b).

Figure 7 .

 7 Figure 7. Cont.

Figure 7 .

 7 Figure 7. Water level time series over different water bodies: (a) Itoiz Reservoir; (b) Irabia Reservoir; (c) Sotonera Reservoir; and (d) Ribarroja Reservoir.

Figure 7 .

 7 Figure 7. Water level time series over different water bodies: (a) Itoiz Reservoir; (b) Irabia Reservoir; (c) Sotonera Reservoir; and (d) Ribarroja Reservoir.

Figure 7 . 1 :

 71 Figure 7.1: Illustration of the forcing dams.

. 2 )•

 2 vol(d) = vol(d -1) + qin sim (d)(7.1)• If dam is more than half-full:qout(d) = [vol(d) -vol target (d)]/τ (7If dam is less than half-full: qout(d) = [vol(d) -vol target (d)] * 30%/τ (7.3)

Figure 7 . 2 :

 72 Figure 7.2: Polynomial function between water level and volume.

Figure 7 . 2 :

 72 Figure 7.2: Results of forcing dams for volumes (left) and stream flows (right) with historical data from 2003 to 2005, where vol obs is the observed volumes, target is the target volume, vol linear is the modelled volumes, max is the maximum volume, overfill is the absolute volume of the reservoir, qin sim is the simulated inflow, q out obs is the observed outflow, q out mod linear is the modelled outflow.

Figure 7 . 3 :

 73 Figure 7.3: Results of forcing dams for volumes (left) and stream flows (right) using satellite data over Itoiz Reservoir (2016-2018). Legend explained in Figure 7.2.

Figure 7 . 4 :

 74 Figure 7.4: Results of forcing dams using satellite data as target (upper) and using observed volumes as target (bottom) over Itoiz Reservoir (2016-2018).

Contents 8 . 1

 81 Main conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 135 8.2 Future research lines . . . . . . . . . . . . . . . . . . . . . . . . . . . . 136

AcronymsAEMET

  Agencia Estatal de Meteorología ALOS Advanced Land Observation Satellite AMSR-E Advanced Microwave Scanning Radiometer for the Earth observing system ASCAT Advanced SCATterometer CD Change Detection COG Centre of Gravity CSIRO Commonwealth Scientific and Industrial Research Organisation DEM digital elevation model EASAC European Academies Science Advisory Council EC European Commission EM Electromagnetic EO Earth observation EOS Earth Observing System ERS European Remote Sensing Satellite ET Evapotranspiration EU European Union EW Extra-Wide swath FAO Food and Agriculture Organization FFT fast Fourier trasform GCM Global Climate Model GEE Google Earth Engine GIAM global irrigated area mapping GLAS the Geoscience Laser Altimeter System GRD Ground Range Detected 139 Acronyms GWC gravimetric water content ICESat Ice, Cloud,and land Elevation Satellite IPCC Intergovernmental Panel on Climate Change IW Interferometric Wide swath LRM low-resolution mode LSM Land Surface Model MAP Mediterranean Action Plan MetOp Meteorological Operational MIRAS Microwave Imaging Radiometer using Aperture Synthesis MODIS Moderate Resolution Imaging Spectroradiometer MSI Multi Spectral Instrument MSS mean square slope MWR MicroWave Radiometer NASA National Aeronautics and Space Administration NDVI Normalized Difference Vegetation Index NIR Near Infrared NN Neural Network NOAA National Oceanic and Atmospheric Administration NRT Near Real Time NTC Non-Time Critical OCN Ocean OCOG Offset Centre of Gravity OLCI Ocean and Land Color Instrument PDGS Payload Data Ground Segment PRF Pulse Repetition Frequency radar RAdio Detection And Ranging RAR real aperture radar RCP Representative Concentration Pathway RCCI Regional Climate Change Index RF Random Forest Acronyms 141 RMSE Root Mean Square Error S2A Sentinel-2A S2B Sentinel-2B SAFRAN Système d'analyse fournissant des renseignements atmosphériques à la neige SAIH Sistema Automático de Información Hidrológica SAR synthetic aperture radar SG Segarra-Garrigues SIGPAC Geographic Information System for Agricultural Parcels SLC Single Look Complex SLSTR Sea and Land Surface Temperature Radiometer SM Stripmap SMAP Soil Moisture Active Passive SMOS Soil Moisture and Ocean Salinity SNR signal-to-noise ratio SRAL Radar Altimeter SRTM Shuttle Radar Topography Mission SSM Surface Soil Moisture SSO Sun-synchronous Orbit STC Slow Time Critical SURFEX SURFace EXternalisée SVM Support Vector Machine SWAT Soil and Water Assessment Tool SWH significant wave height UNESCO United Nations Educational, Scientific and Cultural Organization UNEP United Nations Environment Programme VWC volumetric water content WCM Water Cloud Model WCRP The World Climate Research Programme WV Wave WWF World Wide Fund List of Figures 1.1 Change in average surface temperature (a) and change in average precipitation (b) based on multi-model mean projections for 2081-2100 relative to 1986-2005 under the RCP2.6 (left) and RCP8.5 (right) scenarios. RCP is short for Representative Concentration Pathway, which is a greenhouse gas concentration (not emissions) trajectory adopted by the IPCC for its fifth Assessment Report in 2014 (Pachauri and Meyer , 2014). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 Regional Climate Change Index (RCCI) over 26 land regions of the world calculated from 20 coupled Atmospheric-Ocean General Circulation Models and 3 IPCC emission scenarios (A1B, A2, B1) (Giorgi , 2006; PAGES , 2018). . . . . . 1.3 The global map of water risk index considering both water stress and flood risk (Gassert et al., 2015). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4 Processes and pathways of the water cycle (The Drainage Basin Hydrological Cycle, 2018). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.1 The electromagnetic spectrum. Modified from GSP (2018). . . . . . . . . . . . . 2.2 The EM wave illustration. Modified from Jensen (2005). . . . . . . . . . . . . . . 2.3 One-way transmission rate (%) of microwave through vapor clouds, ice clouds, and rain as a function of frequency (and wavelength) (Ouchi , 2013). . . . . . . . 2.4 Geometry of radar equation (Ulaby et al., 1986b). . . . . . . . . . . . . . . . . . . 2.5 SAR geometry (Moreira et al., 2013). . . . . . . . . . . . . . . . . . . . . . . . . . 2.6 Summary of SAR processing steps. In the first step, the range of compressed data results from a convolution of the raw data with the range reference function. In a second step, the azimuth compression is performed through convolution with the azimuth reference function, which changes from near to far range. The "*" represents the convolution operation (Moreira et al., 2013). . . . . . . . . . . . . 2.7 Scattering Mechanisms (ESA, 2007). . . . . . . . . . . . . . . . . . . . . . . . . . 2.8 Radar topographic effects (Braun and Hochschild , 2017). . . . . . . . . . . . . . . 2.9 Illustration of altimetry principle (Fricker , 2009). . . . . . . . . . . . . . . . . . . 2.10 Waveform retracking. Created based on Deng and Featherstone (2006). . . . . . 2.11 Beam-limited, pulse-limited and Doppler-limited altimetry (Sørensen, 2016). ..3.1 The Ebro river basin. Modified from Ebro Basin (2007). . . . . . . . . . . . . . . 3.2 Study area located in Urgell, Catalonia. The bigger red square in the left figure shows the 60 km by 60 km square for soil moisture study at 1 km resolution; the inside smaller red square shows the 20 km by 20 km square for soil moisture study at 100m resolution and irrigation study at field scale. The blue polygon shows the Segarra-Garrigues (SG) system being developed. . . . . . . . . . . . . . . . . . . 3.3 Sentinel-3 tracking modes (open loop and closed loop) per observed surface (ESA, 2018b). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

  al., 2014; Puma and Cook , 2010; Lee et al., 2011; Douglas et al., 2009) found that irrigation increased evaporation both regionally and globally, and decreased temperature in the local area. Lo et al.

;

  Crawford and Linsley, 1966; Dawdy and O'Donnell , 1965; Sugawara, 1969; Burnash et al., 1973; Burnash, 1995; Arnold et al., 1998; Ma and Cheng, 1998), and to physically based models (Abbott et al., 1986a,b; Beven and Kirkby, 1976, 1979; Beven et al., 1987; Grayson et al., 1992; Fortin et al., 2001a,b; Qu, 2004; Qu and Duffy, 2007; Kumar , 2009; Shi et al., 2013;

  al., 2011; Baghdadi et al., 2002, 2012; Srivastava et al., 2009; Zribi et al., 2007; Zribi and Dechambre, 2003; Zribi et al., 2005; El Hajj et al., 2014). SAR data in the C-band sensor has demonstrated its ability to retrieve soil characteristics over vegetation-covered surfaces (Prevot et al., 1993; De Roo et al., 2001; Sikdar and Cumming, 2001; Gherboudj et al., 2011; Wang et al., 2011; Yu and Zhao, 2011; Zribi et al., 2011; Yang et al., 2012

  Analysis of Retrackers' Performances and Water Level Retrieval over the Ebro River Basin Using Sentinel-3. Remote Sen. 2019, 11, 718. DOI: 10.3390/rs11060718

	Chapter 2
	• Chapter 4, that presents the improvements of soil moisture retrieval, and based partly on published article (Gao et al., 2017) Remote Sensing Principles
	• Chapter 5, presenting the irrigation mapping methodology, based on published article (Gao et al., 2018)
	Contents

• Chapter 6, demonstrating the improved water level retrieval, based on published article (Gao et al., 2019) • Chapter 7, as an expansion of Chapter 6, showing the dam simulation that evaluates the buffering effect of dams in hydrological models • Chapter 8, offering general conclusions and future perspectives The three articles on which Chapter 4, 5, and 6 are based are listed below with details: • Gao, Q.; Zribi, M.; Escorihuela, M.J.; Baghdadi, N. Synergetic Use of Sentinel-1 and Sentinel-2 Data for Soil Moisture Mapping at 100 m Resolution. Sensors 2017, 17, 1966. DOI: 10.3390/s17091966 • Gao, Q.; Makhoul, E.; Escorihuela, M.J.; Zribi, M.; Quintana Seguí, P.; García, P.; Roca, M.

Table 2 .

 2 1: The frequency range and wavelength range of microwave bands (du Preez and Sinha, 2016).

	CHAPTER 2. REMOTE SENSING PRINCIPLES

Ulaby et al., 1986a; Dobson and Ulaby, 1981).

  Study area . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 3.1.1 Ebro basin . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 3.1.2 Urgell site . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31 3.2 Remote sensing dataset . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 3.2.1 Sentinel data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 3.2.2 MODIS data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 3.2.3 DEM data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 3.2.4 SMOS data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37

	Chapter 3
	Study Area and Datasets
	Contents
	3.1 3.3 In situ data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .	37
	3.3.1 Ground soil moisture measurements . . . . . . . . . . . . . . . . . . . . 37
	3.3.2 Meteorological data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
	3.3.3 SIGPAC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
	3.3.4 SAIH Ebro . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
	3.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .	40
	3.1 Study area
	3.1.1 Ebro basin

Table 3 .

 3

	1: Sentinel-1 satellite parameters (Kramer , 2018a)
	Launch date	Sentinel-1A: April 03, 2014 Sentinel-1B: April 22, 2016
	Orbit type	SSO (Sun-synchronous Orbit)
	Revisit time	12 days
	Orbital altitude	693 km
	Sensor complement C-SAR (C-band Synthetic Aperture Radar)
		Stripmap (SM)
	Acquisition mode	Interferometric Wide swath (IW)
		Extra-Wide swath (EW)
		Wave (WV)
	Table 3.2: Sentinel-1 operational modes (ESA, 2018a)
		IW	WV	SM	EW
		Dual	Single	Dual	Dual
	Polarization	(HH+HV,	(HH,	(HH+HV, (HH+HV,
		VV+VH)	VV)	VV+VH) VV+VH)
	Azimuth resolution	20 m	5 m	5 m	40 m
	Ground range resolution	5 m	5 m	5 m	20 m
		L0 RAW	L0 RAW L0 RAW
	Product	L1 SLC L1 GRD L1 GRD L1 SLC	L1 SLC L1 GRD	L1 SLC
		L2 OCN	L2 OCN	L2 OCN	L2 OCN
	-Extra Wide-swath mode (EW)

Table 3 .

 3 3: Sentinel-2 satellite parameters[START_REF] Kramer | Copernicus: Sentinel-2 -the optical imaging mission for land services[END_REF]

	Launch date	Sentinel-2A: June 23, 2015 Sentinel-2B: March 7, 2017
	Orbit type	SSO
	Revisit time	10 days
	Orbital altitude	786 km
	Sensor complement MSI (Multi Spectral Instrument)
	Product Types	Level-1C Level-1A

Table 3 .

 3 4: Sentinel-2 bands[START_REF] Kramer | Copernicus: Sentinel-2 -the optical imaging mission for land services[END_REF]

	Name Scale Resolution	Wavelength	Description
	B1	0.0001 60 meters	443.9nm (S2A) / 442.3nm (S2B)	Aerosols
	B2	0.0001 10 meters	496.6nm (S2A) / 492.1nm (S2B)	Blue
	B3	0.0001 10 meters	560nm (S2A) / 559nm (S2B)	Green
	B4	0.0001 10 meters	664.5nm (S2A) / 665nm (S2B)	Red
	B5	0.0001 20 meters	703.9nm (S2A) / 703.8nm (S2B)	Red Edge 1
	B6	0.0001 20 meters	740.2nm (S2A) / 739.1nm (S2B)	Red Edge 2
	B7	0.0001 20 meters	782.5nm (S2A) / 779.7nm (S2B)	Red Edge 3
	B8	0.0001 10 meters	835.1nm (S2A) / 833nm (S2B)	NIR
	B8a	0.0001 20 meters	864.8nm (S2A) / 864nm (S2B)	Red Edge 4
	B9	0.0001 60 meters	945nm (S2A) / 943.2nm (S2B)	Water vapor
	B10	0.0001 60 meters 1373.5nm (S2A) / 1376.9nm (S2B)	Cirrus
	B11	0.0001 20 meters 1613.7nm (S2A) / 1610.4nm (S2B)	SWIR 1
	B12	0.0001 20 meters 2202.4nm (S2A) / 2185.7nm (S2B)	SWIR 2
	QA10	10 meters		Always empty
	QA20	20 meters		Always empty
	QA60	60 meters		Cloud mask

Two modes of operation are available with the Sentinel-3 mission: High-Resolution Mode/SAR mode and Low-Resolution Mode (LRM) The SRAL mission is always operated at SAR mode. LRM is back-up mode only.

Table 3 .

 3 6: Ground soil moisture measurements in two demonstration fields at Foradada and Agramunt.

	Site	Foradada	Agramunt
	Coordinates	41.866 • N, 1.015 • E	41.782 • N, 1.089 • E
	Field area in hectare	20	2.5
	At soil depths in cm	3, 9, 10, 20	5, 10, 20, 40
	Measurements period	May-August 2015 June-October 2015	February-October 2016 July-November 2016
	Sand, silt, clay in %	41.5, 42.3, 16.2	52.1, 35.3, 12.6
	Irrigation method	Sprinklers	Subsurface drippers
	SSM (min, max) in m3/m3	(0.08, 0.45)	(0.04, 0.28)
	Meteorological station	Baldomar station	Tornabous station

Table 1 .

 1 Ground soil moisture measurements in two demonstration fields at Foradada and Agramunt.

	Site	Foradada	Agramunt
	Coordinates At soil depths in cm	41.866 • N, 1.015 • E 3, 9, 10, 20	41.782 • N, 1.089 • E 5, 10, 20, 40
	Period of ground measurements Sand, silt, clay in %	May-August 2015 February-October 2016 41.5, 42.3, 16.2	June-October 2015 July-November 2016 52.1, 35.3, 12.6
	Irrigation method Surface soil moisture (min, max) in m 3 /m 3 Meteorological station	Sprinklers (0.08, 0.45) Baldomar station	Subsurface drippers (0.04, 0.28) Tornabous station
	2.2.2. Satellite Data
	(1) Sentinel-1 data

Table 2 .

 2 Sentinel-1 database.

	Date	Date	Date	Date	Date	Date
	16 July 2015	13 Nov. 2015	05 Feb. 2016	29 Apr. 2016	03 Aug. 2016	26 Oct. 2016
	28 July 2015	25 Nov. 2015	17 Feb. 2016	11 May 2016	15 Aug. 2016	07 Nov. 2016
	09 Aug. 2015	07 Dec. 2015	29 Feb. 2016	23 May 2016	27 Aug. 2016	19 Nov. 2016
	21Aug. 2015	19 Dec. 2015	12 Mar. 2016	04 June 2016	08 Sept. 2016	-
	02 Sept. 2015	31 Dec. 2015	24 Mar. 2016	28 June 2016	20 Sept. 2016	-
	14 Sept. 2015	12 Jan. 2016	05 Apr. 2016	10 July 2016	02 Oct. 2016	-
	26 Sept.2015	24 Jan. 2016	17 Apr. 2016	22 July 2016	14 Oct. 2016	-
	(2) Sentinel-2 data

Table 2 .

 2 Sentinel-1 database.

	Date	Date	Date	Date	Date	Date
	16 July 2015	13 Nov. 2015	05 Feb. 2016	29 Apr. 2016	03 Aug. 2016	26 Oct. 2016
	28 July 2015	25 Nov. 2015	17 Feb. 2016	11 May 2016	15 Aug. 2016	07 Nov. 2016
	09 Aug. 2015	07 Dec. 2015	29 Feb. 2016	23 May 2016	27 Aug. 2016	19 Nov. 2016
	21Aug. 2015	19 Dec. 2015	12 Mar. 2016	04 June 2016	08 Sept. 2016
	02 Sept. 2015	31 Dec. 2015	24 Mar. 2016	28 June 2016	20 Sept. 2016
	14 Sept. 2015	12 Jan. 2016	05 Apr. 2016	10 July 2016	02 Oct. 2016
	26 Sept.2015	24 Jan. 2016	17 Apr. 2016	22 July 2016	14 Oct. 2016

Table 3 .

 3 Sentinel-2 database.

	Date	Date	Date	Date	Date	Date
	06 July 2015	21 Oct. 2015	19 Mar. 2016	21 May 2016	30 July 2016	28 Sept. 2016
	16 July 2015	20 Nov. 2015	22 Mar. 2016	28 May 2016	06 Aug. 2016	05 Oct. 2016
	02 Aug. 2015	30 Nov. 2015	29 Mar. 2016	07 June 2016	09 Aug. 2016	15 Oct. 2016
	05 Aug. 2015	03 Dec. 2015	01 Apr. 2016	10 June 2016	16 Aug. 2016	18 Oct. 2016
	12 Aug. 2015	23 Dec. 2015	08 Apr. 2016	20 June 2016	19 Aug. 2016	25 Oct. 2016
	15 Aug. 2015	30 Dec. 2015	11 Apr. 2016	27 June 2016	26 Aug. 2016	28 Oct. 2016
	22 Aug. 2015	12 Jan. 2016	18 Apr. 2016	30 June 2016	29 Aug. 2016	04 Nov. 2016
	25 Aug. 2015	19 Jan. 2016	28 Apr. 2016	07 July 2016	05 Sept. 2016	07 Nov. 2016
	11 Sept. 2015	29 Jan. 2016	01 May 2016	10 July 2016	08 Sept. 2016	14 Nov. 2016
	14 Sept. 2015	18 Feb. 2016	08 May 2016	17 July 2016	15 Sept. 2016	17 Nov. 2016
	24 Sept. 2015	09 Mar. 2016	11 May 2016	20 July 2016	18 Sept. 2016	24 Nov. 2016
	01 Oct. 2015	12 Mar. 2016	18 May 2016	27 July 2016	25 Sept. 2016	27 Nov. 2016

  Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

	Chapter 5
	Irrigation Mapping
	Contents
	5.1

73 5.2 Methodology . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74 5.3 Published paper . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75 5.4 Summary and Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . 94

Table 1 .

 1 Confusion matrix.

	Classes	Irrigated Crops	Irrigated Trees	Non-Irrigated Fields
	Irrigated crops	77.5%	1.9%	20.6%
	Irrigated trees	22.3%	73.5%	4.2%
	Non-irrigated fields	14.8%	1.9%	83.3%

Table 2 .

 2 The accuracy with different input metrics.

	Input Metrics	Mean (VV/VH)	Mean (VV/VH) + Variance (VV/VH)	Mean (VV/VH) + Variance (VV/VH) + Correlation Length (VV/VH)
	Accuracy	80.9%	81.9%	82.2%

Table 1 .

 1 Satellite radar altimetry missions.

	Satellite Mission	Mission Period	Inclination (deg)	Revisit Time (Days)	Along-Track Resolution (km)	Ground Track Separation at Equator (km)
	Geosat 1 ERS 2 1 and 2 TOPEX/Poseidon 3 GFO 4 Jason 1, 2, and 3 Envisat 5 CryoSat-2 6 HY-2 7 SARAL 8 Sentinel-3 SWOT 9 JASON-CS 10 / SENTINEL-6	1985-1990 1991-2011 1992-2006 1998-2008 2001-present 2002-2012 2010-present 2011-present 2013-present 2016-present Planned on 2021 Planned on 2022	108 98.5 66 108 66 98.55 92 99.3 98.55 98.6 77.6 66	17 35 10 17 10 30-35 369 Two phases (14 and 168) 35 27 21 10	1.7 1.7 (ocean mode) 3.4 (ice mode) 2.2 1.7 km 2.2 1.7 0.25 (SAR and SARIN) 1.6 (LRM) 1.9 1.4 0.3 0.1 0.3	16 80 315 164 315 80 7.7 100 75 104 0 315

Table 2 .

 2 Selected water bodies.

	Remote Sens. 2019, 11, x FOR PEER REVIEW				6 of 26
	Water Bodies	Coordinates	Width	Satellite Tracks	Tracking Mode	Gauging Station Distance	Average Slope (5 km)
	Ebro Reservoir Itoiz Reservoir Irabia Reservoir Sotonera Reservoir Ribarroja Reservoir Mequinenza Reservoir Cavallers Reservoir San Salvador Reservoir	(43.0 • N, 3.96 • W) (42.81 • N, 1.37 • W) (42.99 • N, 1.15 • W) (42.12 • N, 0.68 • W) (41.24 • N, 0.40 • E) (41.26 • N, 0.04 • W) (42.59 • N, 0.86 • E) (41.78 • N, 0.20 • E)	1.8 km 400 m-2.7 km depends on satellite tracks 130 m 4.5 km 400 m 600 m 800 m 1.2 km	014 165 186 222 242 279 299 242	Closed loop Closed loop Closed loop Closed loop Open loop Open loop Open loop Open loop	8 km 2 km 450 m 1.5 km 3.5 km 30 km 500 m 2.5 km	4% 19% 20% 3% 24% 3.5% 27% 4.5%

Table 2 .

 2 Selected water bodies.

	104			CHAPTER 6. WATER LEVEL RETRIEVAL
	Water Bodies	Coordinates	Width	Satellite Tracks	Tracking Mode	Gauging Station Distance	Average Slope (5 km)
	Ebro Reservoir	(43.0°N, 3.96°W)	1.8 km	014	Closed loop	8 km	4%
	Itoiz Reservoir	(42.81°N, 1.37°W)	400 m-2.7 km tracks depends on satellite	165	Closed loop	2 km	19%
	Irabia Reservoir	(42.99°N, 1.15°W)	130 m	186	Closed loop	450 m	20%
	Sotonera Reservoir	(42.12°N, 0.68°W)	4.5 km	222	Closed loop	1.5 km	3%
	Ribarroja Reservoir	(41.24°N, 0.40°E)	400 m	242	Open loop	3.5 km	24%
	Mequinenza Reservoir	(41.26°N, 0.04°W)	600 m	279	Open loop	30 km	3.5%
	Cavallers Reservoir	(42.59°N, 0.86°E)	800 m	299	Open loop	500 m	27%
	San Salvador Reservoir	(41.78°N, 0.20°E)	1.2 km	242	Open loop	2.5 km	4.5%
		(a)			(b)
		(c)			(d)
			Figure 2. Cont.

Remote Sens. 2019, 11, x FOR PEER REVIEW 7 of 26

Table 3 .

 3 Geophysical corrections from the Sentinel-3 Level-2 product.

	Correction	Model	Variable of Level-2 Product	Range of Correction
	Dry troposphere Wet troposphere	European Center for Medium-Range Weather Forecasts (ECMWF) model [82] ECMWF model [82]	Mod_dry_tropo_cor_meas_altitude_01 Mod_wet_tropo_cor_meas_altitude_01	1.7-2.5 m 0-50 cm
	Ionosphere	Global Ionospheric Map (GIM) [83]	Iono_cor_gim_01_ku	6-12 cm
	Solid earth tide Geocentric polar tide Ocean loading tide	Cartwright model [84] Historical pole location [85] GOT00.2 model [86]	Solid_earth_tide_01 Pole_tide_01 Ocean_tide_sol1_01	-30 to +30 cm -2 to +2 cm -2 to +2 cm

Table 4 .

 4 Definition of the main parameters of the SAR ocean model[65].

	Parameter	Definition	Formulation
	BW 1 PRF	Receive Bandwidth Pulse Repetition Frequency	--
	N p f c c 0 R H orb H s	Number of Pulses in a burst Carrier frequency Speed of light Range to the surface (distance from satellite to the surface) Orbit Height SWH	------
	σ z

Table 5 .

 5 The best result was obtained for the Ribarroja Reservoir, which had an ubRMSE of approximately 16 cm. The results over the Ebro Reservoir, Sotonera and Ribarroja Reservoirs were also satisfactory. The RMSEs for the Itoiz Reservoir and Irabia Reservoir were relatively large with values greater than 1 m. The results using the three different retrackers do not show large differences except for the Sotonera Reservoir and Ribarroja Reservoir. The OCOG retracker gave the smallest ubRMSE for the Sotonera Reservoir with a value of approximately 38 cm. The two-step physical-based retracker and threshold retracker performed better over the Ribarroja Reservoir than the OCOG retracker, which had an ubRMSE of approximately 16 cm. All retrackers gave better results than the Level-2 ocean retracker without the waveform portion selection method.

Table 5 .

 5 The best result was obtained for the Ribarroja Reservoir, which had an ubRMSE of approximately 16 cm. The results over the Ebro Reservoir, Sotonera and Ribarroja Reservoirs were also satisfactory. The RMSEs for the Itoiz Reservoir and Irabia Reservoir were relatively large with values greater than 1 m. The results using the three different retrackers do not show large differences except for the Sotonera Reservoir and Ribarroja Reservoir. The OCOG retracker gave the smallest ubRMSE for the Sotonera Reservoir with a value of approximately 38 cm. The two-step physical-based retracker and threshold retracker performed better over the Ribarroja Reservoir than the OCOG retracker, which had an ubRMSE of approximately 16 cm. All retrackers gave better results than the Level-2 ocean retracker without the waveform portion selection method.

Table 5 .

 5 Comparison of the water level validation results over all monitored water bodies.

	Water Bodies	Width	Tracking Mode	Two-Step Physical	RMSE/ubRMSE (m) OCOG Threshold	Level-2 Ocean
	Ebro Reservoir	1.8 km	Closed loop	0.32/0.29	0.30/0.28	0.30/0.28	2.18/1.76
	Itoiz Reservoir	400 m-2.7 km	Closed loop	1.18/1.02	1.10/1.00	1.10/1.00	1.41/1.03
	Irabia Reservoir	130 m	Closed loop	1.39/1.39	1.39/1.38	1.39/1.38	1.44/1.39
	Sotonera Reservoir	4.5 km	Closed loop	0.60/0.43	0.49/0.38	0.48/0.44	1.65/1.19
	Ribarroja Reservoir	400 m	Open loop	0.18/0.16	0.29/0.28	0.31/0.16	0.44/0.20
	Mequinenza Reservoir	600 m	Open loop		Off-track
	Cavallers Reservoir	800 m	Open loop		Off-track
	San Salvador Reservoir	1.2 km	Open loop		Off-track

  al., 2014; Vignudelli et al., 2011; Zhang et al., 2010). These effects are important sources of errors and need to be further considered quantitatively.

  Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 125 7.2 Models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126 7.2.1 SURFEX land surface model . . . . . . . . . . . . . . . . . . . . . . . . 126 7.2.2 SAFRAN meteorological analysis system . . . . . . . . . . . . . . . . . 126

	Chapter 7
	Dam Simulation
	Contents
	7.1

7.2.3 SASER hydrological model . . . . . . . . . . . . . . . . . . . . . . . . . 127 7.3 Methodology . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 127 7.3.1 Forcing model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 127 7.3.2 Conversion from water level to volume . . . . . . . . . . . . . . . . . . . 128 7.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128 7.4.1 Simulation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128 7.4.2 Satellite data forcing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132 7.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 133

  extensively used in France (Quintana-Seguí et al., 2008; Vidal et al., 2010) and, more recently, it has also been applied in Spain by Quintana-Seguí et al. (Quintana-Seguí et al., 2016; Quintana-Seguí et al., 2017). The Spanish application of SAFRAN uses observational data from the Spanish State Meteorological Agency (AEMET) and ERA-Interim (Dee et al., 2011) as the first estimate (Quintana-Seguí et al., 2016). The estimates are on a gridded dataset, at a 5 km resolution. In this study, the SAFRAN forced SURFEX is used to provide simulated inflows in the simulation test when no satellite data is available.

	7.3. METHODOLOGY	127
	SAFRAN has been

Table 7 .

 7 1: RMSE of the outflow compared to the observations and RMSE-reference over different reservoirs Dams RMSE [hm 3 /day] RMSE-reference [hm 3 /day]

	Itoiz	1.23044	1.24529
	Yesa	2.2956	2.7005
	Mequinenza	10.91728	12.82188
	Barasona	1.19124	1.28102
	Rialb	1.81069	2.32936
	Oliana	1.7312	1.78994

Table 7 .

 7 2: RMSEs over different time span for releasing the volumes

	1 week 2 weeks 3 weeks 4 weeks Reference
	RMSE 1.65444 1.60315	1.57517	1.58179	3.25177

Sensors 2017, 17, 1966; doi:10.3390/s17091966 www.mdpi.com/journal/sensors

4.5. PUBLISHED PAPER

5.3 Published paper

Remote Sens. 2018, 10, 1495; doi:10.3390/rs10091495 www.mdpi.com/journal/remotesensing

Remote Sens. 2018,10, 1495

5.3. PUBLISHED PAPER

Remote Sens. 2019, 11, 718; doi:10.3390/rs11060718 www.mdpi.com/journal/remotesensing

Remote Sens. 2019, 11, x FOR PEER REVIEW

6.3. PUBLISHED PAPER

Remote Sens. 2019,11, 718

Caloz, R., and C. Collet (2001), Précis de Télédétection. Volume III : Traitements Numériques d'Images de Télédétection, Presses Universitaires du Québec-Agence Universitaire de la Francophonie, Sainte-Foy, Québec.

ACKNOWLEDGEMENTS

First of all, I would like to express my sincere gratitude to my supervisors, Maria Jose Escorihuela, Mehrez Zribi and Pere Quintana Seguí for their guidance and suggestions throughout Acknowledgments: This study is part of the REC project partially funded by the European Commission Horizon 2020 Programme for Research and Innovation (H2020) in the context of the Marie Skłodowska-Curie Research and Innovation Staff Exchange (RISE) action under grant agreement No. 645642 and TOSCA/CNES ASCAS

Acknowledgments:

The authors wish to thank the technical teams for their support and the SIGPAC team, which provides plenty of ground information.

Acknowledgments:

The authors wish to thank the technical teams for their support, the Copernicus Programme, which provided open access to the Sentinel data, and the SAIH Ebro team, who provided the ground information.

grant No. 645642 from European Commission Horizon 2020 Programme for Research and Innovation (H2020) in the context of the Marie Sk lodowska-Curie Research and Innovation Staff Exchange (RISE) action, which gives me the opportunity to experience a very international PhD study.

Again, a great thanks to all my friends both in Europe and in China, who gave me numerous support and encouragement both in life and study.

data was downloaded from Google Earth Engine (GEE) and processed data was downloaded from SIGPAC

Introduction

Dam construction is the most common way for the human to restrict the flow of water, and can be used for flood control, irrigation, hydropower, etc. In the Ebro river basin, the water withdrawals for irrigation mostly comes from the reservoirs either by inundation with channels leading water into the fields, or by sprinkler or drip whose water also come from the reservoirs or rivers. Therefore, the surface water within reservoirs and rivers is vital and needs to be estimated and monitored in a better way.

Satellite altimeters have proven to be valuable tools for monitoring the water levels within inland water systems, including lakes and rivers (Alsdorf et al., 2007;Calman and Seyler , 2006;Cretaux and Birkett, 2006), from which water storage estimates can be derived. Traditionally, satellite altimeters were designed to monitor homogeneous surfaces such as oceans or ice sheets, resulting in poor performance over small inland water bodies due to the contribution from land contamination in the returned waveforms. Moreover, most studies about inland water bodies used low-resolution mode (LRM) altimeters. With the advent of synthetic aperture radar (SAR) altimetry, the along-track spatial resolution of which has been improved a lot, this has enabled the measurement of inland water levels with better accuracy and an increased spatial resolution.

This chapter aims to develop an operational methodology for water level retrieval from Sentinel-3 SAR mode Level-1B waveforms over small to middle-sized water bodies, whose width varies from hundreds of meters (close to along-track resolution) to a few kilometres, and/or surrounded by rapid changing topography. Three specialized algorithms or retrackers, namely the threshold retracker (Davis, 1997), the OCOG retracker (Davis, 1997;Wingham et al., 1986), and our in-house SAR physical-based retracker (Makhoul et al., 2018), together with a new waveform portion selection method, are evaluated in order to minimize land contamination in the waveforms, and to select the nadir return associated with the water body being overflown. The accuracy is assessed for each water body using in situ measurements from the gauging stations in the SAIH system.

Author Contributions: Q.G., E.M. and M.J.E. conceived, designed and implemented the research. Q.G. and E.M. coded the program. Q.G. performed the analysis of the data and drafted the manuscript. M.J.E. assisted in the data analysis and interpretation. All authors reviewed and improved the manuscript. The study was supervised by M.J.E.

Funding: Qi Gao received grant DI-15-08105 from the Spanish Education Ministry (MICINN) and DI-2016-078 from the Catalan Agency of Research (AGAUR). This work was partially funded by the Spanish Ministry of Science, Innovation and Universities and the European Regional Development Fund through grant CGL2017-85687-R.

Conflicts of Interest:

The authors declare no conflict of interest.

List of Tables