Remote Sensing-Based Automatic Detection of Shoreline Position: A Case Study in Apulia Region

: Remote sensing and satellite imagery have become commonplace in efforts to monitor and model various biological and physical characteristics of the Earth. The land/water interface is a continually evolving landscape of high scientiﬁc and societal interest, making the mapping and monitoring thereof particularly important. This paper aims at describing a new automated method of shoreline position detection through the utilization of Synthetic Aperture Radar (SAR) images derived from European Space Agency satellites, speciﬁcally the operational SENTINEL Series. The resultant delineated shorelines are validated against those derived from video monitoring systems and in situ monitoring; a mean distance of 1 and a maximum of 3.5 pixels is found. This paper presents a new methodology for shoreline detection, based on edge detection techniques of preprocessed SAR images, freely available from ESA and the Copernicus program. The developed procedure of land/water boundary recognition is carried out with a script, which allows to generate both shorelines and error statistics given the availability of in situ data. When compared to the outcomes from in situ and monitoring programs, the obtained shorelines provide reliable results given the spatial resolution of the satellite data. Moreover, in the case of calm sea and high contrast in the imagery, the method shows good performances in detecting shoreline from aerial images.


Introduction
Coastal areas represent complex and highly dynamic environmental systems continuously facing undergoing changes due to natural as well as anthropogenic factors. The monitoring of their evolution is crucial for the safeguarding thereof. Wind, waves, coastal storms, and sea level rise represent some of the natural coastal threats that have been shaping the coastal zones and the consequences of which have been exacerbated in the last century [1,2]. Sea hazards and extreme storm surges are increasingly threatening coastal areas causing among others flooding, coastal erosion, and damages to ecosystems and infrastructures [3]. Urban pressure, shifting land use, and variations in dunes and sea bed stability [4][5][6] have and are also contributing to a further exacerbation of these hazards. Many coasts throughout Europe, and specifically in the Mediterranean [7], are expected to experience increased flooding in response to current coastal and dune management strategies.
In this context, effective management is of paramount importance to address the effect of the high concentration of human population and activities along the coasts. The assessments of coastal vulnerability and hazards can provide decision makers and policy with diagnostic management tools to support sustainable development of coastal zones and face the negative effects of both human activities and extreme events [8,9]. Environmentally friendly interventions and nature based solutions may help in this direction promoting the restoration of the environment and the ecosystem [10,11].
Information on long-term shoreline variability and periodic surveys of emerged and submerged beaches [12] help to quantify erosion or accretion processes. Among the several physical indicators used to assess beach erosion/accretion is the delineation of shoreline trends [13,14].
An idealized definition of shoreline is the physical interface between land and water [15], identified by means of [16] (i) visually discerning coastal features, (ii) the evaluation of the intersection between beach profile and the local tidal datum, and (iii) specific image processing techniques.
Traditionally, the coastline has been detected by manual visual comparison and interpretation of different sources and images derived from all manner of Earth Observation (EO). However, by using subjective visually discernible features, shoreline position may be strongly influenced by the operator as the boundary between water and land can be difficult to identify due to operational constraints and can vary depending on tidal cycle information. Moreover, such manual delineation is time-consuming, resulting in a high degree of inefficiency.
To overcome such limits, several techniques have been proposed to semi-automatically and automatically extract coastline from images acquired by both remote sensing [22,35] and video monitoring systems [28,36]. Remote sensing systems provide a repetitive and consistent view of the Earth, useful for monitoring short-and long-term changes of coastal zones and specifically shoreline position on a global scale with minimal individual user infrastructure requirements. Video monitoring systems generally consist of more than one camera to cover multiple view angles on the area of interest and typically provide data at very high spatial as well temporal resolution. Usually, the spatial scales of such systems vary from decimeters to kilometers based on the distance of the target from the camera. Images are acquired on a specified time sampling period varying from minutes to hours and are generally automatically collected and processed with low operating costs. A further benefit of such systems is their indifference to most atmospheric obstructions, mainly weather features and events such as clouds, allowing for data to be collected continuously. Camera systems also offer a higher degree of customization per local monitoring objective, as different types of images can be collected to obtain information. For example, an averaged image is useful for the submerged and bar topography [37], while other sophisticated operational video analysis techniques enable the quantification of detailed coastal features such as shoreline evolution and wave run-up [38,39]. Moreover, the evaluation of waves, flow, and rip currents empowers the management of dynamic navigational channels [40][41][42], whereas visitor density detection allows the monitoring of defining stressing factor for zones having an important ecological role [43,44]. Hyperspecific features or data requirements from satellites are not as easily customized, though a wide portfolio of private purchasable data series complements the readily and publicly available data.
Satellites are more sensitive to weather conditions than video systems (particularly for passive and optical data acquisition, while radar and active sensors are limited to a lesser degree) and provide imagery on a lower temporal frequency than video camera systems. However, when considering large temporal time frames, the persistence of satellite imagery provides a powerful resource; whereas video monitoring campaigns last years, and satellite imagery spans decades and continues to grow.
Remote sensing instruments fall into two categories: passive and active. Passive sensors measure various radiation naturally reflected or emitted from Earth's surfaces: atmosphere and clouds. Passive sensors are primarily employed in the monitoring of morphological changes of the near-shore zone and for the monitoring of water quality. ESA, NASA, and other national and international agencies operate an array of satellite missions providing large amounts of data on key physical and environmental indicators including total suspended matter concentration, algal pigment concentration, chlorophylla [45,46], and sea surface temperature. Furthermore, postprocessing techniques can derive currents [47], alongshore and offshore significant waves and tides. New ocean color sensors, such as the Plankton, Aerosol, Cloud, ocean Ecosystem (PACE) sensor, will enable the tracking of new ecological indicators and thus contribute to water quality monitoring.
Active remote sensing, which generates its own signals and records the amount of energy backscattered from the terrain, is capable of passing through gases and clouds collecting clear data regardless of any light availability and weather conditions, a great advantage when compared to passive sensors.
The most common methods to detect the shoreline from active sensors are based on segmentation and edge detection. Segmentation is a process of grouping objects (or pixels) having similar characteristics to later be assigned to a specific class. Several methods have been developed for image segmentation, mainly belonging to a supervised or unsupervised method. Supervised segmentation uses a priori training dataset (or input classes) selected based on the experts' knowledge, while unsupervised segmentation is trained during the segmentation process allowing for a fully automatic recognition of the classes. Postsegmentation processes might be required to merge regions that may derived from an over-segmentation [48].
To overcome these limits, Dellepiane et al. [35] proposed the adaptation of fuzzy sets in defining the shoreline contour. Coherence information extracted by the European Remote sensing Satellite (ERS)-1 and ERS-2 with a spatial resolution 20 × 20 m are used as input of the segmentation process. An accuracy of about 3.5 pixels mean offset resulted by a comparison with the shoreline position detected from aerial images.
Edge detection is a process of finding boundaries between different regions and it is generally simpler than segmentation. When pursing edge detection, locally adaptive thresholding algorithms can be employed. Thresholds are analytically determined by fitting a bimodal Gaussian curve, and the accuracy of the results depends upon their reliability and correctness [48].
Previous studies [49][50][51] used an edge detection method to determine the boundary area between land and water. Niederman et al. [51] apply an edge detection method [52,53] on ERS satellites images. At first, all edges above a certain threshold are recognized and then, a local edge selection is performed to refine the determined land/water boundary. Compared with shoreline extracted by visual decision, the authors found a mean offset of 2.5 pixels.
In general, edge detection methods are simpler than segmentation, however edges may be discontinuous leading to an erroneous definition of the shoreline [49,50], a drawback that may require additional postprocessing.
The present work shows a new automatic method based on edge detection aiming at retrieving shoreline position from Synthetic Aperture Radar (SAR) images. To extract the contour on land/water boundary, the method relies on a global threshold value derived by means of Otsu method [54]. This is a key component of the method since determines its simplicity in comparison to locally adaptive thresholding algorithm as it does not require the partitioning of the whole image in sub-imaged or statistic examination of the intensity values of the local neighborhood of each pixel [55]. The method allows for a fully automatic detection of the coastlines in a short amount of time. Sentinel 1 SAR images with a spatial resolution of 20 m are used as input.
Understanding and quantifying the positional accuracy of the satellite-derived shoreline positions is essential. For this reason, the method is tested along Apulian beach in Southern Italy, Torre Canne (Brindisi), facing the Mediterranean Sea, where spatial and temporal variations of shoreline derived from a video system are available. The accuracy of the method in detecting the shoreline position (Satellite Derived Shoreline (SDS)) is quantified with respect to shoreline data collected by a video system. The comparison is quantified and reported in terms of average offset and root mean square difference.
The following Section 1.1 describes the morphological characteristics of the study area. The data employed in this study and the methodology steps are described in Section 2. The results, actual validation, and main findings are discussed in Section 3, while a summary of the main outcomes with a recommendation for further study is addressed in Section 4.

Study Area
The algorithm is tested along Torre Canne beach (Figure 1), in the south of Italy. Located in Apulia region and facing the Adriatic Sea, Torre Canne beaches are characterized by medium grain size. The emerged shore width may vary from a minimum of about 20 m to a maximum of about 50 m with a gentle foreshore and a steep berm, periodically covered by organic deposits of Posidonia oceanica. During the last decades, despite the mean wave climate is moderate (mostly propagating from N-NW and E-SE with an annual significant wave height comprised between 0.5 and 1.4 m and peak period ranging between 3 to 7 s), the study area has suffered remarkable erosive processes. Human activities and overexploitation of the coastal area, mainly due to touristic activities, combined with the increase of extreme storms occurrences have been negatively affecting the sedimentary balance of the area. The long-term erosive trend has imposed an extensive monitoring program of the area, including topographic and bathymetric surveys. In December 2015, a new video coastal monitoring system was deployed and integrated into the Apulian Region Monitoring Network, managed by the local Basin Authority to study both morphodynamic and hydrodynamic processes [28]. Further detail of the system is reported in Section 2.1.2.

Materials and Methods
To support decision-makers and policy with a management tool for coastal zones defense, this study proposes a new methodology for shoreline detection using Copernicus Sentinel-1 data. This Section describes the datasets used to develop and validate the method as well as the methodological steps employed in this research.

COPERNICUS Sentinel-1 Satellite Data
Sentinel-1 is the first satellite mission developed by European Space Agency (ESA) for the Copernicus Programme, the European Union's Earth observation programme, previously known as GMES (Global Monitoring for Environment and Security). This program aims to serve as Europe's eyes on the Earth and to provide high-quality and high-resolution EO data to public and private industries to aid in the monitoring and strategic planning of activities. The Sentinel-1 carries a C-band Synthetic Aperture Radar (SAR) instrument capable of collecting data regardless of weather conditions or light availability. These data are freely available and downloadable via the Copernicus Open Access Hub web portal (https://scihub.copernicus.eu/, accessed on 25 May 2021) or the five European Data and Information Access Services (DIAS) platforms that allow users to discover, manipulate, process, and download Copernicus data and information. In the present work, Sentinel-1 products acquired on Interferometric Wide swath (IW) mode are elected to test the developed methodology. IW mode combines large swath width (250 km) with high geometric resolution (5 × 20 m on the ground) and is implemented as Terrain Observation with Progressive Scans SAR (TOP-SAR) mode [56]. Depending on the electromagnetic wave transmission and reception directions, data are available in single or dual polarization. Partial dual polarization vertical+horizontal (VH) is tested here. Dual polarization HV and VH have been proved to be adequate for the mapping of flood water [57] and to be less sensitive to wind-induced surface roughness that can influence the backscattering coefficient, especially for VV polarization [58,59]. Moreover, VH polarization generally yields to highest contrast between land and open water and thus is suitable for the purpose of this study.
IW products are acquired at Level-0, the raw data, and processed to Level-1 and Level-2 products through filtering techniques and processing algorithms to obtain higher level information. Level-1 Data (L1D) preprocesses include radiometric normalization, terrain correction, and geolocation. L1D products are available at single look complex (SLC) or ground range detected (GRD). SLC products are provided in slant-range geometry where coordinates are defined as the line-of-sight from the radar to each reflecting object and with variable pixel resolution. GRD products are SLC projected using an Earth ellipsoid model. The resulting pixels have approximately square resolutions and square pixel spacing. Depending on the level of multi-look applied, GRD files are available in Full (FR), High (HR), and Medium (MR) Resolutions. Sentinel-1 GRD HR products acquired on IW mode are suitable to test the proposed algorithm for shoreline detection, due to the easy and low-cost data accessibility and usability. Products have a pixel size equal to 10 × 10 m and a spatial resolution of 20 × 22 m as a result of the aggregation of pixels in the postprocessing phase. An example of the downloaded images is shown in Figure 2.

Validation Data
In order to obtain accuracy and precision metrics, optical images acquired by a video monitoring system installed in 2015 at Torre Canne, Italy (Apulia, South Italy) [28,36,60] are used to compare satellite versus in situ monitoring strategies. Video-derived images spatial resolution varies from decimeters to around 13 m, based on the distance of the near-alongshore footprint from the camera (Figure 3). The camera system is configured to acquire images with a sampling interval of 30 min. The acquisition cycle is of 10 min with a frequency of 1 Hz, thus resulting in 600 images. Those images are automatically processed to extract the shoreline position. The averaged shoreline is then returned every 30 min as a collection of points geolocated into a specific UTM projection identified with the code EPSG:6708. Those points can be transformed into a continuous element defining the shoreline. Results have been validated in a previous work [28] by a comparison with shoreline manually recognized by an operator. A mean difference of 1.17 pixels was found out.
The time of acquisition of both satellite and video system is a key element when comparing the two methods as shoreline detection from image processing is highly influenced by both tide and wave run-up levels. Those can indeed influence the water content and thus the saturation degree of the emerged beach and, consequently, the colors and the reflection proprieties of the sand, based on which the land/water interface is distinguished. Video derived shorelines are already processed by considering the effect of the tide but not satellite images. For this reason, in order to ensure a reliable validation of the algorithm developed within this research, only the images acquired with a maximum lapse of time of 20 min are compared (e.g., on the 29th of September, the image acquired at 16:48 by Sentinel-1 is compared with the image acquired at 16:30 by the video system), namely, the shoreline position variations are assumed to be negligible within the time laps of 20 min.

Methods
The algorithm is based on edge detection techniques, mainly consisting of finding the boundary which separates two different regions. Figure 4 below reports the block diagram of the method. It consists of four main processes: despeckling, binarization, morphological operations, and edge detection. Those processes aim to (i) reduce the speckle effects (despeckling), (ii) recognize water/land boundary by intensity values (binarization, morphological operations), and (iii) automatically extract the shoreline (edge detection). The procedural architecture for automatic shoreline edge extraction from preprocessed Level-1 GRD satellite images obtained from Copernicus Services is built on Python Open Source Computer Vision (OpenCV) library [61] and Scikit-image processing toolbox [62] for SciPy.

Despeckling
Due to the constructive/destructive interferences of backscattered electromagnetic waves, a speckling appears as a salt-pepper noise in the image. The first operation aims at reducing speckle effects on SAR images to obtain an image characterized by lower noises, making the output easier to be processed in subsequent functions. Without this step, such disturbances preclude the further recognition of the water/land boundary by intensity values.
The reduction of speckle noise can be obtained by applying a multi-look integration or a post-image formation method [63]. Multi-look processing technique averages, pixel by pixel, several independent or dependent images, corresponding here to different looks of the same image [64]. Each look is obtained by dividing the signal bandwidth into parts and processing each part as a single look image of the same scene [65]. Post-image formation methods are based on the application of standard filters, such as Huang [66], Kuan [67], Frost [68], and Gamma Filter [69]. Usually, those filters perform well on most SAR images even though spatial resolution is often lost [63]. Due to its low computational cost, simplicity, and good performance in preserving edges [70], a median filtering approach is used in the present work to Sentinel-1 images. Median filtering replaces the intensity of each image pixel y(m, n) with the median of the intensities in a defined neighborhood w around the pixel location (m, n) as defined in the following Equation (1): The neighborhood is defined by the user. The median filtering approach is implemented by means of a structuring element, also called kernel, with a size of 5 × 5.

Binarization
Pixels classification is then applied to the denoised image via the Otsu's method [54] to allow the creation of a binary image. The Otsu's method is one of the best-known thresholding methods and has been proved to obtain satisfactory results when applied to satellite images for the land/water classification [71][72][73]. The Otsu's method is an unsupervised and nonparametric method based on the zeroth and first moment of the grayscale histograms.
The binarization automatically calculates the optimal threshold at level k value by using discriminant criterion measures. Pixels are divided into two classes C 0 and C 1 with levels [1, 2, . . . , k] and [k + 1, k + 2, . . . , L], respectively, where L is the grayscale level.
The probability of occurrence ω (4), class means levels µ (5), and class variance σ (6) are calculated from the histogram of the pixels values distribution for the two classes. Then, to evaluate the threshold goodness, the class separability is evaluated using specific discriminant criteria. The adopted criterium η (2) is defined as the product between the probability of occurrence and the first-and second-order cumulative moments of the histogram. The aim is to find the threshold that minimizes η (2), or equivalent, the weighted within class variance given by σ B (3): where Otsu's method can be easily implemented on the denoised images by means of a Python OpenCv function. In the present work, the binary image is obtained by settings the values of the minimum and maximum pixels based on the gray levels and to be assigned, respectively, to classes C 0 and C 1 . The described classification relies on the founded global threshold value which may therefore strongly influence the position of the satellite derived shoreline. Threshold values estimated using a global thresholding method can be, indeed, arbitrary values in the middle of two peaks. Kittler and Illingworth [74] identified a good performance of the Otsu's method when the pixel distribution histogram is bimodal with a sharp valley. In this study, we assess the impact of the size of the despeckling kernel size on the pixel distribution by analyzing the pixel distribution of the despeckled image.

Morphological Operation
Two morphological operations are then applied to the binary image to distinguish meaningful shape information (water/land) from irrelevant elements present in the satellite image, such as ships or small isolated areas. Morphological operations are widely used in image processing as they clarify the image by eliminating irrelevancies and by preserving objects shape.
Mathematical morphologies are based on set theory [75,76]. Operations such as union, intersection, translation, and subtraction can be performed on each image. Dilation and erosion processes are two of the basic morphological operations. Dilation combines two vectors using addition, whereas erosion uses subtraction [76]. Applying iterative dilation and erosion, it is possible to remove details smaller than the kernel size without causing geometric distortion of the objects [77]. The composition of dilation and erosion results in opening and closing processes. Opening smooths contours by removing sharp peaks and small areas, whereas closing operation helps in filling holes and fusing narrow breaks [78]. If A is the binary image and B is the structuring elements (both set in Euclidean N-space), the opening of A by B is defined as follows: whereas the closing of A by B is defined as follow: where Dilation (D) and Erosion (E) are defined respectively as: The defined morphological operations are employed in the method here proposed by means of a structuring element with a size 7 × 7. The binary image obtained from Otsu's method is used as input in the opening operation. The resulting image is then used as input for the closing operation.

Canny Edge Detector
The recognition of the objects and their boundaries is obtained through Edge Detection Techniques (EDTs). EDTs convolve the image with an operator sensitive to discontinuities in grayscale level from one pixel to another.
Operators mainly belong to a gradient-based or Laplacian-based algorithm. Gradient methods look for the maximum and minimum of the first derivative of the image. The Laplacian-based methods search, instead, for the zero crossing of the second derivate of the image. Based on its geometry, each operator recognizes edges along a fixed direction (vertical, horizontal, or diagonal). Local operators can be used to estimate vertical and horizontal edges. Among the local operators, the most used are the Sobel's and Prewitt's ones that consist of two 3 × 3 convolution kernels. However, the latter as well as the Laplacian operators are sensitive to image noises [79]. To overcome these limits, Marr-Hildreth [80] proposed an algorithm based on the convolution of the image with the Laplacian of the Gaussian function able to smooth the image and, consequently, reduce the noises influence on the edges extraction.
The Canny edge detector algorithm [81] also allows for the reduction of noises, localization of edges, and suppressing of false edges. The Canny edge detector includes three steps: smoothing, differentiation, and labeling, and shows good performance when compared to other filters [82]. The Canny method recognizes edges by convolving the image with an operator which is the first derivative of a two-dimensional Gaussian. Edges are determined by an optimization process and marked at the maxima in gradient magnitude of a Gaussian-smoothed image [81].
In the present work, in the beginning, the Sobel operator with a default size of 3 × 3 is used to find gradient values in both x and y direction (11). If A is the binary image set in Euclidean N-space, the gradients in each direction are defined by the convolution of A and Sobel operators (S x , S y ). The edge gradient magnitude (G) and direction (θ) are then found as in (12,13): where: Then, by using a minimum and a maximum threshold values, here set, respectively, as 100 and 200, pixels are treated as edges if their intensity gradient is higher than the maximum threshold and discarded if lower. Figure 5 shows the outcomes of each of the aforementioned methodology steps.  Figure 5a shows an example of the input images used for testing the proposed methodology. At first, the salt-pepper noise effect visible in the image is reduced by a medianfiltering approach executed with a kernel of 5 × 5 size, thus allowing for a better recognition of the water/land boundary (Figure 5b). The image binarization, namely, the classification of the pixels for water and land, is automatically executed on the denoised image by means of the Otsu's method which allows for the calculation of the optimal threshold value for water/land recognition. The output is a binary image (Figure 5c) where 0 is assigned to the pixels below the threshold, namely, water, whereas 1 is assigned to those pixels in exceedance, corresponding to land.

Sentinel-1 Retrieved Shoreline
It is observed that the application of both opening ( Figure 5d) and closing (Figure 5e) operations contributes in reducing uncertainties and improving the edges detection. Figure 5e reports an example of the effects of such operations, showing the reduction of the number of black dots visible in the mainland comparing with Figure 5c. Figure 5f shows that by applying the Canny edge filter, edges in the binary image are successfully recognized. The pixels ranging between the two thresholds are considered edges only if connected to pixels above the maximum threshold. Figure 5f also shows that the Canny method ensures the filtering of other small pixels noises, under the assumption that edges are long lines. The maximum distance is equal to 44.74 m, which corresponds to an accuracy of 3.51 pixels (Table 1).

The Effect of the Despeckling Kernel Size on the Shoreline Position
In order to study the effect of different despeckling kernel sizes on the final threshold based on which the shoreline position is defined, different kernel sizes have been applied to the same images and the variation of the threshold has been evaluated. Figure 7 shows an example of the variation of the pixels distribution histogram obtained by increasing the kernel size from 5 × 5 to 25 × 25. It can be observed that the distribution is characterized by a bimodal distribution, which varies with the kernel size. In fact, by replacing the pixel with the median of it and its adjacent neighboring values, the despeckling filters influence the histogram distribution and, thus, the threshold. Particularly, it can be seen that by increasing the kernel size, two different picks and a sharp valley are better recognized. The calculated threshold value for each kernel size varies in a range of 10 pixels. In the meantime, as shown in Figure 8, it is observed that the spatial resolution decreases as despeckling kernel size increases. Smaller kernels preserve better spatial resolution, even if the roughness of the extracted shoreline increases and several edges far from the shoreline are detected close to the coasts. On the contrary, the use of bigger kernel sizes results in a smooth and discontinuous line. For this reason and considering that the final shoreline position does not change by varying threshold on such a range of 10 pixels, the 5 × 5 structural element was chosen as the optimal kernel to reduce the speckle effects. Result obtained by using a median-filtering approach by varying the kernel size in order to reduce speckle noise: (a) kernel size 3 × 3; (b) kernel size 5 × 5; (c) kernel size 7 × 7; (d) kernel size 9 × 9.

Versatility of the Presented Methodology
The validated method is applied on other test sites to investigate its versatility.

Torre Lapillo
Torre Lapillo, hamlet of Porto Cesareo (Lecce), located in Apulian region facing the Ionian Sea (see Figure 9), is selected for its interesting morphology. It is constituted by the typical sub-environment of low-lying coasts and narrow sandy beaches with the emerged beach width typically ranging from less than 10 m to about 50 m. As shown in Figure 10a, the morphology of Torre Lapillo results in a limiting factor for the method. Within the study area a basin, called "Bacino dell'Uomo Morto", is situated just behind the dunes. The basin has a maximum width of about 190 m. Beaches width comparable to pixel spacing coupled with the proximity of the basin to the sea leads to some difficulties in detecting the shoreline position. However, the use of VV polarization on acquired images enables to overpass such a problem (Figure 10b).

Aerial Images: Giovinazzo
The algorithm is also applied to optical aerial images to test its versatility on different images. Ortho-images are photogrammetric images obtained from cameras fixed on airplanes. Several images are captured through photogrammetric flights. The resultant image is then retrieved by overlapping pictures in the orders of 60% forward and 20-40% in lateral direction. Freely ortho-images with a spatial resolution of 0.5 × 0.5 m and 0.2 × 0.2 m available as Web Map Services (WMS) for Apulian region (http://www.sit. puglia.it/portal/portale_cartografie_tecniche_tematiche/WMS, accessed on 25 May 2021) are eligible for our purpose.
The method shows good performances when applied to aerial images. Figures 11  and 12a show the shoreline extracted along Giovinazzo (nothern Apulia) beaches where breakwaters, rocky, and sand beaches come in succession. Nevertheless, beach slopes, shadows, and wave foam may constitute limiting factors for using this method on images acquired by using optical sensors as shown in Figure 12b.

Discussion
In this study, the feasibility of an automated method for shoreline position delineation is demonstrated. The data sets used are freely available and easily accessible, making the methodology reproducible in any geographical area; the methodology yields satisfactory results when applied to Sentinel-1 data. The number of processing steps applied to the images is relatively small, making the method relatively straightforward and allowing for a fully automated procedure capable of long-term shoreline evolution assessment. Sentinel-1 satellite-derived shorelines have been extracted for the year 2017. The retrieved shorelines have been validated against video monitoring system-derived shorelines at Torre Canne beaches, and deterministic verification measure such as RMSE have been computed. Results are very promising with an average RMSE of 12.48 m reported, indicating a sub-pixel accuracy. The minimum distance suggests a good overlapping between video monitoring and satellite based shorelines, except for the 6th of October, during which the maximum distance between the two shorelines also occurs. This suggests a lower geolocation accuracy, and thus the importance of accurate georeferenced and orthorectified products for the model.
Tests have been executed in order to evaluate the influence of the kernel size on pixel distribution, and thus on the optimum threshold for water/land classification. Significant improvement of the shoreline position has not been observed; on the contrary a degradation in the spatial resolution has been observed suggesting a size of 5 × 5 as optimum for the kernel.
The study also attempted to explore the versatility of the proposed methodology. To this end, Torre Lapillo, situated along the Apulian coast, has been selected for its peculiar morphology as constituted by narrow sandy beaches with emerged beach width typically ranging from from less than 10 m to about 50 m, comparable to the pixel size, and the presence of a basin behind the dunes. As shown in Figure 10a, the morphology of the area may represent a limiting factor of the methodology. However, imagery collected using a different polarization combination can help overcome this issue (see Figure 10b) allowing a more accurate water/land classification [83].
Depending on the transmitted and received polarizations, the radiation interacts with and is back-scattered differently by the surface. Therefore, different polarizations cause variations in radar response, providing different and complementary information about the targets in the area. This suggests that the combination of different polarization channels can be a potential candidate for further investigation.
A first attempt to employ the methodology on optical images has also been successfully performed. Ortho-images with a spatial resolution of 20 cm have been employed for this purpose. Preliminary results are promising but also show that parameters such as wave foam, shadows, and shallow beaches can influence the retrieved edges and thus the shoreline position. Additional filters may help overcome this issue and need to be further explored.
Areas of methodological improvements include utilizing segmentation algorithms such as locally adaptive thresholding algorithms which can enhance the land/water boundary recognition and thus reduce discontinuity of coastal edges that can occur in low contrast areas in the image [48]. The proposed method of shoreline detection might face difficulties within inter-tidal areas, where porous medium is characterized by a higher saturation degree which induces some uncertainties in the detection of the shoreline. Indeed, the higher the moisture content is, the higher is the probability that its reflectance behavior is assimilated to the water one. In fact, at saturation, the optical path length in water is at its maximum and certain wavelengths may be completely absorbed [84]. For these areas, a multilevel Otsu thresholding algorithm can be used, and three classes may be defined as suggested in [85]. It is therefore recommended to further explore the method performances in different type of coastal landform and different hydrodynamic conditions, and to analyze its dependence on the tide and wave run up and on the different grain sizes of sandy beaches.
Small changes of the coastal areas over short periods may be difficult to detect due to the coarse spatial resolution often obtained via satellite imagery. To enhance the spatial resolution of the employed satellite images, image fusion techniques may be explored. Shoreline mapping can be indeed performed by using a single satellite sensor, as done in this study, or by the combination of multiple images derived by different sensors. Image fusion allows for the combination of different images to enhance the classification accuracy [86][87][88] and has been applied to shoreline mapping studies and to identify dynamic shoreline trends [89,90].
Additional significant improvements of the spatial and temporal resolution of the satellite imagery are also promised by emerging micro-and nanosatellite constellations (PACE, Geosynchronous Littoral Imaging and Monitoring Radiometer (GLIMR)) and commercial satellites (e.g., Radarsat, WorldView, and Pleiades). Those promise to compensate for the limit imposed by the currently implemented satellite's temporal and spatial resolution providing data at a daily temporal resolution and at a higher spatial resolution. The capability of those satellites to retrieve the shoreline position has already been proved [91][92][93].

Conclusions
This paper presents a new methodology for shoreline detection, based on edge detection techniques of preprocessed SAR images, freely available from ESA and the Copernicus program. The developed procedure of land/water boundary recognition is carried out with a script, which allows to generate both shorelines and error statistics given the availability of in situ data. When compared to the outcomes from in situ and monitoring programs, the obtained shorelines provide reliable results given the spatial resolution of the satellite data. Moreover, in the case of calm sea and high contrast in the imagery, the method shows good performances in detecting shoreline from aerial images.
Long-term variability of shoreline trends is widely used in risk management and planning activities on coastal areas and represents one of the most used indicators when accounting for assessing defense strategies for preventing beach erosion or flooding. In such a context, the availability of data with high both temporal and spatial resolution is crucial. Shoreline detection from satellite images could represent a valid alternative to the traditional field survey. The detection method proposed in the present work has the advantage of using freely-available SAR images which are processed automatically to retrieve shoreline position. If applied continuously, it could allow for effective and timely monitoring of long-term changes of coastal area, supporting decision-makers and administrations to manage defense interventions against beach erosion and flooding. The free and ready availability of Sentinel-1 images, the functionality of the algorithm, its operationalization, and its computational efficiency may indeed represent an efficient instrument to store and process a large number of images and data useful for monitoring coastal zones and detecting erosion, assisting in the reporting efforts of EU member countries, and potentially for identifying flooding events. Acknowledgments: Special thanks to Giorgio Santinelli, colleague at Deltares, for his support at various stage of this project.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: