Surface Currents Derived from SAR Doppler Processing: An Analysis over the Naples Coastal Region in South Italy

Abstract

Several studies have shown the capabilities of Synthetic Aperture Radar to map sea currents in ocean regions mainly characterized by large flows. We consider the well known method based on the analysis of the Doppler Centroid. The Doppler, as, in general, the scattering from the sea, is sensitive to several phenomena, occurring between the upper ocean and atmospheric boundary layers. To investigate such phenomena, we considered the combined use of both ENVISAT calibrated amplitude and Doppler data in conjunction with hindcast wind information provided by atmospheric models as well as Wind and Doppler Geophysical Model Functions (W/D-GMF) developed, in the literature, for C-Band systems. This integrated analysis for the interpretation of the Doppler surface currents measurements was carried out on a case study located in the Mediterranean Sea which is characterized by a general low circulation regime: specifically, the coastal region around the city of Naples. In this case study, we show that generally, wind plays a direct significant role in the observed Doppler surface current. The availability of an oceanographic numerical model for one of the analyzed cases also allowed us to attempt to interpret the effect of the typical thermohaline circulation pattern on the Doppler anomaly.

1. Introduction

Driven by the development of geophysical applications, investments in the construction and launch of satellites equipped with high-resolution imaging radars, specifically Synthetic Aperture Radars (SAR), have been multiplied in recent years. In addition to operating at different frequencies, SAR sensors are characterized by the availability of different operating, modes thus giving the possibility to trade-off the resolution and coverage. On one hand, very high-resolution—up to a sub-metric level at X-Band (e.g., COSMO/Skymed and TerraSAR-X/Tandem-X)—data are available with a limited coverage of a few tens of kilometers. On the other hand, by exploiting scanning beams with electronic steering of the antenna, high-resolution modes with decametric resolution at C-Band allow to achieve coverage approaching (ERS and ENVISAT) and even far exceeding (Sentinel-1) hundreds of kilometers. Above all, the increase in the number of available SAR sensors is significantly reducing the access area and the revisiting time, thus expanding the SAR data archives that frequently can go back up to 20 years for many regions of the Earth’s globe [1,2,3].

Sea surface scattering shows large variation over repeated orbits, depending on the sea state conditions. When compared to the return from terrestrial regions, the signal scattered back from the sea shows a larger (relative) dynamic, thus providing a potentially valuable source of information for understanding several geophysical processes. However, the observed scenario is, in this case, rapidly changing and influenced, much more than on the land, by the presence of processes occurring in the atmosphere, thus making the interpretation of the signal more complex.

SAR image formation mainly depends on the small gravity waves, also referred to as ripples: the occurrence of Bragg conditions coupled with the wavelength of the illuminating (microwave) radiation offers favorable conditions for an increase of backscattering. Ripples are generated by the wind and are modulated by longer waves thus providing the possibility to analyze the images to derive information on both the wind and the sea surface waves. [4]. The large sensitivity of the sea scattering makes it possible for SAR images to distinguish oceanographic features related to wind waves and swell as well as ocean internal waves and sea currents [5,6,7,8,9]. These phenomena have different spatial scales and spatial layouts and often, are present together in a single SAR image covering swaths of several tens of kilometers.

Modulation of the radar backscatter due to the wind reflects the spatial properties of the wind itself, which depends on factors such as wind speed, air stability conditions, and interaction with orography in coastal regions. The extraction of the spatial properties of the wind from SAR images has been extensively studied [10,11] and well demonstrated.

Another “hot-point” in the area of extraction of information from sea SAR images, which is the main interest of this work, is related to sea currents [12,13]. The process of hydrodynamic modulation of the sea-surface roughness allows, for example, sea surface currents signatures to be observable in SAR images, as well established in the literature, see, for instance [4,12,13,14]. The extraction of quantitative parameters on currents from the amplitude image alone is, however, a rather complex issue. Many variables are involved for the backscattering determination, such as wind vector, bathymetry, etc. [15]. For this reason, similarly to other SAR applications, the phase information may provide valuable information [4,12].

The possibility to observe the surface ocean currents by SAR images dates back to the eighties, when scientists associated some of the large scale features imaged by SAR with mesoscale oceanic eddies [16] and ocean current fronts [17].

Since then, a quantitative measurement of the radial component of the sea surface currents can be achieved from a single SAR image using spectral analysis. Specifically, it can be obtained via the estimation of the residual Doppler shifts of the received radar echoes, the Doppler Centroid Anomaly (DCA) approach [4], or by means of interferometric acquisitions with two antennas displaced in the along-track direction, Along-Track Interferometry (ATI) approach [13,18]. Both approaches have been attempted to get the radial ocean current, especially over regions of well known strong currents, as in the Gulf Stream [12,19] offshore of the USA and in the Agulhas Greater Current offshore of South Africa [20]. In [4], other experiments have been presented in coastal areas of the Gulf of Normandy, which is characterized by intense tidal currents.

Based on many successful experiments, steps aimed at providing operational services for coastal ocean monitoring have been proposed in recent years [21] and in some cases, implemented. However, regarding currents, the existing examples of applications are mainly limited to the above-cited sites characterized by significant water circulation phenomena.

This work focuses on the DCA, which, similarly to the phase offset for ATI, is induced by the movement of the scattering centers, thus providing a direct measurement of the radial component of the sea surface currents [4].

In particular, we applied the technique of retrieving the radial surface sea current from the DCA in the Mediterranean Sea, a semi-enclosed sea where the sea currents are weak and the geophysical phenomena in coastal areas are strong, due to the interactions between the airflows and the orography. In this respect, this work acts as a feasibility study.

At present, the topic of deriving sea surface currents from SAR involves many issues, both from the engineering and geophysical points of view. The former because of the typically large non-geophysical velocities associated with to the satellite attitude and the SAR instrument imaging characteristics; the latter because of the difficulty to split the different phenomena contributing to the DCA, such as the wind, the waves and the currents [22]. For the investigation of the above open issues, we considered the DCA method and we focused on a case study involving the Naples bay and the adjacent coastal region. Among the Italian coastal zones, the Gulf of Naples and the surrounding coastal regions are particularly interesting as first of all they are areas influenced by numerous environmental, socio-economic, and interacting cultural factors (strong anthropogenic impact, intense maritime traffic, the presence of polluted rivers, tourist and economic activities and the presence of four protected areas) [23]. The area is also interesting from an oceanographic standpoint as both coastal wind-driven and open-sea currents influence the local hydrodynamics [24,25,26]. In particular, several factors make this an area of general oceanographic interest, among them: (i) the presence of the Gulf of Naples, characterized by dimensions comparable to the internal Rossby deformation radius; (ii) the presence of extended shallow water areas at the edge of the continental slope; (iii) a coastal circulation shaped by both the interaction between the larger scale flows with the coastal topography and by the local wind forcing [27].

The peculiar dynamics of water circulation in this area qualify this region as a natural challenging laboratory both for investigating oceanographic phenomenology [27] and for testing the skill of ocean and atmospheric observational platforms as well as numerical models [23]. Furthermore, the choice of the zone is driven by the availability of a large historical data archive regularly acquired over the years. As also better specified in the following sections, we focus in particular on the use of the ENVISAT SAR archive due to the available larger azimuth bandwidth compared to the most recent data gathered by the Sentinel-1 sensor.

Wind fields have been estimated from calibrated SAR data by using the empirical geophysical model function CMOD developed at Institut Francais de Recherche pour l’Exploitation de la Mer (IFREMER), which relates the SAR backscatter to the actual wind vector. SAR wind inversion used the hindcast information provided by the WRF-ARW atmospheric model as reference wind vector. Moreover, the SAR estimated velocity fields were compared with the Doppler model developed by [28].

The paper is organized as described in the following section. The detailed description of the technique for measuring sea radial velocity based on the DCA estimation is presented in Section 2. In Section 3, an overview of the Bayesian approach used to extract the wind information from the SAR data is provided. In Section 4, we describe the Doppler modeling of the wind component. Section 5 is devoted to the description of the overall data analysis method. Section 6 provides details on the analyzed area and motivates the dataset choice. Section 7 is devoted to the description of the experimental results over all the different considered dates. Conclusions are finally addressed in Section 8.

2. Doppler Anomaly Extraction

The DCA technique, described in [4], is based on the measurement of the Doppler shifts induced on the SAR image by the line-of-sight component of the sea surface velocity. In our case, the measurement of DCA has been carried out by using the method based on the estimation of the azimuth auto-correlation function proposed in [29]. In particular, to avoid any bias induced by possible filtering of spectral components in the application of the SAR transfer function at the stage of generation of the Single Look Complex products, the method is applied starting from the raw data.

The block diagram that describes the overall processing chain is shown in Figure 1. More specifically, starting from raw data and after the SAR focusing operation, the corresponding Doppler Centroid frequency, referred to in the following as $f_{DC}$, is estimated over small sliding patches: details on the patch sized choice are provided in Section 7. Subsequently, the DCA, referred to as $f_{DCA}$, is extracted by subtracting the Doppler centroid induced by the antenna pointing and Earth rotation. The DCA and the radial velocity $v_r$ of the backscattering cells are related by the following simple equation:

$$f_{DCA}\left(x,r\right)=f_{DC}\left(x,r\right)-f_{DC0}\left(x,r\right)=\frac{2}{\lambda }{v}_r\left(x,r\right)$$

(1)

where $x$ and $r$ are the azimuth and range, respectively and $\lambda$ is the wavelength. The term $f_{DC0}$ is the Doppler Centroid corresponding to a “stationary” scene, more precisely to a scene with all scatterers moving at the same velocity of the associate one with the Earth rotation.

The radial component of the target velocity $v_r$ is then evaluated by the Doppler Anomaly as follows:

$$v_r=\frac{\lambda \hspace{0.17em}{f}_{DCA}}{2}$$

(2)

Under the white Gaussian assumption for the signal and noise, it has been shown that the accuracy of such estimator is given by [30]

$${\sigma }_{f_{DC}}=0.3407\frac{PRF}{\sqrt{M}}$$

(3)

where $M$ is the number of independent samples used for the spectral estimation and PRF is the pulse repetition frequency of the illuminating sensor.

It is important to underline that the “background” DC, $f_{DC0}$ in its polynomial approximation form could be extracted directly from the ancillary information in the ENVISAT data. In our case, to account also for additional deviations of the beam that might be present in the acquisition, we preferred to perform an estimation directly from the data by a factorization this term in an azimuth $f_{DC0x}(x)$ and a range $f_{DC0r}(r)$ component:

$$f_{DC0}(x,r)=f_{DC0x}(x)f_{DC0r}(r)$$

(4)

Based on the observed experimental patterns for the ENVISAT sensor we modeled $f_{DC0x}(x)$ as a sinusoidal contribution:

$$f_{DC0x}(x)=a\sin (2\pi {\xi }_{0}x)$$

(5)

whereas for $f_{DC0r}(r)$ is expressed by the use of a polynomial form as follows:

$$f_{DC0r}(r)=a_0+a_{1}r+a_2{r}^2+.....={\displaystyle \sum _{n=0}^L{a}_n{r}^n}$$

(6)

The azimuth component was modeled as in Equation (5) because such a behavior was observed in the entire scene processed in the experimental results discussed in Section 7 and takes place at a constant frequency of 31.72 Hz. This component, whose origin is assumed to be related to the sensor, likely due to a systematic oscillation of the carrier frequency, was compensated simply via a notch filtering.

As for the range component in $f_{DC0}$, being the scene located in a coastal region involving both land and sea, we used an approach based on separate estimation in two regions starting from the sea mask. We therefore distinguished $f_{DC0r}$ between a term relative to the land $f_{DC0r_L}$ and the term relative to the sea $f_{DC0r_S}$:

$$f_{DC0r_L}(r_L)={\displaystyle \sum _{n=0}^N{a}_{L_n}\hspace{0.17em}{r}_{\hspace{0.17em}L}^n}$$

(7)

$$f_{DC0r_S}(r_S)={\displaystyle \sum _{n=0}^N{a}_{S_n}\hspace{0.17em}{r}_S{}^n}$$

(8)

where we used the subscript L and R to distinguish the terms related to the land and to the sea, respectively; and N is the polynomial order. Our tests also show that a range expansion to a fourth-order, i.e., N = 4, provides a satisfactory trade-off between the removal of the regular pattern related to the background variation of the Doppler Centroid and the preservation of components associated to the sea surface velocity.

The polynomial coefficients are estimated through a classical Least Square approach, and the resulting coefficients are combined afterward as explained in the following.

More specifically, focusing for example on the land component, we assume to have P measurement points on the land:

$$\{\begin{array}{c}{f_{DC0}}_{r_L}(r_{L_1})=f_1\\ \begin{array}{l}.\\ .\\ .\end{array}\\ {f_{DC0}}_{r_L}(r_{{}_{L_P}})=f_P\end{array}$$

(9)

thus allowing the following matrix problem formulation:

$$\left[\begin{array}{cccc}1& r_{L_1}& \dots & r_{L_1}{}^N\\ .& .& \dots & .\\ .& .& \dots & .\\ .& .& \dots & .\\ 1& r_{L_P}& \dots & r_{L_P}{}^N\end{array}\right]\left[\begin{array}{c}{a}_{L_0}\\ .\\ .\\ .\\ a_{L_N}\end{array}\right]=\left[\begin{array}{c}{f}_1\\ .\\ .\\ .\\ f_P\end{array}\right]$$

(10)

that is:

$$R_L{a}_L=f_L$$

(11)

where $R_L$ is the Px (N + 1) range matrix expressed in Equation (10) and $a_L$ is the (N+1) × 1 vector of coefficients on the land area that we need to find. Hence, now, it is possible to estimate the coefficients $a_L$ as follows:

$$a_L={\left(R_L{R}_L^{\mathrm{T}}\right)}^{-1}{R}_L^{\mathrm{T}}{f}_L$$

(12)

In the same way, it is possible the estimate the coefficients relative to the say Q points located in the sea region:

$$a_S={\left(R_S{R}_S^{\mathrm{T}}\right)}^{-1}{R}_S^{\mathrm{T}}{f}_S$$

(13)

The estimated coefficients relevant to the sea and land areas are then combined to achieve the final estimation of the range component: in this process, only the coefficients from order 1 to N are averaged in order to decouple the zero order contribution on the two different regions. Assuming, as an approximation, that the variance of the estimated coefficients is proportional to the inverse of the number of measurements in two areas, i.e,. P and Q, we applied the following combination to find t:

$${\left[a_1.....a_N\right]}^{\mathrm{T}}={\left[a_{L_1}.....a_{L_N}\right]}^{\mathrm{T}}\frac{P}{P+Q}+{\left[a_{S_1}.....a_{S_N}\right]}^{\mathrm{T}}\frac{Q}{P+Q}$$

(14)

The range component to be subtracted from the data at each azimuth line, i.e., ${f_{DC0}}_r(r)$, is achieved by substituting $a_0=a_{L_0}$, i.e., the land offset, and $a_1.....a_N$ resulting from Equation (14) in Equation (6). The final DCA estimation, is evaluated by subtracting ${f_{DC0}}_r(r)$ from $f_{DC}\left(x,r\right)$ estimated over the patches with the Madsen approach described in [29], and then applying the previously mentioned azimuth notch filter.

The final DCA is geocoded to provide a product loadable in GIS environments. During the geocoding stage, the line of sight (LOS) director cosine is evaluated and provided as an output together with the sea mask.

3. Wind Retrieval from SAR

The dominant contribution to the radar backscatter comes from the wind generated surface gravity-capillary waves through the resonant Bragg wave scattering mechanism [31,32]. Since the caps of these waves tend to align perpendicular to the local wind, the radar backscatter is also wind direction-dependent. Therefore, SAR backscatter can be used in principle to derive sea-surface wind information. Several empirically derived backscatter models, known as Wind Geophysical Model Functions (W-GMF), have been proposed in the literature [33,34,35] to relate the 10 m wind vector to SAR observations. They are generally defined as

$${\sigma }^0=\mathrm{B}0{\left[1+\mathrm{B}1\cos \varphi +\mathrm{B}2\cos \left(2\varphi \right)\right]}^n$$

(15)

where ${\sigma }^0$ is the backscattering coefficient measured on the SAR image, $\varphi$ is the wind direction relative to the antenna pointing ($\varphi =0$ indicates an upwind, i.e., a wind aligned to the range and blowing towards the radar), and the coefficients B0, B1 and B2 depend on the wind speed $w$, the local incidence angle, and the polarization and frequency of the radar beam. The Bi parameters were determined statistically after drawing comparisons with wind data from ECMWF atmospheric modeling outputs, or from operational buoy measurements. The exponent n could be 1.6 or 1, depending on the model considered. In this work, the GMF used was developed at IFREMER and is referred to as CMOD-IFREMER [36].

Unfortunately, the inverse relationship of Equation (15) is an inherently underdetermined problem for single-view measurement instruments, as happens for SAR imaging. Indeed, Equation (15) establishes that the backscatter is a monotonically increasing function of the wind speed and that is sinusoidal with $\varphi$. In particular, fixing the wind speed, ${\sigma }^0$ is maximum for the wind blowing towards (or away from) the radar ($\varphi =0°,180°$) and minimum for the wind blowing perpendicular to the radar, see Figure 2.

As a consequence, both wind speed and direction cannot be estimated from a single ${\sigma }^0$ measurement so that some a priori wind information is needed. For this research, a downscaling of the 10 m wind fields obtained from the non-hydrostatic mesoscale model Weather Research and Forecast-Advanced Research (WRF-ARW) version 3.3.1 [37] were considered as a priori wind vector information. The choice fell on this wind database since it has been already assessed for the Mediterranean basin. Indeed, [33] studied the performance of the third-generation wave model WAVEWATCH III (WWIII) in the Mediterranean by using as a wind forcing term the downscaled WRF data and validating their simulation results through buoys data provided by the Rete Ondametrica Nazionale (RON) [38], as better specified in Section 6.

Following [39], the SAR retrieved wind vector components are computed by minimizing the cost function:

$$J={\left(\frac{{\sigma }^0\left(u,v\right)-{\sigma }_m^0}{\mathsf{\Delta }\sigma }\right)}^2+{\left(\frac{u_{ww}-u}{\mathsf{\Delta }u}\right)}^2+{\left(\frac{v_{ww}-v}{\mathsf{\Delta }v}\right)}^2$$

(16)

where ${\sigma }_m^0$ is the SAR measured backscatter, ${\sigma }^0$ is the backscatter value computed by solving Equation (15) with wind vector of components $u$ and $v$, $\Delta \sigma$ = 0.078 ${\sigma }_m^0$, and $\Delta u=\Delta v=\sqrt{3 } m/s$. Finally, $u_{ww}$ and $v_{ww}$ are the reference wind components provided by the WWIII hindcast. Cost function minimization is performed by looking at the global minimum of Equation (16) computed for velocities from 0 m/s to 25 m/s with a 0.05 m/s spacing and wind direction from 0° to 360° and spacing of 1°.

4. Doppler Modeling of the Wind Component

The topic of evaluating the influence of the wind influence on the Doppler centroid has been already addressed in the literature, see f.i. [4,40]. We focus here on the investigations carried out in [28], which deal with the problem of using the SAR DCA as an additional source of information for wind retrieval. In particular, the goal in [28] was to integrate the DCA measurement in classical approaches, such as the one described in the previous section, for wind inversion of models for the Normalized Radar Cross Section (NRCS) [36]. In doing such an integration, the authors proposed and empirically derived Doppler Geophysical Model Function (D-GMF) relating the wind speed and direction to the Doppler anomaly trough parameters associated with the SAR observation geometry and characteristics. This empirical function, which is summarized in the following, was derived by analyzing a match-up database between the wind speed and direction provided by ECMWF and Doppler measurements derived by ASAR acquisitions for both VV and HH polarizations. The analysis emphasizes that the Doppler shift due to sea surface wind can be described by the wind speed ($w$), the wind direction ($\varphi$) with respect to the antenna look angle in degrees, the incidence angle of the impinging radiation in degree ($\vartheta$) and the radar polarization ($pol$), i.e.,

$$f_{W\_DCA}=CDOP\left(\varphi ,w,\vartheta ,pol\right)$$

(17)

Expression (17) is developed in [28] as

$$f_{W\_DCA}={\alpha }_{pol}F\left[X\left(\varphi ,w,\vartheta ,pol\right)\right]+{\beta }_{pol}$$

(18)

where ${\alpha }_{pol}$ and ${\beta }_{pol}$ are two coefficients depending on polarization, and

$$F\left(x\right)={\left(1+e^{-x}\right)}^{-1}$$

(19)

the argument of the function $F$ is in Equation (18) is defined as

$$X\left(\varphi ,w,\vartheta ,pol\right)={\gamma }_0^{pol}+{\displaystyle {\sum }_{i=1}^{11}{\gamma }_i^{pol}}F\left[{\mathsf{\Gamma }}_i^{pol}\left(\tilde{\varphi },w,\vartheta \right)\right]$$

(20)

where

$${\mathsf{\Gamma }}_i^{pol}\left(\tilde{\varphi },w,\vartheta \right)={\omega }_{i,0}^{pol}+{\omega }_{i,1}^{pol}{V}_1^{pol}\left(\tilde{\varphi }\right)+{\omega }_{i,2}^{pol}{V}_2^{pol}\left(w\right)+{\omega }_{i,3}^{pol}{V}_3^{pol}\left(\vartheta \right)$$

(21)

$$V^{pol}=\left(\begin{array}{c}{V}_1^{pol}\\ V_2^{pol}\\ V_3^{pol}\end{array}\right)=\left(\begin{array}{c}\tilde{\varphi }{\lambda }_{20}^{pol}+{\lambda }_{21}^{pol}\\ \vartheta {\lambda }_{00}^{pol}+{\lambda }_{01}^{pol}\\ w{\lambda }_{10}^{pol}+{\lambda }_{11}^{pol}\end{array}\right)$$

(22)

with $\tilde{\varphi }=360°-\varphi$ if $\varphi <180°$ and $\tilde{\varphi }=\varphi$ elsewhere. Notice that ${\omega }^{pol}, {\lambda }^{pol}$ and ${\gamma }^{pol}$ are polarization-dependent coefficients whose values can be found in [28].

The model described by [28] and here summarized was validated on two cases of complex meteorological situation and it is found that the high sensitivity of the Doppler to the wind direction is useful to retrieve more realistic wind patterns in the case of complex and rapidly changing meteorological situations. More detail about the results and the areas of study can be found in [28].

Figure 3a shows the absolute value of the Dopper Anomaly ($f_{W\_DCA}$), evaluated according to Equations (18)–(22), as a function of wind speed, the latter ranges from 1 m/s to 18 m/s. For the ENVISAT case, the incidence angle is equal to 23° and the polarization is VV. In particular, in this figure are shown the upwind case (green continuous line), which is the case of wind blowing toward the antenna, and the downwind (red dashed line) case, which is the case of wind blowing away from the antenna. In both cases, the measurements shows that a reduction of the wind leads to decrease of the Doppler frequency; however, the DCA does not fall down beyond 8 Hz, corresponding to about 30 cm/s surface radial velocity. Notice that this figure is comparable to Figure 1b reported in [28]. In Figure 3B it is shown the DCA, as a function of wind direction ($\varphi$), as defined above, for three different values of wind speed: 7 m/s blue line, 5 m/s green line, 2 m/s red line. In this case, the DCA and the surface current both approach zero in the case of orthogonality between the wind direction and the radar LOS.

As our focus was on the analysis of the estimated sea surface currents, we exploited the D-GMF in a direct form, which is the opposite of the method used in [28]. More specifically, given the radar geometry and the wind speed and direction provided estimated by the Bayesian retrieval involving the hindcast of the meteorological ocean model and the measured NRCS a Doppler Centroid shift $f_{W\_DCA}$ were evaluated through the D-GMF in Equations (18)–(22) and converted to surface velocity through Equation (2).

5. The Data Analysis Procedure

The former sections addressed the separate description of the different elements considered in this study to investigate the possibility of exploiting the Doppler information derived from the SAR data to achieve products that can provide an indication about the sea currents. As previously noticed, the sea surface current directly extracted from the SAR data can be heavily impacted by meteorological factors, mainly the wind.

The procedure used in this work for the data analysis is depicted in Figure 4: it involves the joint processing of SAR acquisitions and the integration with the external information about the wind achieved by the external atmospheric model.

Figure 4 highlights that SAR data are exploited following two main streams based on the use of the intensity and of the complex (amplitude and phase) information for the (spectral) Doppler analysis. As for the latter, the DCA is generated according to the procedure described in Section 2 and summarized in Figure 1 starting from the complex information associated with the raw data. The choice to start from the Level-0 product is, as already pointed out in Section 2, related to the need for avoiding any possible bias introduced at the data focusing stage for the generation of Level-1 products supplied by the data provider. The latter may in fact incorporate filtering for image enhancement based on ad hoc trade-off of the resolution and the sidelobes suppression. The achieved DCA, one converted to surface current (point 6 left output), is then stored for the subsequent comparison with the response associated with to meteorological factors carried at the final stage of the processing.

The estimation of the Doppler component induced by the meteorological factors follows the right side flow in Figure 4 and is composed of two sequential steps. In the first step, the wind information provided by the weather forecast model is enriched via the assimilation, based on the Bayesian strategy described in Section 3, of the High Resolution SAR intensity. A High-Resolution wind map is thus also produced (level 3).

The second step exploits the Doppler GMF described in Section 4 to derive the wind related Doppler information from the wind intensity and direction achieved after from both the original wind data by the model and the HR wind results after integration with the intensity SAR data. The retrieved wind Doppler components are converted to surface current, see level 6 at middle and right, according to Equation (2). At this stage, a final subtraction between the DCA derived surface current and the one derived by the wind information could be implemented to extract the residual radial surface current component. As long as wind data derived from the weather forecast are given with a high confidence level, the output produced by the procedure shown in Figure 4 should therefore contain only the information relevant to the current associated with the ocean current and possible effects of the sea state. This subtraction is, however, not included because of the aim of this work was to evaluate the influence of the wind on the surface currents derived from the Doppler analysis.

The above described procedure was adopted in the experimental section for the data analysis of the case studies. The green numbers shown in Figure 4 identify the results of the different processing steps, which are shown and described in the experimental analysis carried out in Section 7. To summarize, they are

(1) The HR detected ground image (SAR intensity).

(2) The Wind vector information from the external atmospheric model.

(3) The Wind vector retrieval from SAR intensity image (point 1 in Figure 4) using as a-priori guess the external atmospheric model (point 2 in in Figure 4).

(4) The Estimated wind driven Doppler component.

(5) The Doppler Centroid Anomaly.

(6) The radial surface current computed from the DCA (point 5 in Figure 4) and the Doppler modeled starting from the two different wind information (point 4 in Figure 4).

6. Test Site and Data Set

The test site is situated in the south part of Italy, see Figure 5 on the left, and covers a coastal area of about 90 km × 100 km (see Figure 5 on the right). More specifically, it starts in the south of the Gulf of Gaeta, includes the Gulf of Naples with the Ischia, Procida and Capri Islands, and part of the Gulf of Salerno.

In Figure 5 left, red placemarks were highlighted to delimitate the corners of the validation area of WRF-ARW wind model. The validation was performed at the sites that are vertices of a quadrangular containing the study area of this paper: Ponza, Cetraro, Palermo and Capo Comino. Moreover, [41] developed hindcast covering 32 years (from 1979 to 2010) and containing, among the other, the wind downscaled data used in the present study.

In this work, we focused our study on the coast described above because it is an area densely urbanized and characterized by heavy maritime traffic. Hence, the analysis carried out could be functional for the study of the sea state and have an interesting outcome for the activities oriented toward the definition of a long-term strategy in support of sustainable growth in the marine and maritime sectors, i.e., those known as Blue Economy or Blue Growth. Moreover, the results on the surface current obtained in this work could also be useful for a meteorological analysis.

Among the several available SAR sensors, in this work, we decided to use ENVISAT ASAR for its characteristics relative to the swath. As a matter of fact, ENVISAT, in case of the (standard) stripmap acquisition, provides seamless and continuous coverage over a single swath width up to 100 km.

The swath analyzed in this work is shown on the right side of Figure 5 and delimitated by the yellow placemarks and red line. For the case of COSMO-SkyMed the swath width, still considering the stripmap acquisition, would have been approximatively halved in both directions. A more recent SAR sensor, which offers a very wide swath and operates in C-Band, is Sentinel-1 (S1), which can boast swath of about 250 km. S1 captures three sub-swaths using Terrain Observation with Progressive Scans SAR (TOPSAR) mode; this mode allows covering wide coastal areas for each single orbit. The operation mode provides, however, additional complexity in the estimation of Doppler centroid for the azimuth beam steering reasons, and would lead from burst to burst to the presence of discontinuities in the estimation Doppler Anomaly due to the change on the borders of the antenna pointing: the discontinuities impacts also the final radial velocity product [28].

Following the above considerations, as the best choice for this area, we found 47 images acquired by ENVISAT over ascending orbits between 2002 and 2010 (i.e., almost the whole operative period of the sensor). The ENVISAT sensor operates at C-band (wavelength equal to 5.6 cm), in VV polarization. Such a sensor is characterized by a spatial resolution of 5 m by 20 m in azimuth and ground range, respectively.

The dataset was screened for the presence of potential amplitude patterns related to the sea-air interface, in particular from the wind or from mass displacements of marine currents. Among them, a total of five images were selected for the subsequent numerical analysis on the Doppler characteristics, in particular, they span from 2004 to 2010 and are listed in

Table 1. Acquisitions dates and used models.

Acquisition Date	WRF-ARW Model	Oceanografic Numerical Model
17 November 2004	✔
22 December 2004	✔
02 November 2005	✔
11 January 2006	✔
22 September 2010	✔	✔

together with the indication of the available model/models.

It is worth noting that the dataset presented had the aim of discussing the open issues, especially concerning the influence of wind on surface currents and hence, represents a case study, rather than a statistically significant analysis. We underline that the choice of the period is related to the availability of ENVISAT data that, as discussed above, are characterized by an operational mode more suitable to the scope of the work.

7. Results

In this section, we discuss the results of the numerical analysis carried out on the selected images of Table 1.

DCA was evaluated, in particular, according to the method described in Section 2 with patches sliding along the whole SAR image. We considered patches of width 512 azimuth by 128 range pixels corresponding to a 2.5 km- by 2.5 km-wide window. The size of the patch was set in such a way to guarantee a good balancing between spatial and spectral resolutions. In particular, this latter is equal to 2.85 Hz. The stationary component $f_{DC0}\left(r\right)$ was evaluated by removing azimuth independent range polynomial components of the fourth order estimated on land and sea regions. The resulting radial velocity was estimated from (1): the notch filter for the removal of the undulations of the Doppler Centroid modeled as in Equation (5) and present on all the processed images.

On the other hand, the wind vector retrieval through CMOD model was performed by multi-looking the SAR data, i.e., averaging the number of cells, in order to reduce the speckle noise. As a trade-off between spatial resolution and speckle noise, the removal a multi-looked pixel of 1 km × 1 km was considered. The wind a-priori guess, provided by the 10 m wind fields obtained using the non-hydrostatic mesoscale model WRF-ARW, is instead downscaled to a computational domain model, covering the whole Mediterranean, with a ~10 km resolution Lambert conformal grid. Initial and boundary conditions for the atmospheric simulations with the WRF-ARW model were provided from the CFSR (Climate Forecast System Reanalysis) database [42].

Figure 6, Figure 7, Figure 8, Figure 9 and Figure 10 show the results obtained for the five dates listed in the previous section, organized in a panel comprising eight images. Specifically, for each date are shown:

• the ${\sigma }^0$ image (point 1 in Figure 4);
• the estimated Doppler Anomaly maps converted to the sea surface (line-of-sight component) velocity (point 6 left arrow in Figure 4);
• the wind direction obtained from the (Bayesian) CMOD and from the WRF-ARW model (point 3 in Figure 4);
• the wind speed relative to the CMOD model and to the WRF-ARW model (point 3 in Figure 4)
• the surface radial velocity due to the wind component, evaluated by the D-GMF as described in Section 4 for both winds estimated by WRF-ARW and refined by the CMOD-based Bayesian model (point 6 in Figure 4, middle and right products).

It can be noticed that the estimated velocity ranges between ± 100 cm/s, red and blue colors (redshift and blue shift) of the estimated sea surface velocity component indicate a motion toward and away from the radar antenna, respectively.

Figure 6 refers to the acquisition of 17 November 2004; in this case, the amplitude image (Figure 6A) shows an evident pattern in the north of the Ischia Island area. This pattern is observed also in the image corresponding to the surface currents (Figure 6B) with a redshift located in the same position. The image in Figure 6B also shows another pattern with a component of the surface current toward the sensor in the south-most region of the swath: no significant variation, however, is observed in the same area in Figure 6A. Both patterns can be explained, although partially, by looking at the wind data of the WRF-ARW model wind: the arrows’ direction (Figure 6D) indicates in fact the presence of a wind blowing towards the sensor (upwind condition). However, the north pattern has a sharper transition than the southern one in terms of wind direction. In the most westerly part of both Figure 6C,D, a wind blowing almost parallel to the azimuth direction, i.e., in the direction orthogonal in the LOS, can be noticed: a feature that seems to be correctly mapped in the DCA measurement on the left part of Figure 6B. Such a transition is even more emphasized in the refined wind measurements following the MAP based on CMOD estimator for the wind described in Section 3 (see Figure 6C,E,F for the intensities). The surface radial velocity obtained via the Doppler model described in Section 4 shows a rather good match with the Doppler measurement achieved by processing the complex SAR data (Figure 6B) for both the WRF-ARW (Figure 6H) and the CMOD wind (Figure 6G). However, the CMOD derived surface velocity (Figure 6G) shows a higher spatial match with Figure 6B when compared to the WRF-ARW derived surface velocity in Figure 6H. Finally, in Figure 6G it is noticeable only a partial mismatch on the northern pattern, likely due to a bias of the WRF-ARW wind source on the estimation of the wind direction in CMOD MAP based estimator.

In Figure 7 are presented the results relative to the acquisition on 22 December 2004, also in this case, in the amplitude image (Figure 7A), there is a pattern in the northern part of Ischia Island (less evident than the case of Figure 6A), which is also present in the surface radial velocity of Figure 7B. In the south-east area of the DCA map (Figure 7B) a redshift is again observed; it represents a surface current that moves toward the sensor. Wind directions of Figure 7C,D can explain such a behavior because they show, also in this case, the presence of an upwind: the wind magnitude (see Figure 7E,F) confirms this evidence too. Finally, Figure 7G,H highlight that in the upwind condition the radial velocity obtained via the D-GMF described in Section 4 tends to overestimate the radial velocity, i.e., shows a redshift also in areas of with low wind magnitude. This behavior is in agreement with the plots in Figure 3, which show that the in the upwind condition the D-GMF provides radial velocities different from zero even in the case of very low wind magnitude (Figure 3 left), whereas with wind directions orthogonal to the LOS (Figure 3 right), it does not provide contributions. Despite the bias for low wind, there are many congruencies between the patterns shown in Figure 7B,C, especially in areas with substantial wind conditions.

Figure 8 shows the results for the acquisition on 2 November 2005. Moreover, in this case, the presence of a strong backscattering pattern in the northern area and in the south part of Figure 8A it is evident and the surface radial velocity shows a similar pattern, confirming a good correspondence between the amplitude and the product derived by the DCA and hence, likely an effect of the wind. Observing Figure 8C, and similarly Figure 8D, it is possible to notice a change of the wind direction: it starts blowing from land to sea then it rotates orthogonally to the LOS. Thiss change is mapped the red to green shift present in Figure 8B. Wind intensities evaluated from the WRF-ARW model (Figure 8F) and CMOD MAP based estimator (Figure 8E) look slightly different. Anyway, the surface radial velocity due to wind component (Figure 8G,H) shows that the redshift in Figure 8B can be due to the presence of wind, which is in a quasi-upwind condition. It is also noticed that, in agreement with the plot in Figure 3 right, the change of wind direction in the western part is mapped in a decrease of the surface current signal. In this case the CMOD derived surface velocity (Figure 8G) is in a better agreement with the behavior of Figure 8B than the WRF-ARW derived surface velocity (Figure 8H).

The results for the acquisition of 11 January 2006 are presented in Figure 9: they again involve specific patterns present in the amplitude image (Figure 9A) and in the surface radial velocity (Figure 9B) in the Gulf of Naples, as well as in the north and south regions. It can be noticed that CMOD wind and WRF-ARW show differences mainly in the wind intensity (Figure 9E,F), which also impacts the direction (Figure 9C,D) especially in the region in the north of the Gulf of Naples. Here, it is rather evident in both the amplitude (Figure 9A) and surface velocity (Figure 9B) images in the north region that the decrease of the surface current in the north-east area, with reduced wind intensity, is impacted by the change of direction. The intensity increase of surface velocity visible in the south part of Figure 9B is on the other hand due to the quasi-upwind condition, particularly evident in Figure 9D and to the high wind intensity shown in Figure 9F. The CMOD derived surface velocity (Figure 9G) is in agreement with the surface radial velocity (Figure 9B), while for the WRF-ARW derived surface velocity (Figure 9H), a distributed overestimation of the model essentially due to the alignment of the wind to the LOS is present. The comparison of Figure 9B,G shows also a decrease of the surface velocity estimated with the DCA: this reduction could be originated by the subtraction polynomial range component discussed in Section 2 at the stage of generation of the DCA.

Finally, in Figure 10 are shown the results for the acquisition of 22 September 2010. This case study is interesting because, as in some other cases, there is a change of the wind direction going offshore, but the wind is rather intense all over the image, see Figure 10D. In the amplitude image (Figure 10A), a pattern in the intensity image in the northern region is evident, which is also mapped by the surface velocity derived by the DCA in Figure 10B with a redshift corresponding to a current toward the sensor. This pattern, however, is not mapped in the intensity of the WRF-ARW model (Figure 10F). Apart from the area shown in the south-east part of the image, the surface current derived by the WRF-ARW model (Figure 10H) does not show a satisfactory match with the measured DCA shown in Figure 10B. The pattern in the north of the Gulf of Naples in Figure 10B, is only roughly mapped in Figure 10H. The situation drastically changes with the wind derived by CMOD. The pattern with large intensity in the north of the Gulf of Naples in Figure 10A provides, in fact, a refinement mainly affecting the wind intensity, see Figure 10E. This refinement originates from a better matching of the surface current estimated by the wind in Figure 10G and the surface current derived by the DCA in Figure 10B with respect to Figure 10H.

In order to a comprehensive interpretation of the results with respect to the circulation currents, the surface current in the area has been investigated by using an oceanographic numerical model for the date of 22 September 2010. Specifically, a state-of-the-art ocean model (SHYFEM, [43]) was applied to compute the circulation patterns in the area during the whole of September 2010. The model has been widely used to reproduce with and predict with success the three-dimensional circulation patterns induced by wind, tides and thermohaline forcing, at both coastal and shelf scale in the Mediterranean basin [44,45,46,47,48,49]. A simulation run was carried out considering the main meteo-marine forcing affecting the local circulation as the heat flows to the surface, the direct wind action and thermohaline forcing. Atmospheric and oceanographic data needed as model surface and boundary conditions were provided by the operational meteorological prediction system SKIRON [50] and by the European Earth Observation Programme Copernicus, respectively. The same model parameters and simulation setup used in [43] was adopted, unless the model domain that, in this application, included the Naples Gulf and the surrounding coastal and shelf areas.

In Figure 11, upper panels, the surface current fields obtained from simulation results between the 0 and 3 m depths, corresponding to the model first vertical layer, are reported for the interested domain at the time between the 20:00 and 21:00 UTC on 22 September 2010. In both panels the flow pattern is characterized by the presence of an intense coastal flow directed northwestward with maximum current speeds higher than 0.3 m/s.

In the nearby Gulf of Naples, an HF radar system is operationally providing the water current at the surface but limited to the bay area [23,51]. A comparison between the measurements derived with the Doppler Centroid Anomaly from SAR and the measurements of the HF radar system available in the bay is reported in [52], and shows reasonable agreement although the limitation of covered area.

On the other hand, several studies based on the use of oceanographic models have been carried out to reproduce the water circulation in both the Gulf of Naples and the surrounding coastal zone [23,27]. In particular, in [27] a similar modeling application was carried out applying an ocean model to investigate the water circulation in both the shelf and coastal area. Even if in the two applications the temporal coverages differ, there is an evidence that in [27] results obtained for both Spring and Autumn periods reproduce exactly the same hydrodynamic patterns reported in Figure 11. In particular in [27], both the model results and the analysis of the geostrophic currents depicted a northward flow leaning on the shelf break that interacting with the coast geometry generates cyclonic circulation patterns on the adjacent coastal area (see Figure 3, 13 and 14 in [27]). The obtained results are in line with previous studies in the area indicating the northward coastal flow as a typical circulation structure characterizing the area [27].

Changes in the flow directions are detected near the coastlines and on the leeside of the main Capes and Islands. Specifically, the interaction between the coastal flow, the coastline geometry and the wind action give rise to an edge separating the main external flow from the circulation cells in the shallow water part of the domain, north to the Gulf.

Similarly, in Figure 10B, inversions of the signal direction are detected for the same areas, thus confirming the presence of potential oceanographic heterogeneities as fronts or edges. Considering the wind role in modulating the pattern flow in this area, in the off-shore part of the domain, due to the higher water depths, the Coriolis force allows a direct contribution of the wind action to the main northwestward flow. In the shallow part of the domain, the off-shore wind direction and the complicated coastal geometry generate highly variable oceanographic structures that cannot be accurately solved by the adopted numerical model application. Within this context, thought the model results are very helpful in detecting the edge between open ocean and coastal circulation patterns, they cannot provide accurate information on the wind induced transport processes characterizing the skin layers of the coastal waters. Consequently, in the lower panels of Figure 11, the decomposition of the modeled flow field into the radial component is not favoring the detection of a direct correspondence between the wind speed and direction and the effects on the computed surface circulation.

8. Discussion and Conclusions

SAR sensors have significant applications in the marine sector due to the high sensitivity to different phenomena occurring under and above the sea surface.

SAR data can be used to extract the sea radial surface velocity, a key variable in maritime scenario. In the present work, we investigated the potentialities and the limitation in the extraction from ENVISAT data of information about the sea currents in the coastal area of the Gulf of Naples (South Italy) in the Mediterranean Sea.

In the cases examined in the previous section, rather typical of the area considered, we have seen that the wind shaped the Doppler Anomaly derived from the analysis of the SAR images. One of the reasons of that is the weakness of the surface currents, which, according to the ocean model hindcast (Figure 11), roughly range between ± 5 cm/s, which is similar to the resolution of the Doppler Anomaly derived from SAR. Despite the inability to quantify the effect of circulation currents on the Doppler measurements, the study does provide evidence of the effects in terms of the observed Doppler patterns.

Doppler Anomaly measured form SAR image shows some differences with the model-derived one, a fact that may be ascribed to the coarse spatial resolution of the wind’s model, which was adopted both to force the Doppler Geophysical Model Functions [28] (D-GMF) and as a first guess to the retrieve wind field from SAR (D-GMF).

The results presented in this paper suggest that this well known issue [53,54] could be partly overcome. The coarse spatial resolution of external data could indeed be refined by wind/Doppler approach that highlights changes in wind direction, as shown in this paper. By combining D-GMF and W-GMF with an atmospheric wind model, a wind vector fine map could be thus obtained.

Once all the above problems are solved, the Doppler contribution produced by the D-GMF should be subtracted to that derived from SAR, to perform a validation of the residual Doppler Anomaly. Future investigations of the effects of wind and currents on SAR derived surface currents should thus consider coastal regions for which current data (e.g., from moored Acoustic Doppler Current Profiler velocity data [55], High-Frequency Radar currents (e.g., [56]), are available for a direct quantitative comparison. It should be remarked that over large areas (as those studied here) that cannot be sensed by the coastal HF radars, there is the need of a more detailed ocean surface modelling assimilating the available in-situ data, possible nowadays with the enhanced computer facilities and the data assimilation techniques. Along this line, altimetry-derived surface currents (e.g., [57,58]) can be exploited also in coastal regions [59] as an additional source of information.

The work along these lines would provide a contribution to the exploitation of SAR derived currents products to the advance of coastal ocean observing systems [60].