Global Analysis of Coastal Gradients of Sea Surface Salinity

: Sea surface salinity (SSS) is a key variable for ocean–atmosphere interactions and the water cycle. Due to its climatic importance, increasing efforts have been made for its global in situ observation, and dedicated satellite missions have been launched more recently to allow homogeneous coverage at higher resolution. Cross-shore SSS gradients can bear the signature of different coastal processes such as river plumes, upwelling or boundary currents, as we illustrate in a few regions. However, satellites performances are questionable in coastal regions. Here, we assess the skill of four gridded products derived from the Soil Moisture Ocean Salinity (SMOS) and Soil Moisture Active Passive (SMAP) satellites and the GLORYS global model reanalysis at capturing cross-shore SSS gradients in coastal bands up to 300 km wide. These products are compared with thermosalinography (TSG) measurements, which provide continuous data from the open ocean to the coast along ship tracks. The comparison shows various skills from one product to the other, decreasing as the coast gets closer. The bias in reproducing coastal SSS gradients is unrelated to how the SSS biases evolve with the distance to the coast. Despite limited skill, satellite products generally agree better with collocated TSG data than a global reanalysis and show a large range of coastal SSS gradients with different signs. Moreover, satellites reveal a global dominance of coastal freshening, primarily related to river runoff over shelves. This work shows a great potential of SSS remote sensing to monitor coastal processes, which would, however, require a jump in the resolution of future SSS satellite missions to be fully exploited. Author Contributions: Conceptualization, G.A.; Data curation, G.A.; Formal analysis, A.N.D. and A.M.D.; Funding acquisition, G.A., A.C.d.S. and A.B.; Investigation, A.N.D. and A.M.D.; Methodology, G.A.; Project administration, A.C.d.S. and A.B.; Software, A.N.D. and A.M.D.; Supervision, G.A., A.C.d.S. and A.B.; Visualization, A.N.D. A.M.D.; Writing—original A.N.D.; Writing—


Introduction
Sea surface salinity (SSS) is an essential climate variable that bears the signature of the water cycle at the ocean-atmosphere interface, where most water fluxes occur [1]. The large-scale spatial patterns of SSS mostly result from the balance between evaporation and precipitation. High SSS is found in the subtropical gyres, where evaporation dominates [2,3], while lower SSS is found around the intertropical convergence zones and at high latitudes, where precipitation dominates. SSS is also affected at high latitudes by the formation or melting of sea ice, which has much smaller salt content than seawater [4]. Ocean circulation and mesoscale eddies also contribute to the transport of water properties and thus to the spatial distribution of SSS (e.g., [5]).

In Situ Data
In situ SSS data used in the present study come from underway TSG measurements available in the global Ocean from 1993 to the present (Figure 1a). The main source is the French TSG database, including measurements from VOS, managed by the SSS Observation Service [17], from research vessels (RV) [37] and sailing ships [38], representing together about half of the global TSG network [39]. Other sources include TSG data from the US-managed Shipboard Automated Meteorological and Oceanographic System initiative [40], VOS managed by the US National Oceanic and Atmospheric Administration (NOAA), by the Japanese National Institute for Environmental Studies and VOS Nippon organization, RV operated by NOAA, the Australian Commonwealth Scientific and Industrial Research Organization and the Japan Agency for Marine-Earth Science and Technology. These measurements are mostly collected from the SBE-21 TSG manufactured by Seabird Electronics. The TSG is generally installed in the engine room of VOS. The depth at which the seawater is collected depends on the ship and its load but is around 5-10 m for VOS and a bit shallower for RV. The temporal resolution is 1 and 5 min for RV and VOS, respectively, which corresponds roughly to a 2-3 km along-track resolution for VOS, typically cruising at 15 knots, and even finer for RV. TSG data density is high in the north Atlantic and western Pacific, and relatively low in the southern hemisphere ( Figure 1a). The present study is based on TSG data available in the nearest 300 km from the coast, with a focus on a selection of five coastal regions: the Amazon River plume, northeast These measurements are mostly collected from the SBE-21 TSG manufactured by Seabird Electronics. The TSG is generally installed in the engine room of VOS. The depth at which the seawater is collected depends on the ship and its load but is around 5-10 m for VOS and a bit shallower for RV. The temporal resolution is 1 and 5 min for RV and VOS, respectively, which corresponds roughly to a 2-3 km along-track resolution for VOS, typically cruising at 15 knots, and even finer for RV. TSG data density is high in the north Atlantic and western Pacific, and relatively low in the southern hemisphere ( Figure 1a). The present study is based on TSG data available in the nearest 300 km from the coast, with a focus on a selection of five coastal regions: the Amazon River plume, northeast Brazil, the Bay of Bengal, the California upwelling system and the Great Australian Bight (Figure 1b; Table 1).  Figure 1a).

Regions Name
Amazon River Plume (AMZ) The Soil Moisture and Ocean Salinity mission is conducted by the European Space Agency (ESA) in collaboration with the Centre National d'Etudes Spatiales (CNES) in France and the Centro para el Desarrollo Technoligica Industrial (CDTI) in Spain [19]. The satellite is equipped with a radiometer that measures the microwaves of electromagnetic radiation emitted by the Earth's surface at L-band (1.4 GHz) using an interferometric radiometer called the Microwave Image Radiometer using Aperture Synthesis (MIRAS) designed by ESA. MIRAS measures the brightness temperature (TB) for the recovery of SSS over the oceans and soil moisture on land. SMOS satellite follows the sun-synchronous polar orbit at an altitude of 763 km. The satellite achieves global coverage every three days with an incidence range of 0-60 • and a swath of~1000 km at full polarization.
In this study, we used two SMOS SSS products. The first is the de-biased SMOS Level 3 (L3) gridded product version 5 made by the Laboratoire d'Oceanographie et du Climat (LOCEAN) and the Ocean Salinity Expertise Center (CEC-OS) of Centre Aval de Traitements des Données SMOS (CATDS). Compared to previous versions, this product uses an improved adjustment of land-sea biases close to the coast and ice filtering at high latitudes [35]. The 9-day running maps are available on a 0.25 • × 0.25 • grid at 4-day time scale for the 2010-2019 period. This product is referred to as SMOS LOCEAN hereafter. The second SMOS product is the so-called L3 version 2 product from Barcelona Expert Center (BEC) [41]. This product was derived from the debiased non-Bayesian retrieval that reduces the spatial biases induced by land-sea and radio frequency interference contamination close to the coast. The 9-day running maps , the so-called L3 product are daily produced and are available on a 0.25 • × 0.25 • grid at daily time scale for the 2011-2019 period. This product is referred to as SMOS BEC hereafter.

SMAP Satellite Products
The National Aeronautics and Space Administration (NASA) leads the Soil Moisture Active Passive (SMAP) mission that also monitor the SSS. The SMAP satellite is in a nearpolar orbit, sun-synchronous, at an inclination of 98 • and an altitude of 685 km. With its swath (~1000 km), the satellite reaches its global coverage in approximatively 3 days with a spatial resolution of approximatively 40 km. The sensors onboard the mission are passive microwave radiometers operating at 1.4 GHz (L-band). SSS is retrieved from the TB measurement.
Here, we used two SMAP SSS maps produced by different centers. One was the SMAP L3 version 4.3 product developed by Jet Propulsion Laboratory (JPL). This product is based on the retrieval algorithm that was developed in the context of the previous Aquarius mission and updated to SMAP. It has an inherent feature resolution of 60 km and is mapped on a 0.25 • × 0.25 • grid, at an 8-day resolution (based on temporal running means). This product includes land-sea bias corrections [42]. It is referred hereafter as SMAP JPL. The other product is the SMAP L3 version 4.0 product developed by Remote Sensing system (RSS, www.remss.com/missions/smap, accessed on 25 June 2021). The product is derived from a 40-km spatial resolution product using Backus Gilbert Optimum Interpolation Remote Sens. 2021, 13, 2507 5 of 23 (OI). It was then smoothed with an approximatively 70-km spatial resolution using simple next neighbor averaging at each grid cell that allows noise reduction. It is available on a 0.25 • × 0.25 • , 8-day grid. This product includes land correction with retrieval of salinity within 30-40 km off the coast. It is referred to hereafter as SMAP RSS. These two SMAP SSS products are available from 1 April 2015 to the present and are distributed by the NASA Physical Oceanography Distributed Active Archive Center (PO.DAAC).

Reanalysis Product
The Global 1/12 • Ocean Reanalysis version 1 (GLORYS12V1), developed by Mercator and distributed by the European Copernicus Marine Environment Monitoring Service (CMEMS), aims at producing state of the art eddy resolving ocean simulations constrained by oceanic observations by means of data assimilation [43]. The ocean model component is NEMO, forced at the surface by ECMWF-ERA interim reanalysis including precipitation and radiative fluxes corrections. The assimilation scheme is 3D-VAR, providing a correction for the slowly evolving large-scale biases in temperature and salinity. The assimilated observations are Reynolds 0.25 AVHRR-only SST, delayed time SLA from the satellites, in situ T/S profiles from the Coriolis Ocean dataset for Reanalysis version 4.1 (CORAV4.1) database [44]. The ETOPO1 Earth topography at 1 arc-minute resolution was used for the deep ocean topography while the Global Earth Bathymetric Chart of the Oceans (GEBCO) 1-resolution dataset was used for the coast and continental shelf topography. The GLORYS12V1 product includes daily and monthly outputs of salinity, temperature, currents, sea level, mixed layer depth and ice parameters available on a regular grid at 1/12 • and 50 standard vertical levels. The product spans from January 1993 to June 2019. In this study, we used the daily salinity outputs estimated at the ocean surface (~0.5 m depth). This product is referred to as GLORYS hereafter.

Methods
The computation of coastal SSS gradients from TSG, satellite and reanalysis data is based on a linear least square fit method. Gradients were estimated in bands ranging from 50 to 150, 200 and 250 km from the coast, along cross-shore sections. Data in the first 50 km from the coast were excluded because of the poor skill of the satellite data in this band. Due to the spatial scales considered here, coastal data around islands smaller than 150 km in diameter were excluded from the present study.

SSS Gradient Calculation from TSG Data
SSS gradient from TSG data were calculated along cross-shore ship tracks data available in coastal regions ( Figure 2). The extractions of cross-shore ship tracks work as follows. For a given ship trip (Figure 2a), tracks are first extracted in each coastal region whatever their orientation to the coast (red ship tracks in Figure 2b). To estimate the orientation of the tracks, we use the nearest coastline (green contour in Figure 2b). As a first step, a track section is considered as cross-shore when the angle between the ship trajectory and the coastline, linearized by a least square fit in a large box enclosing the whole coastal part of the track (black thick line in Figure 2b), is in a 90 • ± 40 • range. As a second step, after excluding other parts of the track, for each individual cross-shore section, the coastline is linearized in a smaller box enclosing the section ( Figure 2c) and only sections with angle of 90 • ± 20 • relative to the coastline are kept. Sections with length less than 25 km are excluded. Finally, the SSS cross-shore gradient is estimated by linear least square fit applied on the SSS data collected along each selected cross-shore section.

SSS Gradient Calculation from Gridded Products
Reanalysis and satellite products are available on a global regular grid. Therefore, the estimation of SSS gradients is different from the TSG. First, every 25 km along the global coastline, we compute cross-shore axes, which local direction is estimated relatively to a linear least-square fit of the local coastline (within a 200 km radius of the coastal point considered, Figure 3b). Second, for each gridded product, coastal SSS data are selected ( Figure 3a) and those found in a 50-km wide swath following each cross-shore axis (blue rectangle in Figure 3b) are interpolated at points every 25 km along this axis. Third, SSS cross-shore gradients are estimated at each time step by a linear least-square fit of these interpolated data. A minimum of three data points available on each cross-shore axis is required to estimate the gradient. The number of gradients computed with this method varies from one product to another because they have different coverage in coastal regions.

SSS Gradient Calculation from Gridded Products
Reanalysis and satellite products are available on a global regular grid. Therefore, the estimation of SSS gradients is different from the TSG. First, every 25 km along the global coastline, we compute cross-shore axes, which local direction is estimated relatively to a linear least-square fit of the local coastline (within a 200 km radius of the coastal point considered, Figure 3b). Second, for each gridded product, coastal SSS data are selected ( Figure 3a) and those found in a 50-km wide swath following each cross-shore axis (blue rectangle in Figure 3b) are interpolated at points every 25 km along this axis. Third, SSS cross-shore gradients are estimated at each time step by a linear least-square fit of these interpolated data. A minimum of three data points available on each cross-shore axis is required to estimate the gradient. The number of gradients computed with this method varies from one product to another because they have different coverage in coastal regions.

SSS Gradient Calculation from Gridded Products
Reanalysis and satellite products are available on a global regular grid. Therefore, the estimation of SSS gradients is different from the TSG. First, every 25 km along the global coastline, we compute cross-shore axes, which local direction is estimated relatively to a linear least-square fit of the local coastline (within a 200 km radius of the coastal point considered, Figure 3b). Second, for each gridded product, coastal SSS data are selected ( Figure 3a) and those found in a 50-km wide swath following each cross-shore axis (blue rectangle in Figure 3b) are interpolated at points every 25 km along this axis. Third, SSS cross-shore gradients are estimated at each time step by a linear least-square fit of these interpolated data. A minimum of three data points available on each cross-shore axis is required to estimate the gradient. The number of gradients computed with this method varies from one product to another because they have different coverage in coastal regions.

Global Statistics on Coastal SSS and SSS Gradient
In addition to computing the SSS gradients, we evaluated the difference between SSS measurements from TSG and estimation from satellites and reanalysis products in the nearest 300 km from the coast, divided in 50 km width coastal bands. SSS from TSG and from gridded products were co-located in time and space for the common period from April 2015 to June 2019. Each TSG measurement was associated with SSS from a given product at the nearest grid point, and then all TSG measurements associated with the same grid point were averaged to create pairs of TSG and product SSS values. Based on these pairs, global statistics of differences between the SSS from TSG and from each gridded product can be quantified. We computed the bias and standard deviation (STD) of the differences as follows: where SSS product corresponds to the SSS from a satellite or reanalysis product, SSS TSG corresponds to the SSS from TSG measurements and N corresponds to the total number of co-located pairs (number of pixels). The same methodology was applied to investigate satellites and reanalysis performance in SSS gradient estimation. In this case, we first evaluated if the gradients were of the same sign. Then, for the gradients of the same sign, each TSG SSS gradient was associated with the nearest gradient estimated from a given product, if any available within a radius of 50 km and 4 days. Bias and STD were computed by replacing SSS with SSS gradient in the formulas above.

Comparison of SSS from TSG Measurement and Satellites and Reanalysis Products
When comparing the satellite and reanalysis SSS products with the TSG measurements, we first investigated the global statistics. Then, we focused on five regions of interest (Table 1, Figure 1b) selected based on their particularly strong or remarkable coastal SSS gradients observed with the TSG data and the different driving physical processes.

Global Statistics in Coastal SSS
The differences between the coastal TSG measurements and the satellite or reanalysis SSS products vary according to the product ( Figure 4). They are relatively small between 100 and 300 km from the coast. In this band, SMOS LOCEAN, SMOS BEC and SMAP RSS show a generally positive mean bias of 0.09 pss (0.07 to 0.14 pss), 0.08 pss (0.01 to 0.12 pss) and 0.15 pss (−0.02 to 0.25 pss), respectively. It is smaller and negative for SMAP JPL and GLORYS with values of −0.07 pss (−0.09 to −0.03 pss) and −0.02 pss (−0.03 to 0.03 pss), respectively. In terms of STD, farther than 100 km from the coast, the ranking model skill is very dependent on the 50-km wide band considered. However, SMOS BEC shows the overall largest STD, ranging from 0.76 to 1.32 pss, and GLORYS shows the smallest STD, ranging from 0.38 to 0.85 pss, with the other products in between. Except for the bias of SMAP RSS and the STD of SMOS LOCEAN, there is no obvious loss in data quality when approaching the coast farther than 100 km offshore. density of ship traffic, but the decreases in the last 50 km in all products but GLORYS are due to less available satellite data. The number of pixels also varies from one satellite product to another, reflecting more or less data availability. We could have selected only grid points where all satellite products are available for the statistical comparison with TSG data, but we adopted the point of view of the coastal data user, exploiting all the data available in the product distribution. In addition, GLORYS has one order of magnitude more pixels than the satellites products, even relatively far from the coast.

Regions of Interest
In order to validate the accuracy of the satellites and reanalysis products, we selected five regions of interest due to their strong or remarkable coastal SSS gradients observed with TSG data, representative of specific physical processes. We provide a qualitative comparison of gridded products with ship TSG sections to illustrate regional mesoscale features. The selected regions are the Amazon River plume (AMZ), northeast Brazil (NEB), the Bay of Bengal (BoB), the California coast (CAL) and the Great Australian Bight (GAB) ( Figure 1b; Table 1). The AMZ region is an area where significant surface freshening takes place due to the discharge of the Amazon River [45][46][47]. Moreover, it is a region where the North Brazil Current (NBC), a strong western boundary surface current, flows in a northwestern direction along the Brazilian coast [48]. The NBC brings salty waters, with a maximum salinity core located around 100 m depth that can be traced back to the SSS maximum in the South Atlantic subtropical gyre in boreal winter, which creates an SSS contrast with the low salinity of the Amazon river plume [49][50][51]. At about 6-8 • N, 48-50 • W, the NBC separates from the coast and partly retroflects eastward to connect with the North Equatorial Counter-Current (NECC), generating large anticyclonic rings that can reach the Caribbean Sea (e.g., [52,53]). The upper NBC transports low-salinity waters from the Amazon plume along the coast, which can later feed the NECC or swirl around the NBC rings [54,55].
In this region, all tested products reproduce the alongshore extension of the Amazon water lens at about 54 • W ( Figure 5). In addition, Amazon water lenses were also observed in all products extending along the continental shelf. GLORYS, which provides the finest resolution, shows a ring centered at 10 • N, 57 • W that is not represented on satellites maps, probably due to modelling artifact. Differences in coastal flagging from one product to another can be seen and correspond to the different number of pixels available in the 0-50 km coastal band ( Figure 4c). All products represent well the strong cross-shore salinity gradient observed in the TSG transect with SSS increasing from more than 37 pss offshore to less than 32 pss at the coast. However, the distance from the coast of the haline front (around 300 km) is underestimated in all products, which are also too salty until 48 • W along the TSG track. The different vertical sampling of salinity at 0-10 m in TSG data vs. 5 m in GLORYS vs. skin layer in satellite products [17,34], potentially important in this highly stratified region [56], cannot explain these differences. Temporal sampling (8-day running mean vs. 3-day transect) may contribute to them. Amazon River Plume The AMZ region is an area where significant surface freshening takes place due to the discharge of the Amazon River [45][46][47]. Moreover, it is a region where the North Brazil Current (NBC), a strong western boundary surface current, flows in a northwestern direction along the Brazilian coast [48]. The NBC brings salty waters, with a maximum salinity core located around 100 m depth that can be traced back to the SSS maximum in the South Atlantic subtropical gyre in boreal winter, which creates an SSS contrast with the low salinity of the Amazon river plume [49][50][51]. At about 6-8°N, 48-50°W, the NBC separates from the coast and partly retroflects eastward to connect with the North Equatorial Counter-Current (NECC), generating large anticyclonic rings that can reach the Caribbean Sea (e.g., [52,53]). The upper NBC transports low-salinity waters from the Amazon plume along the coast, which can later feed the NECC or swirl around the NBC rings [54,55].
In this region, all tested products reproduce the alongshore extension of the Amazon water lens at about 54°W ( Figure 5). In addition, Amazon water lenses were also observed in all products extending along the continental shelf. GLORYS, which provides the finest resolution, shows a ring centered at 10°N, 57°W that is not represented on satellites maps, probably due to modelling artifact. Differences in coastal flagging from one product to another can be seen and correspond to the different number of pixels available in the 0-50 km coastal band ( Figure 4c). All products represent well the strong cross-shore salinity gradient observed in the TSG transect with SSS increasing from more than 37 pss offshore to less than 32 pss at the coast. However, the distance from the coast of the haline front (around 300 km) is underestimated in all products, which are also too salty until 48°W along the TSG track. The different vertical sampling of salinity at 0-10 m in TSG data vs. 5 m in GLORYS vs. skin layer in satellite products [17,34], potentially important in this highly stratified region [56], cannot explain these differences. Temporal sampling (8-day running mean vs. 3-day transect) may contribute to them.

Northeast Brazil
Northeast Brazil, located 500-1000 km south of the Amazon River mouth, is the source region for the NBC. Between 10 • S and 5 • S, the North Brazil Undercurrent (NBUC) flows equatorward along the coast [57][58][59]. This western boundary current originates from the northward bifurcation of the westward southern South Equatorial Current (sSEC) at the Brazilian coast further south. As the coast orientation changes at around 5 • S, it is joined by the central South Equatorial Current (cSEC) to form the NBC. The NBUC has a weak surface expression due to the southward wind-driven surface drift in this region [60]. It transports salty waters from subtropical origin with maximum salinity found around 100 m depth close to the coast [57].
It is obvious from the SSS maps that the NBUC salty core also has a surface signature [61], with SSS increasing over 36.5 pss toward the coast ( Figure 6). This may be favored by vertical mixing due to the wind-induced shear stress or friction with the shelf slope. SSS is particularly high over the narrow shelf, possibly due to evaporation, as seen in TSG and GLORYS data that extend closer to the coast than satellite data. A direct comparison with TSG data is difficult for SMAP RSS between 6 • S and 8 • S due to the coastal flagging, while SMOS BEC seems particularly noisy. However, all products represent the high salinity extending farther offshore south of 7.5 • S well, although it is too high in the GLORYS and SMAP products. This pattern is probably associated with the orographic effect here that shifted the path of the NBUC away from the coast [58]. From the coast to the islands located around 4 • S, GLORYS reproduces the TSG observations better than the satellite products. Northeast Brazil, located 500-1000 km south of the Amazon River mouth, is the source region for the NBC. Between 10°S and 5°S, the North Brazil Undercurrent (NBUC) flows equatorward along the coast [57][58][59]. This western boundary current originates from the northward bifurcation of the westward southern South Equatorial Current (sSEC) at the Brazilian coast further south. As the coast orientation changes at around 5°S, it is joined by the central South Equatorial Current (cSEC) to form the NBC. The NBUC has a weak surface expression due to the southward wind-driven surface drift in this region [60]. It transports salty waters from subtropical origin with maximum salinity found around 100 m depth close to the coast [57].
It is obvious from the SSS maps that the NBUC salty core also has a surface signature [61], with SSS increasing over 36.5 pss toward the coast (Figure 6). This may be favored by vertical mixing due to the wind-induced shear stress or friction with the shelf slope. SSS is particularly high over the narrow shelf, possibly due to evaporation, as seen in TSG and GLORYS data that extend closer to the coast than satellite data. A direct comparison with TSG data is difficult for SMAP RSS between 6°S and 8°S due to the coastal flagging, while SMOS BEC seems particularly noisy. Error! Bookmark not defined.However, all products represent the high salinity extending farther offshore south of 7.5°S well, although it is too high in the GLORYS and SMAP products. This pattern is probably associated with the orographic effect here that shifted the path of the NBUC away from the coast [58]. From the coast to the islands located around 4°S, GLORYS reproduces the TSG observations better than the satellite products.

Bay of Bengal
The Bay of Bengal (BoB) is a region of high SSS contrasts due to strong seasonal freshwater inputs by summer monsoon rains, which feed large rivers. Oceanic rainfall and river discharge contribute equally to freshwater inputs in the northern BoB, while the Ganga-Brahmaputra River located at its the northern end represents more than 60% of the total discharge in the BoB [15,[62][63][64]. After monsoons, the East Indian Coastal Current (EICC) reverses and flows southward, transporting this freshwater along what has been

Bay of Bengal
The Bay of Bengal (BoB) is a region of high SSS contrasts due to strong seasonal freshwater inputs by summer monsoon rains, which feed large rivers. Oceanic rainfall and river discharge contribute equally to freshwater inputs in the northern BoB, while the Ganga-Brahmaputra River located at its the northern end represents more than 60% of the total discharge in the BoB [15,[62][63][64]. After monsoons, the East Indian Coastal Current (EICC) reverses and flows southward, transporting this freshwater along what has been called a River in the Sea (RIS), which is about 100 km wide and can extend up to 3000 km from the main river mouth [15,65,66].
Some satellite SSS products have been tested in the BoB, with available in situ data. The semi-enclosed configuration of the BoB induces strong land-sea contamination and early SMOS products had very limited skill in this region [67], but improved versions perform much better and allow the study of interannual SSS variations thanks to the decade-long record [68]. SMAP would even work better. It compares well with in situ SSS cross-shore sections in the BoB and reproduces a large part of the seasonal SSS variations measured at coastal stations along the RIS [18,66,69]. In October 2015 (Figure 7), SMOS LOCEAN, SMAP JPL and GLORYS captured the coastal freshening (SSS < 32 pss) seen by the TSG data at 13 • N, clearly associated with the northern BoB freshening transported southward by the RIS. SMOS BEC and SMAP RSS lack coastal data to enable the comparison. The RIS southward extension is impacted by interannual and intraseasonal anomalies [69]. It was relatively limited in 2015 due to a negative Indian Ocean Dipole and, according to GLORYS, by a mesoscale cyclonic eddy that diverts its path [66][67][68][69].
called a River in the Sea (RIS), which is about 100 km wide and can extend up to 3000 km from the main river mouth [15,65,66].
Some satellite SSS products have been tested in the BoB, with available in situ data. The semi-enclosed configuration of the BoB induces strong land-sea contamination and early SMOS products had very limited skill in this region [67], but improved versions perform much better and allow the study of interannual SSS variations thanks to the decade-long record [68]. SMAP would even work better. It compares well with in situ SSS cross-shore sections in the BoB and reproduces a large part of the seasonal SSS variations measured at coastal stations along the RIS [18,66,69]. In October 2015 (Figure 7), SMOS LOCEAN, SMAP JPL and GLORYS captured the coastal freshening (SSS < 32 pss) seen by the TSG data at 13°N, clearly associated with the northern BoB freshening transported southward by the RIS. SMOS BEC and SMAP RSS lack coastal data to enable the comparison. The RIS southward extension is impacted by interannual and intraseasonal anomalies [69]. It was relatively limited in 2015 due to a negative Indian Ocean Dipole and, according to GLORYS, by a mesoscale cyclonic eddy that diverts its path [66][67][68][69].

California Upwelling System
One of the major Eastern Boundary Upwelling Systems (EBUS) takes place along the California coast. This wind-driven upwelling generally uplifts salty subsurface waters to the surface, while the California Current, farther offshore, transports relative fresh waters equatorward, which creates a cross-shore SSS gradient [70][71][72][73].
The upwelling is particularly strong along the central California coast in summer ( Figure 8). TSG data show an SSS contrast of about 1 pss with an SSS over 33 pss in a coastal band that gets wider (around 100 km) towards the south. GLORYS reproduces this cross-shore SSS gradient remarkably well. Among satellite products, only SMOS BEC reproduces this gradient, except for a few grid points around 36°N. It is also the satellite product with the best coastal coverage in this region. Other satellite products do not capture the coastal SSS increase. On the contrary, SMAP JPL shows a strong contrast between a relatively high SSS (>34) offshore and a low SSS (<32.5) at the coast, south of 36°N. This product, however, reproduced the in situ SSS gradients much better at monthly to yearly timescales in this region [74]. For the SMAP RSS, which has a data gap at the time of the TSG observations (17-25 June 2019), we instead show data from 14-15 June 2019, which may partly explain the differences.
Note that, while some upwelling regions are associated with coastal SSS increase [75,76], others are associated with coastal freshening (e.g., [77]). The mean vertical salinity

California Upwelling System
One of the major Eastern Boundary Upwelling Systems (EBUS) takes place along the California coast. This wind-driven upwelling generally uplifts salty subsurface waters to the surface, while the California Current, farther offshore, transports relative fresh waters equatorward, which creates a cross-shore SSS gradient [70][71][72][73].
The upwelling is particularly strong along the central California coast in summer ( Figure 8). TSG data show an SSS contrast of about 1 pss with an SSS over 33 pss in a coastal band that gets wider (around 100 km) towards the south. GLORYS reproduces this cross-shore SSS gradient remarkably well. Among satellite products, only SMOS BEC reproduces this gradient, except for a few grid points around 36 • N. It is also the satellite product with the best coastal coverage in this region. Other satellite products do not capture the coastal SSS increase. On the contrary, SMAP JPL shows a strong contrast between a relatively high SSS (>34) offshore and a low SSS (<32.5) at the coast, south of 36 • N. This product, however, reproduced the in situ SSS gradients much better at monthly to yearly timescales in this region [74]. For the SMAP RSS, which has a data gap at the time of the TSG observations (17-25 June 2019), we instead show data from 14-15 June 2019, which may partly explain the differences.
gradients are indeed different and can be of opposite sign in the four major EBUS [78]. This depends on the upwelled subsurface water masses, that can vary locally within the same upwelling system, at seasonal or interannual timescales, and can be influenced by other processes such as river plumes [79][80][81].

Great Australian Bight
In the Great Australian Bight (GAB), the Spencer Gulf is a remarkable inverse estuary, where strong evaporation over shallowing waters produces a year-round meridional SSS gradient, with maximum SSS at its head [82]. This phenomenon occurs in other hot and dry regions of the world [83][84][85]. The salinity increase is particularly strong during summer and is followed by winter cooling. This results in a density increase that generates gravity currents out of the gulf mouth and along the shelf slope, compensated by a surface inflow of fresher waters [86,87]. It was later found that this process occurs at larger scale in the whole GAB that has been called a mega inverse estuary [87,88].
The available TSG observations do not penetrate far inside the GAB but clearly show that the SSS increases toward the entrance of the Spencer Gulf, where it already reaches 36 pss compared to values lower than 35.5 pss farther offshore (Figure 9). SMAP JPL is too fresh compared to the TSG data south of 36°S and shows a large scale SSS increase in the GAB, less pronounced than in other products. Moreover, the SMAP JPL coastal freshening is unexpected in this dry region with very little freshwater runoff. SMOS LOCEAN is generally too salty north of 36°S while the few grid points in the Spencer Gulf do not capture the expected local SSS maximum. On the contrary, GLORYS, SMOS BEC and SMAP RSS reproduce the SSS along the ship track well. The SSS increase in the Spencer Gulf captured by GLORYS is consistent with the few grid points available here for the latter satellite products. Except for probably spurious local freshening at a few points along the coast in SMAP RSS, SMOS BEC and SMAP RSS largely agree with GLORYS on the large-scale meridional gradient between 40°S and 32°S. Note that, while some upwelling regions are associated with coastal SSS increase [75,76], others are associated with coastal freshening (e.g., [77]). The mean vertical salinity gradients are indeed different and can be of opposite sign in the four major EBUS [78]. This depends on the upwelled subsurface water masses, that can vary locally within the same upwelling system, at seasonal or interannual timescales, and can be influenced by other processes such as river plumes [79][80][81].

Great Australian Bight
In the Great Australian Bight (GAB), the Spencer Gulf is a remarkable inverse estuary, where strong evaporation over shallowing waters produces a year-round meridional SSS gradient, with maximum SSS at its head [82]. This phenomenon occurs in other hot and dry regions of the world [83][84][85]. The salinity increase is particularly strong during summer and is followed by winter cooling. This results in a density increase that generates gravity currents out of the gulf mouth and along the shelf slope, compensated by a surface inflow of fresher waters [86,87]. It was later found that this process occurs at larger scale in the whole GAB that has been called a mega inverse estuary [87,88].
The available TSG observations do not penetrate far inside the GAB but clearly show that the SSS increases toward the entrance of the Spencer Gulf, where it already reaches 36 pss compared to values lower than 35.5 pss farther offshore (Figure 9). SMAP JPL is too fresh compared to the TSG data south of 36 • S and shows a large scale SSS increase in the GAB, less pronounced than in other products. Moreover, the SMAP JPL coastal freshening is unexpected in this dry region with very little freshwater runoff. SMOS LOCEAN is generally too salty north of 36 • S while the few grid points in the Spencer Gulf do not capture the expected local SSS maximum. On the contrary, GLORYS, SMOS BEC and SMAP RSS reproduce the SSS along the ship track well. The SSS increase in the Spencer Gulf captured by GLORYS is consistent with the few grid points available here for the latter satellite products. Except for probably spurious local freshening at a few points along the coast in SMAP RSS, SMOS BEC and SMAP RSS largely agree with GLORYS on the large-scale meridional gradient between 40 • S and 32 • S.

Comparison of SSS Gradient from TSG Measurement and Satellites and Reanalysis Products
In this section, we examine the satellites and the reanalysis's skill at reproducing the observed SSS gradients computed from TSG. We first present the global distribution of the SSS gradients computed from the TSG, satellites and reanalysis. Then, we present the global statistics on the SSS gradients comparison between gridded products and TSG. In the regions of interest, as it is also the case globally (Figure 4c), the satellite products have limited coverage near the coast, while reanalysis products are available up to the coastline. We therefore only estimated the gradients at distances higher than 50 km from the coast.

Global Distribution of SSS Gradients
Collocating the available TSG SSS gradients computed in the 50-250 km coastal band with those of other products, as described in Section 2.2, allows the most direct comparison. Indeed, each TSG SSS gradient is compared with the nearest gradient estimated from a given product within a radius of 50 km and 4 days. However, if several TSG SSS gradients are available within the collocation radius, their average is used for the comparison. As observed in Figure 4c for coastal SSS, the number of cross-shore sections used to compare with TSG varies from one product to another. All products agree with TSG on SSS gradients with typical order of ±0.1 pss over 100 km (Figure 10), with a distribution skewed toward negative values. However, only SMAP JPL and GLORYS reproduce the observed mean SSS gradient. Other products overestimate the observed SSS gradient. In terms of the STD of the SSS gradient, all products show lower values than the TSG observations, although SMOS BEC shows weaker difference with the TSG 0.7 pss over 100 km. Other products exhibited differences higher than 1.1 pss over 100 km.
The spatial distribution of the collocated gradients shows that all products reproduce to some extent the sign of the gradients (Figure 11). In northeast Brazil, all products capture the positive SSS gradient observed by TSG (see also Figure 6). In the Amazon River plume region, all products but SMAP RSS also reproduce the negative gradient due to the river discharge, as seen in Figure 5. In the Bay of Bengal, the in situ TSG SSS gradient is positive, which seems contradictory with the negative gradient seen in Figure 7. This discrepancy may be because the gradient estimated from TSG corresponds to the mean gradient of several sections and the one corresponding to the case shown in Figure 7 may be flooded in. The seasonal variation of the extension of the RIS due to the mesoscale activity and the negative Indian Ocean dipole may also favor this discrepancy.

Comparison of SSS Gradient from TSG Measurement and Satellites and Reanalysis Products
In this section, we examine the satellites and the reanalysis's skill at reproducing the observed SSS gradients computed from TSG. We first present the global distribution of the SSS gradients computed from the TSG, satellites and reanalysis. Then, we present the global statistics on the SSS gradients comparison between gridded products and TSG. In the regions of interest, as it is also the case globally (Figure 4c), the satellite products have limited coverage near the coast, while reanalysis products are available up to the coastline. We therefore only estimated the gradients at distances higher than 50 km from the coast.

Global Distribution of SSS Gradients
Collocating the available TSG SSS gradients computed in the 50-250 km coastal band with those of other products, as described in Section 2.2, allows the most direct comparison. Indeed, each TSG SSS gradient is compared with the nearest gradient estimated from a given product within a radius of 50 km and 4 days. However, if several TSG SSS gradients are available within the collocation radius, their average is used for the comparison. As observed in Figure 4c for coastal SSS, the number of cross-shore sections used to compare with TSG varies from one product to another. All products agree with TSG on SSS gradients with typical order of ±0.1 pss over 100 km (Figure 10), with a distribution skewed toward negative values. However, only SMAP JPL and GLORYS reproduce the observed mean SSS gradient. Other products overestimate the observed SSS gradient. In terms of the STD of the SSS gradient, all products show lower values than the TSG observations, although SMOS BEC shows weaker difference with the TSG 0.7 pss over 100 km. Other products exhibited differences higher than 1.1 pss over 100 km.
The spatial distribution of the collocated gradients shows that all products reproduce to some extent the sign of the gradients (Figure 11). In northeast Brazil, all products capture the positive SSS gradient observed by TSG (see also Figure 6). In the Amazon River plume region, all products but SMAP RSS also reproduce the negative gradient due to the river discharge, as seen in Figure 5. In the Bay of Bengal, the in situ TSG SSS gradient is positive, which seems contradictory with the negative gradient seen in Figure 7. This discrepancy may be because the gradient estimated from TSG corresponds to the mean gradient of several sections and the one corresponding to the case shown in Figure 7 may be flooded in. The seasonal variation of the extension of the RIS due to the mesoscale activity and the negative Indian Ocean dipole may also favor this discrepancy.  After collocating the SSS gradients from the TSG and gridded products, we now analyze the distribution of all SSS gradients available from TSG, to take advantage of the long period of observations, and all the SSS gradients available from the gridded products on their common period, to take advantage of their global homogeneous coverage. Note  After collocating the SSS gradients from the TSG and gridded products, we now analyze the distribution of all SSS gradients available from TSG, to take advantage of the long period of observations, and all the SSS gradients available from the gridded products on their common period, to take advantage of their global homogeneous coverage. Note After collocating the SSS gradients from the TSG and gridded products, we now analyze the distribution of all SSS gradients available from TSG, to take advantage of the long period of observations, and all the SSS gradients available from the gridded products on their common period, to take advantage of their global homogeneous coverage. Note that the SSS gradients shown for satellites and reanalysis products correspond to the temporal average of SSS gradients computed every 4 days over the satellites' common period. The distribution skewness toward negative cross-shore gradients is confirmed for TSG and all products ( Figure 12). Although the TSG data are relatively sparse in time and space, they display the largest range of gradients, almost reproduced at global scale by GLORYS, which shows a slightly smaller STD. However, the gradient STD is around three times smaller for the SMAP products, and five times smaller for the SMOS products, probably due to the lower resolution of satellite products. Among the gridded products, GLORYS and SMAP show the most negative mean gradients (−0.35 pss/100 km) while it is less than half for SMOS products (−0.15 pss/100 km). space, they display the largest range of gradients, almost reproduced at global scale by GLORYS, which shows a slightly smaller STD. However, the gradient STD is around three times smaller for the SMAP products, and five times smaller for the SMOS products, probably due to the lower resolution of satellite products. Among the gridded products, GLORYS and SMAP show the most negative mean gradients (−0.35 pss/100 km) while it is less than half for SMOS products (−0.15 pss/100 km).
At global scale, the geographical distribution of SSS gradients reveals a mosaic of negative and positive gradients from all products ( Figure 13). TSG data show a relatively balanced occurrence of positive and negative gradients. On the contrary, satellite products are clearly dominated by negative gradients. The SMOS LOCEAN product shows more positive gradients than other satellite products, although they are not found in the regions where we highlighted such occurrence. Compared to the satellite products, GLORYS shows more regions of well-pronounced positive gradients, including the Great Australian Bight, northeast Brazil and California upwelling regions. Although GLORYS seems to perform better than satellites products on these maps, global statistics based on comparison with collocated TSG data are required to objectively and quantitatively evaluate products performance in the estimation of SSS coastal gradients.  At global scale, the geographical distribution of SSS gradients reveals a mosaic of negative and positive gradients from all products ( Figure 13). TSG data show a relatively balanced occurrence of positive and negative gradients. On the contrary, satellite products are clearly dominated by negative gradients. The SMOS LOCEAN product shows more positive gradients than other satellite products, although they are not found in the regions where we highlighted such occurrence. Compared to the satellite products, GLORYS shows more regions of well-pronounced positive gradients, including the Great Australian Bight, northeast Brazil and California upwelling regions. Although GLORYS seems to perform better than satellites products on these maps, global statistics based on comparison with collocated TSG data are required to objectively and quantitatively evaluate products performance in the estimation of SSS coastal gradients.

Global Statistics of SSS Gradients Comparison
Differences between the TSG measurements and the satellite or reanalysis SSS products in terms of cross-shore SSS gradient vary according to the product ( Figure 14). All products failed in estimating the sign of the gradient in a proportion higher to 30%, whatever the distance from the shelf (Figure 14a). SMOC BEC is the product with the highest error (ranging from 43 to 48%) while SMAC JPL has the lowest error (ranging from 30 to 32%).
When looking at the differences considering only the gradients of the same sign (Figure 14b,c), in the three selected coastal bands, SMOS LOCEAN, SMAP JPL, SMAP RSS and GLORYS overestimate the observed gradients (positive bias) while SMOS BEC underestimates them (negative bias). Overall, all products except GLORYS show a decreasing bias (in absolute value) from the largest (50-250 km) to the narrowest (50-150 km) coastal bands. The SMOS products show the highest STD whatever the coastal band. A minimum criterion on STD rather than bias is more relevant to identify the best product when STD is one order larger than the bias, as is the case here. GLORYS shows an increasing bias and STD from the largest to the narrowest coastal bands, which are on average between those found for the SMAP products. The SMAP products show on average lower bias (in absolute value) and STD than the SMOS products. It must be noted that the number of cross-shore sections also varies according to the product, especially in the 50-150 km band. These variations are consistent with the different numbers of SSS pixels in the coastal zone ( Figure 4).

Global Statistics of SSS Gradients Comparison
Differences between the TSG measurements and the satellite or reanalysis SSS products in terms of cross-shore SSS gradient vary according to the product ( Figure 14). All products failed in estimating the sign of the gradient in a proportion higher to 30%, whatever the distance from the shelf (Figure 14a). SMOC BEC is the product with the highest error (ranging from 43 to 48%) while SMAC JPL has the lowest error (ranging from 30 to 32%).

Discussion
In the present study, we investigated the skill of different satellite and reanalysis products at reproducing SSS in the global coastal ocean (<300 km), as observed with TSG measurements. Four gridded SMOS or SMAP satellite products and a model reanalysis were compared with TSG measurements in terms of coastal SSS and cross-shore SSS When looking at the differences considering only the gradients of the same sign (Figure 14b,c), in the three selected coastal bands, SMOS LOCEAN, SMAP JPL, SMAP RSS and GLORYS overestimate the observed gradients (positive bias) while SMOS BEC underestimates them (negative bias). Overall, all products except GLORYS show a decreasing bias (in absolute value) from the largest (50-250 km) to the narrowest (50-150 km) coastal bands. The SMOS products show the highest STD whatever the coastal band. A minimum criterion on STD rather than bias is more relevant to identify the best product when STD is one order larger than the bias, as is the case here. GLORYS shows an increasing bias and STD from the largest to the narrowest coastal bands, which are on average between those found for the SMAP products. The SMAP products show on average lower bias (in absolute value) and STD than the SMOS products. It must be noted that the number of cross-shore sections also varies according to the product, especially in the 50-150 km band. These variations are consistent with the different numbers of SSS pixels in the coastal zone ( Figure 4).

Discussion
In the present study, we investigated the skill of different satellite and reanalysis products at reproducing SSS in the global coastal ocean (<300 km), as observed with TSG measurements. Four gridded SMOS or SMAP satellite products and a model reanalysis were compared with TSG measurements in terms of coastal SSS and cross-shore SSS gradient.
Our study evaluated gridded (Level 3-L3) satellite products rather than along-track (Level 2-L2) satellite products because of their wider use by the scientific community. In addition, comparing L2 products with TSG measurements would greatly reduce the number of collocated data, which would limit the analyses on a global scale.
Through the global statistics on SSS and SSS gradients, our results revealed that the four SMOS and SMAP L3 gridded products deviated from each other and from the TSG measurement in the first 300 km from the coast. Several factors can explain the discrepancy between satellites and TSG data. First, spatiotemporal scales are different. At a typical merchant-ship speed of 15 knots, the TSG data resolution is about 2-3 km [17], and a 100 km coastal section is sailed along in about 3 h. The satellite gridded products have a typical 25-km and 8-day resolution [33,36]. The difference in spatial scales is not so much a concern in the cross-shore direction when averaging SSS in 50-km bands and estimating linear SSS by least square fit was applied in 100-km bands. However, differences in the temporal and along-shore spatial sampling mean that mesoscale and sub-mesoscale processes, which are particularly strong in coastal regions [34,89], can affect the comparison between satellite and TSG data. Second, there are also differences in vertical levels as satellites estimate salinity in the upper few centimeters, while TSG data are acquired at 5-10 m depth [17,34]. However, given the dominance of coastal freshening (Figure 11), probably due to river discharges or coastal rains, which are associated with strong stratification [54,[90][91][92], the different depth of acquisition would induce a negative bias in satellite SSS at the coast. The fact that we observed the opposite (Figure 4) suggests a small influence of the depth difference. The increasingly positive bias found in most satellite products when approaching the coast may be associated with Radio frequency interference (RFI) pollution. For instance, SMOS is mostly disturbed by RFI over coastal areas of Europe, southern Asia and the Middle East [93]. It is also consistent with an earlier study for the SMOS LOCEAN product at global scale [35]. Compared to both, the bias and standard deviation are larger in our study due to less spatiotemporal filtering of data. However, the increase in STD when approaching the coast may be favored by the increase of SSS variability as observed by [35].
Nevertheless, the bias in the coastal SSS absolute value is largely independent from the bias in reproducing coastal SSS gradients. For instance, the SMAP RSS bias in SSS increases of~1 pss toward the coast in the 50-150 km band (Figure 4a). One could expect an induced positive bias in SSS gradient of 0.02 pss km −1 (1 pss over 50 km) in this band, whereas it is only 0.002 pss km −1 (Figure 14b). With the same reasoning and still in the 50-150 km band, one could expect an induced bias in SSS gradient that would be negative for SMOS LOCEAN and positive for SMOS BEC, while the opposite is found. The fact that there is no linear relation between SSS bias and SSS gradient bias is probably due to the standard deviation difference (STD) in coastal SSS that is much larger than the bias (Figure 4). Moreover, the cross-shore SSS gradients arguably provide a better monitoring of coastal processes than coastal SSS.
To summarize, the satellite SMOS and SMAP gridded products still need to be improved in coastal regions, particularly in the nearest 100 km from the coast for SSS retrievals. Based on the comparison with collocated TSG data in this coastal band (Figure 4), The SMAP JPL product performs better with a clear minimum in bias and STD, followed by the SMOS LOCEAN, GLORYS reanalysis, SMOS BEC and SMAP RSS products. More SMAP JPL data are also available near the coast. There are fewer differences between products farther than 100 km offshore.
All products failed in estimating the sign of the gradient in more than 30% of cases. SMOC BEC performs while SMAP products perform better (Figure 14a). When considering only the cross-shore SSS gradients of the same sign, from 250, 200 or 150 km to 50 km from the coast (Figure 14b-d), SMAP JPL performs better, followed by SMAP RSS, then (with STD as the main criteria) GLORYS, SMOS LOCEAN and SMOS BEC. This seems to reveal a trade-off between data availability and quality. Indeed, SMAP RSS has the lowest number of cross-shore sections collocated with TSG data, while SMOS BEC has the largest. In this respect, the SMAP JPL product shows the best compromise. It is interesting to note that most satellite products tend to show a positive bias corresponding to an underestimation of coastal freshening. The fact that two satellite products are better than GLORYS means that they globally provide more information in the near-coastal band than a state-of-the-art model fed with Argo profiles, satellite SST and SLA observations. Despite their differences, all evaluated products agree with the dominance of coastal freshening at a global scale (Figures 11 and 13). This was suggested by ship TSG sections (Figures 12 and 14), which are relatively scarce and could be biased due to the frequent location of harbors in estuaries. The satellites and reanalysis products extend this important result globally, thanks to their homogeneous coverage. This can be attributed to river discharges as illustrated here in the Amazon River plume and the Bay of Bengal regions (Figures 6 and 8) and previously documented in regional studies [15,46,47,69,[94][95][96]. More recently, SMOS and SMAP have been shown to reproduce well seasonal and interannual SSS variation in the mouths of major rivers [97]. Strong coastal precipitations in the tropics could also induce such freshening [13,45] but their potential contribution could not be quantified here. Indeed, splitting the respective role of river runoff and precipitation on coastal freshwater plumes requires dedicated regional studies (e.g., [91]).
The SSS increase at the coast also occurs, as illustrated here in the regions of northeast Brazil, the California upwelling system and the Great Australian Bight (Figures 7, 9 and 10), and can be induced by different processes. In northeast Brazil, the increase in coastal salinity is induced by the rise of the salty core of the NBUC that flows along the coast [58,59,61]. In the upwelling system of California, the increase in coastal salinity can be related to the coastal upwelling, which brings higher-salinity waters to the surface [71,72,81]. Finally, the Great Australian Bight, known as inverse estuary, is a region where dominant evaporation promotes high salinity along the coast [87,98,99] In conclusion, the present study shows a great potential of SSS remote sensing to monitor coastal processes, somewhat limited with existing SMOS and SMAP sensors and associated products. The future Surface water and ocean Topography (SWOT) satellite mission will give access to mesoscale geostrophic currents close to the coast [100], which will allow to better understand how coastal SSS anomalies are advected in areas of strong SSS gradient such as river plumes. Ship-drift currents estimated from Automated Information Systems (AIS) data [101] will also be useful for that in the coastal zone. Nevertheless, a more direct and global access to coastal SSS dynamics would require a jump in the resolution of future SSS satellite missions [102].