Di ﬀ erential SAR Interferometry Using Sentinel-1 Imagery-Limitations in Monitoring Fast Moving Landslides: The Case Study of Cyprus

: Cyprus, being located on the Mediterranean fault zone, exhibits a unique geodynamic regime since its tectonic evolution is driven by the interaction of the Eurasian and African plate. Besides its seismological interest, many active landslides and slope instabilities in areas of steep topography occur in Cyprus, having substantial impact on the built environment, by posing an imminent threat for entire settlements and critical infrastructure. Moreover, extreme meteorological events occur rarely, like severe rainfall and thunderstorms, that combined with the geological properties in some areas and the seismically stressed ground, can lead to landslides, causing severe damages to critical infrastructure. In the present study, the DInSAR methodology is applied for the detection of two individual landslide events that were triggered by heavy rainfall in Limassol and Paphos Districts in February 2019. Six co-event interferometric Synthetic Aperture Radar (SAR) pairs were used to produce displacement maps in vertical and east-west directions to study the resulting slope deformations. The above are carried out using Sentinel-1 imagery that are freely provided under the Copernicus umbrella. The limitations that arise from the speed and complexity of the deformations under study and the adverse residing meteorological conditions that caused these phenomena are investigated, as found in literature. Indeed, the sparse vegetation at the slopes a ﬀ ected by the landslides, the residing meteorological conditions, the heavy rainfall that triggered the two landslides, and the temporal phase aliasing e ﬀ ect due to the speed of the ground deformation were found to be the main limitations for the application of DInSAR methodology, resulting in the underestimation of the ground deformation that occurred. displacement fast-moving descending images for the case study near Pissouri, dated 5–17 February 2019 and 11–23 February 2019 respectively, there is a clear phase shift at the west / northwest part of the area a ﬀ ected by the landslide. Author Contributions: Conceptualization, M.T.; methodology, M.T., C.D. and D.G.H.; software, M.T.; validation, M.T. and C.D.; formal analysis, M.T., C.D. and D.G.H.; investigation, M.T.; resources, M.T., C.D. and D.G.H.; data curation, M.T. and C.D.; writing—original draft preparation, M.T.; writing—review and editing, M.T., C.D. and D.G.H.; visualization, M.T.; supervision, C.D. and D.G.H.


Introduction
Landslides can be caused by geological factors (soil stability, lithology, faults, etc.); geomorphological factors, such as soil erosion in the case of Cyprus, vegetation removal due to fire or draught and volcanic or tectonic uplift; physical factors, such as topography, tectonic activity, intense rainfall and snow melt; and factors associated with human activity, such as excavations, agricultural activities, mining and deforestation [1]. The area geology plays an important role on the behavior of land during a seismic event or a landslide. This tendency, especially in Cyprus, is intensified by the long history of powerful seismic activity in the region.
Cyprus is in the earthquake zone of the Alps-Himalayas, where 15% of earthquakes occur worldwide, and exhibits a unique geodynamic regime since its tectonic evolution is driven by the interaction of the Eurasian and the African plate. In fact, 19 earthquakes of MS >5 occurred in Cyprus Moreover, Cyprus has a semi-arid climate type (BSh), in the north-eastern part of island, and a hot Mediterranean climate type (Csa), in the rest of the island according to Köppen-Geiger climate classification, with very mild winters and warm to hot summers. Summers are generally dry with rain occurring mostly between November and March [4,5]. In fact, based on data obtained from the Department of Meteorology in Cyprus, average precipitation is approximately 470 mm/year. Extreme meteorological events occur rarely in Cyprus, such as heavy rainfalls and severe thunderstorms.
These extraordinary events combined with the geological properties in some areas and the seismically stressed ground, can lead to ground deformation, in terms of land subsidence/uplift or landslides [6], causing severe damages to critical infrastructure and urban areas. The most severe landslides have occurred in the mountainous and hilly areas of Paphos and Limassol district, causing damages to residences; civil infrastructure networks, such as roads and road networks; utility services pipelines, etc. Considering the impact of these phenomena combined with damages due to seismic events, villages such as Statos, Agios Photios, Kivides, and Korfi have been abandoned and relocated at safer areas [6].
There are many techniques for monitoring ground deformation due to landslides. These techniques can be classified as ground-based or satellite-based geodetic, geotechnical or physical, photogrammetric, and Earth Observation (EO)-based techniques [7]. For each one of these categories, different instrumentation and equipment needs to be installed or used for periodic monitoring or for continuous data collection either remotely or in-situ [8].
Area inaccessibility, bad weather conditions and high infrastructure cost can be limiting factors for the use of ground-based techniques. The spatial and temporal characteristics of the active satellite sensors used in remote sensing, with their very short revisit times and high coverage, can overcome the obstacles that can hinder the monitoring of natural hazards via conventional means. Earth Observation, through the exploitation of SAR images that can be acquired day and night and regardless of weather conditions [9], can show the way for a systematic monitoring of land Moreover, Cyprus has a semi-arid climate type (BSh), in the north-eastern part of island, and a hot Mediterranean climate type (Csa), in the rest of the island according to Köppen-Geiger climate classification, with very mild winters and warm to hot summers. Summers are generally dry with rain occurring mostly between November and March [4,5]. In fact, based on data obtained from the Department of Meteorology in Cyprus, average precipitation is approximately 470 mm/year. Extreme meteorological events occur rarely in Cyprus, such as heavy rainfalls and severe thunderstorms.
These extraordinary events combined with the geological properties in some areas and the seismically stressed ground, can lead to ground deformation, in terms of land subsidence/uplift or landslides [6], causing severe damages to critical infrastructure and urban areas. The most severe landslides have occurred in the mountainous and hilly areas of Paphos and Limassol district, causing damages to residences; civil infrastructure networks, such as roads and road networks; utility services pipelines, etc. Considering the impact of these phenomena combined with damages due to seismic events, villages such as Statos, Agios Photios, Kivides, and Korfi have been abandoned and relocated at safer areas [6].
There are many techniques for monitoring ground deformation due to landslides. These techniques can be classified as ground-based or satellite-based geodetic, geotechnical or physical, photogrammetric, and Earth Observation (EO)-based techniques [7]. For each one of these categories, different instrumentation and equipment needs to be installed or used for periodic monitoring or for continuous data collection either remotely or in-situ [8].
Area inaccessibility, bad weather conditions and high infrastructure cost can be limiting factors for the use of ground-based techniques. The spatial and temporal characteristics of the active satellite sensors used in remote sensing, with their very short revisit times and high coverage, can overcome the obstacles that can hinder the monitoring of natural hazards via conventional means. Earth Observation, through the exploitation of SAR images that can be acquired day and night and regardless of weather conditions [9], can show the way for a systematic monitoring of land displacement and assessment of the impact of natural hazards on infrastructure resilience over large areas in a timely manner.
Differential SAR Interferometry (DInSAR) is an active remote sensing technique based on the principle that, due to the very high stability of the satellite orbits, it is possible to exploit the informative contribution carried by the phase difference between two or more SAR images looking at the same scene from comparable geometries [10][11][12]. The availability of a large set of data and the multi-pass DInSAR approach leads to achieving a more accurate estimation of atmospheric contribution and topographic phase component and allows the reduction of error sources [13,14]. In addition, the dense data temporal sampling makes it possible to study the temporal evolution of a phenomenon in order to provide time series and rate of the deformation process.
DInSAR can be used to detect and monitor, under certain conditions, ground deformations on a centimetric scale, with a good spatial resolution (a few meters), offering an effective method for their assessment over large areas [15]. These characteristics are making DInSAR one of the most promising space-geodetic techniques for monitoring Earth's surface deformations [16]. Therefore, DInSAR techniques can be used for various engineering applications, such as land use and urban planning, as they contribute to the better understanding of an area's geomorphology and assist in the impact assessment of natural hazards, such as earthquakes and landslides [17,18].
DInSAR has been used successfully in numerous applications for the detection of the impact of natural hazards in terms of ground deformation caused by earthquakes [19][20][21][22][23] and landslides [24][25][26]. Moreover, DInSAR has proved to be performing efficiently in the detection of damages in civil infrastructure like road networks and dams [27,28].
Based on these studies, the selection of the type of SAR data should be case-specific and should always be based on the properties and characteristics of the area of interest and its surroundings to achieve optimum results. Therefore, SAR acquisitions with a long wavelength, such as L-band and P-band, having enhanced capabilities for penetration, are recommended for areas covered with vegetation [29,30]. In the cases of DInSAR and MT-InSAR, SAR data with shorter wavelengths such as the C-band and X-band are preferred, as there is a need to improve the accuracy of the calculated deformations, however, having certain limitations due to their sensitivity to vegetation coverage and temporal decorrelation. In the case of early warning systems and detailed monitoring, high-resolution SAR data are required [31]. Integration of different SAR interferometry techniques can overcome any obstacles by combining their respective advantages [32].
However, several criteria have to be met, as described before, the most important being the maintenance of coherence, at some extent, during the period of study [33]. Nevertheless, the low spatial resolution characteristics of some satellite imagery impose some limitations. Additionally, in some cases the temporal resolution is also very limited. The increased spatial resolution characteristics of the TerraSAR-X and COSMO-SkyMed satellites and the temporal resolution, as well as the free, full, and open data policy of the Copernicus program, that the Sentinel-1 constellation offers, can be considered as milestones for the scientific progress of radar space technology [34]. The launch of Sentinel-1A in 2014 and Sentinel-1B 2 years later, has made possible the systematic supply of data for multitemporal analyses.
Indeed, the temporal characteristics of the Sentinel-1 mission, such as the 12-day repeat cycle of the Sentinel-1 sensors in a single pass (6-day repeat cycle of the two-satellite constellation at the equator), ascending or descending, and the spatial resolution of Level-1 Single Look Complex (SLC) products can be used for mapping rapid landscape changes. The fact that these Copernicus Sentinel-1 data are freely distributed datasets, provided by the European Space Agency (ESA), along with the availability of open-source software for processing these images, provides us with the opportunity to systematically monitor and assess the impact of natural hazards to built environment worldwide, overcoming area inaccessibility issues and weather conditions, even in regions and countries with limited budgets [19,23,35].
Additionally, although the moderate spatial resolution of the IW swath images for Level-1 SLC products (2.7 × 22 m to 3.5 × 22 m), which are available for Cyprus, has proven to be satisfactory for the detection of land subsidence and/or uplift [22,23,36], the detailed study and monitoring of more complex ground displacement phenomena, especially of fast-moving landslides [1,37,38] and slope deformations, cannot be achieved using Sentinel-1 images due to the limitations of the C-band wavelength and the images' medium spatial resolution [39][40][41].
In the present paper, freely available Sentinel-1 imagery and the free Sentinel's Application Platform (SNAP) software developed by ESA were used. The objective of the present study is the investigation of the limitations in the application of DInSAR methodology for the study of two landslides that occurred in Cyprus in February 2019 due to heavy rainfall, arising from: • the spatial characteristics of the Sentinel-1 datasets, • the speed of the ground deformation, • the residing meteorological conditions and • vegetation In the next section, the areas of study are presented, followed by the proposed methodology and the materials used for its application. In Section 4, the results are presented, and finally, the results from the application of the proposed methodology are discussed and the future steps for overcoming the aforementioned limitations are highlighted.

Area of Study
In the present study, two landslide events that took place at the south-west part of the island, and more specifically in Limassol and Paphos District, were studied, as presented in Figures 2 and 3 respectively. Both landslides affected thousands of people commuting every day between the two districts, using the specific parts of the road network. The areas affected by the two landslides were defined based on directions by the Geological Survey Department of Cyprus.
Geosciences 2020, 10, x FOR PEER REVIEW 4 of 26 complex ground displacement phenomena, especially of fast-moving landslides [1,37,38] and slope deformations, cannot be achieved using Sentinel-1 images due to the limitations of the C-band wavelength and the images' medium spatial resolution [39][40][41].
In the present paper, freely available Sentinel-1 imagery and the free Sentinel's Application Platform (SNAP) software developed by ESA were used. The objective of the present study is the investigation of the limitations in the application of DInSAR methodology for the study of two landslides that occurred in Cyprus in February 2019 due to heavy rainfall, arising from:  the spatial characteristics of the Sentinel-1 datasets,  the speed of the ground deformation,  the residing meteorological conditions and  vegetation In the next section, the areas of study are presented, followed by the proposed methodology and the materials used for its application. In Section 4, the results are presented, and finally, the results from the application of the proposed methodology are discussed and the future steps for overcoming the aforementioned limitations are highlighted.

Area of Study
In the present study, two landslide events that took place at the south-west part of the island, and more specifically in Limassol and Paphos District, were studied, as presented in Figures 2 and 3 respectively. Both landslides affected thousands of people commuting every day between the two districts, using the specific parts of the road network. The areas affected by the two landslides were defined based on directions by the Geological Survey Department of Cyprus.  The first landslide took place on 15 February 2019, following substantial rainfall, next to the Paphos-Limassol motorway A6 (Figure 2), in the direction from Paphos to Limassol between exits of Avdimou and Pissouri. Extensive earthworks were carried out to remove the disturbed soil and the Paphos-Limassol motorway was again open to traffic on 23 March 2019.
Moreover, a landslide that took place on 20 February 2019, following heavy rainfall, was studied ( Figure 3). A 3-km-long part of the old road from Paphos to Limassol was closed, approximately 200 m away from "Petra to Romiou". Structural issues were found in that part of the road, which remains closed until today for security reasons. The access from Paphos to Limassol is still carried out only through the Paphos-Limassol motorway A6. The extensive precipitation and residing soil erosion issues were found to be the main reasons for this phenomenon. In this region, the soils are particularly problematic and prone to landslides, rockslides and other geotechnical related failures. The wider area's geological formations consist of bentonite, lava, limestone, quartz sandstone, argillaceous shale and hornstone, serpentinite, pyroxenite, gabbro, chalks and marls, rocks with different mechanical characteristics in terms of weathering and erosion resistance [42]. The overall geology of the case study area is shown in Figure 4. Moreover, a landslide that took place on 20 February 2019, following heavy rainfall, was studied ( Figure 3). A 3-km-long part of the old road from Paphos to Limassol was closed, approximately 200 m away from "Petra to Romiou". Structural issues were found in that part of the road, which remains closed until today for security reasons. The access from Paphos to Limassol is still carried out only through the Paphos-Limassol motorway A6. The extensive precipitation and residing soil erosion issues were found to be the main reasons for this phenomenon. In this region, the soils are particularly problematic and prone to landslides, rockslides and other geotechnical related failures. The wider area's geological formations consist of bentonite, lava, limestone, quartz sandstone, argillaceous shale and hornstone, serpentinite, pyroxenite, gabbro, chalks and marls, rocks with different mechanical characteristics in terms of weathering and erosion resistance [42]. The overall geology of the case study area is shown in Figure 4. Geosciences 2020, 10, x FOR PEER REVIEW 6 of 26 Pissouri, located to the south and southeast of the case study area a, as seen in Figure 4, is covered mostly by calcarenites, carbonate marls and sandstones of the Nicosia formation, as in the case of the case study area near A6 motorway [43]. Underneath lies the Pissouri marl, which is characterized by sandy interlayers. The marlstone is susceptible to fast weathering and in most sites the fresh marl is covered by a layer of weathered marl [44]. The wider area of Petra tou Romiou, near the village of Kouklia, lies on a complex geology with the formations of the Mamonia Complex, the Troodos Ophiolite Complex, the Kannaviou Formation, the Lefkara Formation, the Pachna Formation, as well as the Plio-Pleistocene deep-sea and alluvial depositions. The limestone of Petra tou Romiou is exposed in detached blocks of various sizes, as the Fasoulla lavas of the Fasoulla formation [42,45]. The case study area located near "Petra tou Romiou", i.e., area b, and the surrounding area between the villages of Kouklia and Pissouri lies mainly on the Pachna Formation, which consists of yellowish marls and chalks [45].
Photographs showing the magnitude and the impact of the two landslide events that are under investigation were provided by the Cyprus Geological Survey Department, as presented in Figure 5. Both landslides took place overnight, resulting in substantial ground displacement. Figure 5a presents the impact of the landslide that took place next to Paphos-Limassol motorway near Pissouri on the 15th of February 2019 and Figure 5b shows the impact of the landslide that occurred on the 20th of February next to Paphos-Limassol old coastal road near to Petra tou Romiou. Pissouri, located to the south and southeast of the case study area a, as seen in Figure 4, is covered mostly by calcarenites, carbonate marls and sandstones of the Nicosia formation, as in the case of the case study area near A6 motorway [43]. Underneath lies the Pissouri marl, which is characterized by sandy interlayers. The marlstone is susceptible to fast weathering and in most sites the fresh marl is covered by a layer of weathered marl [44]. The wider area of Petra tou Romiou, near the village of Kouklia, lies on a complex geology with the formations of the Mamonia Complex, the Troodos Ophiolite Complex, the Kannaviou Formation, the Lefkara Formation, the Pachna Formation, as well as the Plio-Pleistocene deep-sea and alluvial depositions. The limestone of Petra tou Romiou is exposed in detached blocks of various sizes, as the Fasoulla lavas of the Fasoulla formation [42,45]. The case study area located near "Petra tou Romiou", i.e., area b, and the surrounding area between the villages of Kouklia and Pissouri lies mainly on the Pachna Formation, which consists of yellowish marls and chalks [45].
Photographs showing the magnitude and the impact of the two landslide events that are under investigation were provided by the Cyprus Geological Survey Department, as presented in Figure 5. Both landslides took place overnight, resulting in substantial ground displacement. Figure 5a presents the impact of the landslide that took place next to Paphos-Limassol motorway near Pissouri on the 15th of February 2019 and Figure 5b shows the impact of the landslide that occurred on the 20th of February next to Paphos-Limassol old coastal road near to Petra tou Romiou.

Materials and Methods
For the needs of the interferometric analysis, a high coherence between the two image acquisitions is very important. Coherence is used as a quantitative indicator of the quality of the associated interferogram and the DInSAR findings and is an indication of any significant radar signal changes between acquisitions. Changes in the surface reflectivity over time due to vegetation or other natural phenomena often exhibit and introduce sudden changes in coherence in SAR imagery, decorrelating the phase measurements significantly [47]. Thus, for higher accuracy of results, it is suggested that the final displacement maps are masked according to high coherence for land cover such as urban areas.
Therefore, it is recommended that the acquisition times for the pair of SAR images are temporally close in order to minimize the loss of coherence. Small perpendicular baseline interferograms are preferable, as these do not suffer from geometric decorrelation, show good fringes, and are usually easy to unwrap [10]. For DInSAR, the limit of their perpendicular baselines was set to 200 m [48]. Due to the high revisit frequency and the excellent orbital accuracy of Sentinel-1, the perpendicular baseline is usually not a limiting factor, facilitating the application of InSAR in many countries [40]. Interferometric phase is also influenced by other factors, as presented in Equation (1), such as the atmosphere, pixel reflectivity, radar system phase noise, and the uncertainty in knowing the platform position and baseline. topo disp atm pixel noise error Δφ = Δφ + Δφ + Δφ + Δφ + Δφ + Δφ (1) The error associated with the atmospheric effect can be corrected by stacking or averaging many interferograms, since the atmospheric effect has short spatial correlation length with virtually no temporal correlation between data acquisition times [13,14,49].

Materials and Methods
For the needs of the interferometric analysis, a high coherence between the two image acquisitions is very important. Coherence is used as a quantitative indicator of the quality of the associated interferogram and the DInSAR findings and is an indication of any significant radar signal changes between acquisitions. Changes in the surface reflectivity over time due to vegetation or other natural phenomena often exhibit and introduce sudden changes in coherence in SAR imagery, decorrelating the phase measurements significantly [47]. Thus, for higher accuracy of results, it is suggested that the final displacement maps are masked according to high coherence for land cover such as urban areas.
Therefore, it is recommended that the acquisition times for the pair of SAR images are temporally close in order to minimize the loss of coherence. Small perpendicular baseline interferograms are preferable, as these do not suffer from geometric decorrelation, show good fringes, and are usually easy to unwrap [10]. For DInSAR, the limit of their perpendicular baselines was set to 200 m [48]. Due to the high revisit frequency and the excellent orbital accuracy of Sentinel-1, the perpendicular baseline is usually not a limiting factor, facilitating the application of InSAR in many countries [40]. Interferometric phase is also influenced by other factors, as presented in Equation (1), such as the atmosphere, pixel reflectivity, radar system phase noise, and the uncertainty in knowing the platform position and baseline.
Geosciences 2020, 10, 236 The error associated with the atmospheric effect can be corrected by stacking or averaging many interferograms, since the atmospheric effect has short spatial correlation length with virtually no temporal correlation between data acquisition times [13,14,49].
Moreover, a limitation is found in literature, called temporal phase aliasing, which mainly affects the phase unwrapping step of the DInSAR methodology [41,50], when the displacement of a pixel or an area under investigation exceeds the value of λ/4 between two consecutive satellite acquisitions that form an interferometric SAR pair [51,52]. λ is the wavelength of the SAR sensor, which in the case of Sentinel-1 satellite's C-band is approximately 5.55 cm. This is due to the fact that InSAR phase signals are 2π-wrapped, making it impossible to measure phase differences exceeding π unambiguously, that in turn correspond to LoS displacements over λ/4 [51,53]. Therefore, deformations that can be two, three or n times greater than the λ/4 threshold can be underestimated as they produce similar observed phase signals [41,51,53]. Temporal phase aliasing thresholds for the currently available SAR satellites have been calculated based on their nominal wavelengths and revisit times [41,[51][52][53].
For the investigation of the two landslides, 10 Sentinel-1 (6 Sentinel-1A and 4 Sentinel-1B), Interferometric Wide (IW) swath SLC images of VV + VH polarization were downloaded for the period starting from 4 February 2019 until 1 March 2019, to carry out a comprehensive study on the impact of the landslides to the affected areas. A pre-selection of the downloaded images was carried out prior to the interferometric processing. The acquisition of the SAR images was performed through the Copernicus Open Access Hub by specifying the area of interest for each case.
The processing of Sentinel-1 SAR acquisitions was carried out using the open-source software Sentinel Applications Platform (SNAP) (version 7.0, as of 22nd of July 2019) [54]. Within SNAP, the Sentinel-1 Toolbox (S1TBX) was used for the application of the proposed DInSAR processing chain presented in Figure 6.
Geosciences 2020, 10, x FOR PEER REVIEW  8 of 26 Moreover, a limitation is found in literature, called temporal phase aliasing, which mainly affects the phase unwrapping step of the DInSAR methodology [41,50], when the displacement of a pixel or an area under investigation exceeds the value of λ/4 between two consecutive satellite acquisitions that form an interferometric SAR pair [51,52]. λ is the wavelength of the SAR sensor, which in the case of Sentinel-1 satellite's C-band is approximately 5.55cm. This is due to the fact that InSAR phase signals are 2π-wrapped, making it impossible to measure phase differences exceeding π unambiguously, that in turn correspond to LoS displacements over λ/4 [51,53]. Therefore, deformations that can be two, three or n times greater than the λ/4 threshold can be underestimated as they produce similar observed phase signals [41,51,53]. Temporal phase aliasing thresholds for the currently available SAR satellites have been calculated based on their nominal wavelengths and revisit times [41,[51][52][53].
For the investigation of the two landslides, 10 Sentinel-1 (6 Sentinel-1A and 4 Sentinel-1B), Interferometric Wide (IW) swath SLC images of VV + VH polarization were downloaded for the period starting from 4 February 2019 until 1 March 2019, to carry out a comprehensive study on the impact of the landslides to the affected areas. A pre-selection of the downloaded images was carried out prior to the interferometric processing. The acquisition of the SAR images was performed through the Copernicus Open Access Hub by specifying the area of interest for each case.
The processing of Sentinel-1 SAR acquisitions was carried out using the open-source software Sentinel Applications Platform (SNAP) (version 7.0, as of 22nd of July 2019) [54]. Within SNAP, the Sentinel-1 Toolbox (S1TBX) was used for the application of the proposed DInSAR processing chain presented in Figure 6. Figure 6. The proposed DInSAR processing chain. The steps from data acquisition until the production of the land displacement map are presented. This processing chain is suggested for processing Sentinel-1 using SNAP software, although the same can be applied with different satellite images and processing software.
These images were used to form six co-event pairs, with all images in the respective SAR pairs being acquired 12 days apart, to study the two landslide events that took place near Pissouri and Petra tou Romiou on 15 February 2019 and 20 February 2019 respectively. The suitability of all SAR Figure 6. The proposed DInSAR processing chain. The steps from data acquisition until the production of the land displacement map are presented. This processing chain is suggested for processing Sentinel-1 using SNAP software, although the same can be applied with different satellite images and processing software. These images were used to form six co-event pairs, with all images in the respective SAR pairs being acquired 12 days apart, to study the two landslide events that took place near Pissouri and Petra tou Romiou on 15 February 2019 and 20 February 2019 respectively. The suitability of all SAR image pairs was checked via the criteria of the perpendicular baseline and other image pair characteristics, such as the temporal baseline and the modelled coherence. Perpendicular baseline and modelled coherence values were calculated automatically through the InSAR stack overview tool of SNAP software. Modelled coherence is only an estimate of the coherence of the interferometric SAR pair, which is calculated accurately at the step of interferogram formation of the DInSAR methodology. The characteristics of the co-event SAR pairs are presented in Table 1. The perpendicular baselines for all six co-event SAR pairs were shorter than the threshold of 150 m that was set. Moreover, all six pairs showed high modelled coherence equal or greater than 0.90. However, the modelled coherence shown above is only an estimate. The coherence of the interferometric pair will be calculated after the completion of the step of the interferogram formation. Therefore, all six pairs were processed further using the proposed DInSAR processing chain. Interferometric pairs 1, 2, 3, and 4 were used to study the landslide that occurred on 15 February 2019 near Pissouri, next to the A6 motorway, whereas the interferometric pairs 3, 4, 5, and 6 were used to investigate the impact of the landslide that took place near Petra tou Romiou by the old Paphos-Limassol road.
Sentinel Precise Orbits and the SRTM (Shuttle Radar Topography Mission) 1 sec HGT Digital Elevation Model (DEM) were automatically downloaded through SNAP software and were used at the step of co-registration. Interferograms and coherence maps were then produced for the six interferometric pairs. To produce the final products, the topography component was removed by using the provided DEM, and the black line between the bursts of the image was removed using the Deburst tool. Finally, Goldstein phase filter was applied to enhance the developed interferograms and make the interferometric fringes visually clearer. The final interferograms and coherence maps were produced.
The phase unwrapping step of the Interferometric processing of the multiple SAR images was carried out within the environment of SNAP software using the SNAPHU plugin [55]. The unwrapped interferogram that was produced from the previous step was used to convert phase to Line-of-Sight (LoS) displacement. The Range-Doppler Terrain Correction was then applied in order to geo-locate the displacement map products. In this step, the resampling method of bilinear interpolation and the map projection of WGS84 were selected. Areas without elevation were masked out from the geo-referenced displacement map.
Following the geo-location of all displacement maps, the relative LoS displacement values that were calculated at this stage were converted to absolute values by specifying a "stable" point with high coherence in the area, in our case, a nearby building. The pixel value of displacement at the coordinates of the building was subtracted from all the displacement map values in all six cases individually. The displacement products were then stacked together using the first interferometric pair according to the acquisition date as the master image to assist in the comparison of results on the same basis.
As very small values of LoS displacement were observed, these were converted into vertical, i.e., up/down, and horizontal, east-west only, components, to investigate the impact of both landslides further. The north-south part of the horizontal displacement requires different processing with various combinations of ascending and descending SAR images [56][57][58]. The calculation of the vertical (up-down) and east-west displacement components can be achieved using the already developed average LoS displacement products from the ascending and descending image pairs, together with the scene's incidence angle through Equations (2) and (3), as seen in Figure 7.
interpolation and the map projection of WGS84 were selected. Areas without elevation were masked out from the geo-referenced displacement map. Following the geo-location of all displacement maps, the relative LoS displacement values that were calculated at this stage were converted to absolute values by specifying a "stable" point with high coherence in the area, in our case, a nearby building. The pixel value of displacement at the coordinates of the building was subtracted from all the displacement map values in all six cases individually. The displacement products were then stacked together using the first interferometric pair according to the acquisition date as the master image to assist in the comparison of results on the same basis.
As very small values of LoS displacement were observed, these were converted into vertical, i.e., up/down, and horizontal, east-west only, components, to investigate the impact of both landslides further. The north-south part of the horizontal displacement requires different processing with various combinations of ascending and descending SAR images [56][57][58]. The calculation of the vertical (up-down) and east-west displacement components can be achieved using the already developed average LoS displacement products from the ascending and descending image pairs, together with the scene's incidence angle through Equations (2) and (3), as seen in Figure 7.
Coherence maps were also developed to study the effect of coherence loss to the displacement results of the DInSAR methodology. To remove the effect of vegetation in coherence, a NDVI (Normalized Difference Vegetation Index) mask was created through Equation (4), using a Sentinel-2 image, taking into account that for NDVI values greater than 0.2, vegetation is present. NDVI values Coherence maps were also developed to study the effect of coherence loss to the displacement results of the DInSAR methodology. To remove the effect of vegetation in coherence, a NDVI (Normalized Difference Vegetation Index) mask was created through Equation (4), using a Sentinel-2 image, taking into account that for NDVI values greater than 0.2, vegetation is present. NDVI values range from −1 to + 1, with the negative values of NDVI generally corresponding to water, whereas NDVI values between −0.1 and + 0.2 usually represent barren areas of rock, sand, or snow [59].

Results
The final interferograms and coherence maps that were produced after the phase filtering step in the proposed processing chain (Figure 6), geo-located for easier visualization, are presented in Figure 8 below. The produced interferograms were then converted into KMZ (Keyhole Markup language Zipped) files to be overlaid in Google Earth to perform a quick check for any fringes that could potentially indicate the rapid change on the Earth's surface in the two areas of interest. The results are presented for the landslide near Pissouri and Petra tou Romiou in Figure 9a,b respectively. The interferograms produced were filtered for phase noise using the Goldstein phase filter in SNAP. The black lines in all interferograms in both figures are the limits of motorway A6 and the coastline, respectively.
The produced interferograms were then converted into KMZ (Keyhole Markup language Zipped) files to be overlaid in Google Earth to perform a quick check for any fringes that could potentially indicate the rapid change on the Earth's surface in the two areas of interest. The results are presented for the landslide near Pissouri and Petra tou Romiou in Figure 9a and 9b respectively. The interferograms produced were filtered for phase noise using the Goldstein phase filter in SNAP. The black lines in all interferograms in both figures are the limits of motorway A6 and the coastline, respectively.
(a) (b) Figure 9. Differential interferograms after the application of Goldstein phase filtering for the landslides (a) near Pissouri and (b) near Petra tou Romiou. In the case of the interferograms based on S1A and S1B descending images for the case study near Pissouri, dated 5-17 February 2019 and 11-23 February 2019 respectively, there is a clear phase shift at the west/northwest part of the area affected by the landslide.
The fringes, or else the color variations, represent phase changes between two images in terms of radians of a full 2π cycle. Each cycle of colors represents half the sensor's wavelength, i.e., in the case of Sentinel-1, it represents a LoS displacement of ~2.7 cm. In Figure 9a, i.e., the case of the 15 th of February landslide near Pissouri, it can be seen that there is substantial phase noise in both ascending S1A and S1B SAR image pairs at the landslide area and the areas nearby. This is also supported by the relatively low coherence values within the Area of Interest (AOI) that range from 0.286 to 0.535 Figure 9. Differential interferograms after the application of Goldstein phase filtering for the landslides (a) near Pissouri and (b) near Petra tou Romiou. In the case of the interferograms based on S1A and S1B descending images for the case study near Pissouri, dated 5-17 February 2019 and 11-23 February 2019 respectively, there is a clear phase shift at the west/northwest part of the area affected by the landslide.
The fringes, or else the color variations, represent phase changes between two images in terms of radians of a full 2π cycle. Each cycle of colors represents half the sensor's wavelength, i.e., in the case of Sentinel-1, it represents a LoS displacement of~2.7 cm. In Figure 9a, i.e., the case of the 15th of February landslide near Pissouri, it can be seen that there is substantial phase noise in both ascending S1A and S1B SAR image pairs at the landslide area and the areas nearby. This is also supported by the relatively low coherence values within the Area of Interest (AOI) that range from 0.286 to 0.535 with the mean coherence being 0.390. The noise at the north of the landslide outlined area could be due to the vegetation and at the south due to the presence of the motorway (Figure 2). Additionally, as analyzed in the introduction, heavy rainfall and changes in soil moisture can alter the radar signal significantly, leading to loss of coherence.
In Figure 9b, i.e., the case of the 20th of February landslide near Petra tou Romiou, there appears to be phase noise in all interferometric SAR image pairs to the west and northwest from the outlined landslide. It is clear from Figure 3 that this is mainly due to the existing vegetation, which can lead to radar signal coherence reduction. There is also significant phase noise to the south after the black line (coastline), as sea is characterized by very low, close to zero, coherence. The six produced interferograms were then unwrapped within SNAP, using the SNAPHU plugin, and the unwrapped results were then re-imported to the DInSAR processing chain to calculate displacement values.

Displacement
Following, the unwrapped interferograms were used to convert phase to LoS displacement, i.e., deformation along the slant range direction of the satellite. All displacement products were then geo-located using the Range-Doppler Terrain Correction function. Finally, the individual displacement products were stacked in chronological order within SNAP software using Geolocation as the initial offset method. The relative LoS displacement values that were calculated at this stage were then converted to absolute values by specifying a "stable" point of high coherence in the area, in our case a nearby building. The pixel values of displacement at the coordinates of the building were subtracted by all the displacement map values in all six cases individually.
In all cases, the extreme displacement values in intense red and blue color appeared at the sea area around Cyprus, and most of the pixel values of displacement inland are around zero, with some small variations seen in light red and blue, for subsidence and uplift, respectively. In the cases where more intense red and blue colors can be seen inland in the displacement maps, this is since zero is close to the confidence intervals of the histogram for 95% level of confidence. The displacement maps produced in SNAP were inserted into Google Earth as KMZ files and the results are presented for the landslide near Pissouri in Figure 10 and for the landslide near Petra tou Romiou in Figure 11.
as analyzed in the introduction, heavy rainfall and changes in soil moisture can alter the radar signal significantly, leading to loss of coherence.
In Figure 9b, i.e., the case of the 20th of February landslide near Petra tou Romiou, there appears to be phase noise in all interferometric SAR image pairs to the west and northwest from the outlined landslide. It is clear from Figure 3 that this is mainly due to the existing vegetation, which can lead to radar signal coherence reduction. There is also significant phase noise to the south after the black line (coastline), as sea is characterized by very low, close to zero, coherence. The six produced interferograms were then unwrapped within SNAP, using the SNAPHU plugin, and the unwrapped results were then re-imported to the DInSAR processing chain to calculate displacement values.

Displacement
Following, the unwrapped interferograms were used to convert phase to LoS displacement, i.e., deformation along the slant range direction of the satellite. All displacement products were then geolocated using the Range-Doppler Terrain Correction function. Finally, the individual displacement products were stacked in chronological order within SNAP software using Geolocation as the initial offset method. The relative LoS displacement values that were calculated at this stage were then converted to absolute values by specifying a "stable" point of high coherence in the area, in our case a nearby building. The pixel values of displacement at the coordinates of the building were subtracted by all the displacement map values in all six cases individually.
In all cases, the extreme displacement values in intense red and blue color appeared at the sea area around Cyprus, and most of the pixel values of displacement inland are around zero, with some small variations seen in light red and blue, for subsidence and uplift, respectively. In the cases where more intense red and blue colors can be seen inland in the displacement maps, this is since zero is close to the confidence intervals of the histogram for 95% level of confidence. The displacement maps produced in SNAP were inserted into Google Earth as KMZ files and the results are presented for the landslide near Pissouri in Figure 10 and for the landslide near Petra tou Romiou in Figure 11.   In (a,c,d), most of the pixel values of displacement are around zero, with some small variations seen in light red and blue, for subsidence and uplift, respectively. In the case of (b), more intense blue colors can be seen in the displacement maps, as zero is close to the upper confidence interval of 0.02 (95% level of confidence). Figure 11. Absolute LoS displacement maps for landslide by the old Limassol-Paphos road. The boundaries of the area affected by the landslide are marked with a red polygon. In (a), (b) and (c), most of the pixel values of displacement are around zero, with some small variations seen in light red and blue, for subsidence and uplift, respectively. In the case of (d), more intense red colors can be seen in the displacement maps, as zero is close to the lower confidence interval of −0.04 (95% level of confidence).
In both sites, zero and/or near to zero LoS displacements were calculated. As stated in the introduction, the averaging or stacking of many interferometric products can assist in minimizing or even correcting the error associated with the atmospheric effect, since the atmospheric effect has short spatial correlation length with virtually no temporal correlation between data acquisition times. The averaging of the geo-located displacement maps for each landslide event was performed by stacking the various displacement products via the Band Maths module in SNAP.
Thus Moreover, the LoS displacement was decomposed into vertical, i.e., up/down and horizontal, east-west only, components, to investigate both events further, through Equations (2) and (3). The average LoS displacement products from the ascending and descending image pairs were used together with the scene's incidence angle. Due to the magnitude of the incident angle, θ = 39.05°, displacement values were expected to increase for the vertical and east-west components but not significantly, since the denominator in both equations was over 0.5, i.e., 0.777 for the vertical and 0.63 for the east-west displacement.
The displacement maps from SNAP were then inserted in ArcGIS, where the sea was removed and displacement values were classified for the development of the final LoS, vertical and east-west displacement products, as presented in Figure 12 for the case study near Pissouri and Figure 13 for the case study near Petra tou Romiou. Figure 11. Absolute LoS displacement maps for landslide by the old Limassol-Paphos road. The boundaries of the area affected by the landslide are marked with a red polygon. In (a-c), most of the pixel values of displacement are around zero, with some small variations seen in light red and blue, for subsidence and uplift, respectively. In the case of (d), more intense red colors can be seen in the displacement maps, as zero is close to the lower confidence interval of −0.04 (95% level of confidence).
In both sites, zero and/or near to zero LoS displacements were calculated. As stated in the introduction, the averaging or stacking of many interferometric products can assist in minimizing or even correcting the error associated with the atmospheric effect, since the atmospheric effect has short spatial correlation length with virtually no temporal correlation between data acquisition times. The averaging of the geo-located displacement maps for each landslide event was performed by stacking the various displacement products via the Band Maths module in SNAP.
Thus Moreover, the LoS displacement was decomposed into vertical, i.e., up/down and horizontal, east-west only, components, to investigate both events further, through Equations (2) and (3). The average LoS displacement products from the ascending and descending image pairs were used together with the scene's incidence angle. Due to the magnitude of the incident angle, θ = 39.05 • , displacement values were expected to increase for the vertical and east-west components but not significantly, since the denominator in both equations was over 0.5, i.e., 0.777 for the vertical and 0.63 for the east-west displacement.
The displacement maps from SNAP were then inserted in ArcGIS, where the sea was removed and displacement values were classified for the development of the final LoS, vertical and east-west displacement products, as presented in Figure 12 for the case study near Pissouri and Figure 13 for the case study near Petra tou Romiou.  The variations in LoS, vertical and East-West displacement values that are observed outside the boundaries of the area that was affected by the landslide, towards the east and west, northwest and southwest directions are mainly due to the residing vegetation and in some cases lie along the road network. Additionally, due to the meteorological conditions of the specific period and the intense precipitation changes in the surface reflectivity could have led to confusing results. Another interesting point is the fact that displacement values outside the AOI seem generally greater than those within the AOI. As seen in Figure 13, the AOI near Petra tou Romiou, which consists of 79 pixels in total, appears to have moved opposite to the slant range direction between 2.4 mm and 5.2 mm (negative values) with the average value being 4.1 mm. Moreover, the AOI appears to have subsided by 5.3 mm on average with the maximum subsidence (negative values) being 6.7 mm and the minimum 3.1 mm. In It appears, from Figure 12, that in the AOI near Pissouri, there is a maximum LoS displacement along the slant range direction of 5.2 mm with a single pixel moving in the opposite direction at 0.3 mm (negative value). The average LoS displacement is approximately 3 mm along the slant range direction. Moreover, a maximum uplift of 6.7 mm is observed, with the average being 3.9 mm (min + 0.9 mm and max + 6.7 mm). In the east-west direction the entire area seems to have moved to the east with the maximum shift being 8 mm and the minimum 1.1 mm (average eastwards + 4.8 mm). The only exception is a small area (a single pixel) at the northeast part of the AOI, which appears to have marginally subsided (−0.4 mm) and moved to the west (−0.4 mm). As the standard deviation for the 28 pixels that the AOI consists of, is ± 1.4 mm for the LoS displacement, ± 1.8 mm for the vertical displacement component and ± 2.2 mm for the East-West displacement, these opposite to the slant direction, subsiding and westbound movements are considered insignificant.
The variations in LoS, vertical and East-West displacement values that are observed outside the boundaries of the area that was affected by the landslide, towards the east and west, northwest and southwest directions are mainly due to the residing vegetation and in some cases lie along the road network. Additionally, due to the meteorological conditions of the specific period and the intense precipitation changes in the surface reflectivity could have led to confusing results. Another interesting point is the fact that displacement values outside the AOI seem generally greater than those within the AOI.
As seen in Figure 13, the AOI near Petra tou Romiou, which consists of 79 pixels in total, appears to have moved opposite to the slant range direction between 2.4 mm and 5.2 mm (negative values) with the average value being 4.1 mm. Moreover, the AOI appears to have subsided by 5.3 mm on average with the maximum subsidence (negative values) being 6.7 mm and the minimum 3.1 mm. In the east-west direction, the AOI appears to have moved to the west by maximum 8.3 mm and minimum 3.9 mm (average 6.5 mm). Standard deviations for the LoS, up-down and east-west displacements were ± 0.6 mm, ± 0.8 mm and ± 0.9 mm, respectively. As in the case of the other study site, displacement values outside the AOI seem generally greater than those within the AOI. This could be mainly due to the vegetation and the precipitation that took place during the specific period.
Based on the photos of the landslide that were kindly provided by the Geological Survey Department ( Figure 5), as well as from site visits, the DInSAR displacement values significantly underestimated the actual ground movements. To investigate the reason behind the small displacement values calculated through the proposed DInSAR methodology, the coherence values were further studied in both areas.

Coherence Change Analysis
The coherence products that were developed after the completion of the "Development of Interferogram" step of the processing chain were geo-located and stacked, resulting in the development of the average coherence maps. The results for Pissouri and Petra tou Romiou are presented in Figure 14.
In both cases, within the two AOIs, coherence values are lower compared to the entire scenes. In the case of the AOI by the A6 motorway near Pissouri (Figure 14a), the mean coherence of the AOI was 0.390, with coherence values ranging from 0.286 to 0.535 with a standard deviation of ± 0.078. In the case of the AOI near Petra tou Romiou (Figure 14b), coherence values range from 0.378 to 0.771 with a mean coherence of 0.546 and a standard deviation of ± 0.086.
However, there are other areas in both scenes individually that are characterized by similar values. The low coherence values could be due to the residing vegetation surrounding the case study site near Pissouri and the vegetation towards the north and west from the AOI near Petra tou Romiou [47]. Additionally, as both events were induced by heavy rainfall, this is an important factor along with changes in soil moisture that can alter the radar signal, significantly leading to loss of coherence, as analyzed in the introduction of the present paper [33,47,60,61]. To remove the effect of vegetation from the phase decorrelation (coherence loss), NDVI was calculated from Equation (4), based on a Sentinel-2 satellite image, dated 3 February 2019. It was the closest to the landslides cloud-free and atmospherically-corrected image (Level-2A). All areas with NDVI >0.2 were removed from the final coherence maps to exclude the vegetation's impact on phase decorrelation in the two areas of interest. The removal of the areas of vegetation from the coherence To remove the effect of vegetation from the phase decorrelation (coherence loss), NDVI was calculated from Equation (4), based on a Sentinel-2 satellite image, dated 3 February 2019. It was the closest to the landslides cloud-free and atmospherically-corrected image (Level-2A). All areas with NDVI >0.2 were removed from the final coherence maps to exclude the vegetation's impact on phase decorrelation in the two areas of interest. The removal of the areas of vegetation from the coherence maps is presented in Figures 15 and 16 for the area affected by the landslide near Pissouri and Petra tou Romiou respectively. maps is presented in Figures 15 and 16 for the area affected by the landslide near Pissouri and Petra tou Romiou respectively.   Geosciences 2020, 10, x FOR PEER REVIEW 19 of 26 maps is presented in Figures 15 and 16 for the area affected by the landslide near Pissouri and Petra tou Romiou respectively.   Indeed, at the case study of the landslide by the A6 Paphos-Limassol motorway near Pissouri, there seems to be a substantial match of a scale of 85% between the area of landslide and the area of very low coherence (< 0.4), as presented in Figure 15. The areas to the north, northwest and east of the landslide area are covered with vegetation, a matter that could explain the loss of coherence between the SAR acquisitions. These were removed through the application of the NDVI mask. However, due to the fact that there is very limited vegetation at the specific slope, with the exception of the part at the top of AOI, within the determined area of landslide, other factors such as the heavy rainfall or a sudden change in the topography could be the reasons for such a substantial phase decorrelation.
At the second case study, near Petra tou Romiou (Figure 16), there are also indications that a sudden change in the topography has occurred within the predefined landslide area. However, there is some sparse vegetation that resides within the defined AOI and in close proximity to the area. The vast majority of the vegetation was successfully removed with the application of the NDVI mask, especially to the north and northwest part of the scene. Apart from the AOI, there seems to be an area (yellow circle) where there are low coherence values (< 0.4).
From the first overview of the coherence maps, there seems to be valuable information regarding the extents of the two landslide events that took place on the 15th and 20th of February 2019. Based on the findings from the analysis of the coherence values obtained from the processing of the co-event interferometric SAR pairs, and the sudden coherence loss between the 12 days of the temporal resolution of the SAR images, there are clear indications that changes to the Earth's surface occurred. DInSAR failed to detect these changes in the area's topography. In fact, there is a substantial underestimation of the slope deformation that took place at both areas under study.
Based on the initial findings from the coherence analysis, a comprehensive study should be carried out in the future to differentiate between coherence loss due to the landslides and due to the natural decrease in coherence due to the temporal baseline of the SAR images and the residing vegetation and meteorological conditions.

Discussion-Conclusions
Areas of rough topography, such as mountainous areas and areas with steep slopes, can be affected by geometric distortions, i.e., shadowing and layover, due to satellite acquisition geometries and local topography, cannot be monitored at all. Therefore, the selection of imagery of appropriate acquisition geometry and the knowledge of the residing meteorological conditions at the moment of image acquisitions are two very important factors to obtain accurate results.
A proper selection of the suitable images for the case study was performed, with minimum temporal baselines and avoiding, at the same time, geometric distortion issues that could arise from the area's topography. Moreover, a stack of interferograms was created with six individual co-event interferometric pairs in order to minimize the atmospheric effect on the displacement calculations and provide us with displacement in the vertical and east-west directions.
Based on the information obtained through photos of the landslides that were kindly provided by the Geological Survey Department, as well as from site visits, the DInSAR displacement results are significantly underestimated. In fact, as supported from the review that was carried out in the introduction, the moderate spatial resolution of the IW swath images for Level-1 SLC products (2.7 × 22 m to 3.5 × 22 m), which are available for Cyprus, has proven to be inadequate for the detailed study and monitoring of complex ground displacement phenomena, such as fast-moving landslides and slope deformations [39][40][41]52]. This is mainly due to a limitation called temporal phase aliasing, which affects the accuracy of DInSAR results due to issues created during the phase unwrapping step of the DInSAR methodology [41,50], when the displacement of a pixel or an area under investigation exceeds the value of λ/4 between two consecutive satellite acquisitions that form an interferometric SAR pair [41,51,52]. This is due to the fact that InSAR phase signals are 2π-wrapped, making it impossible to measure phase differences exceeding π unambiguously, that in turn correspond to LoS displacements over λ/4 [51,53]. Therefore, deformations that can be two, three or n times greater than the λ/4 threshold can be underestimated as they produce similar observed phase signals [41,51,53].
This was supported by the findings of our study and more specifically by the fact that displacement values were generally greater outside the extent of the area that was affected by the landslide than those within the AOI. As both landslides resulted in fast slope deformations, the phase signal that was produced was similar to that which would have been produced by an area with minimum displacement.
Moreover, due to the fact that both landslide events were triggered by heavy rainfall, this could have introduced additional errors caused by the adverse meteorological conditions and the phase decorrelation that can be caused by changes in soil moisture and rainfall. Additionally, the slopes that were affected by the landslides were sparsely vegetated, a matter that can cause even further phase decorrelation [47]. As stated earlier, one of the fundamental conditions for DInSAR methodology to be effective is that the coherence at a certain amount should remain unchanged between two subsequent SAR images [33].
To investigate the reason behind the small displacement values calculated through the proposed DInSAR methodology, the coherence values were studied in both areas of interest. The coherence maps from the six co-event interferometric pairs were stacked together to minimize the impact of the atmosphere effect, and a NDVI mask was also applied to remove the effect of vegetation. In the future, further investigation of the NDVI mask in situations with different vegetation is planned. The results from coherence analysis were promising as they showed the extents of the areas affected by landslides in both cases. However, no quantitative information can be obtained by the study of coherence variations.
The complementarity of DInSAR with coherence change analysis needs to be further investigated. A comparison of pre-event and co-event coherence maps would need to be conducted to study the coherence loss and differentiate between coherence loss due to the landslides and due to the temporal baseline of SAR images. Coherent Change Detection has been widely used for monitoring land use/land cover changes [62][63][64]; agricultural studies [65][66][67]; forestry applications [68][69][70]; and lately earthquake impact assessment [71][72][73][74]. In fact, the Coherent Change Detection methodology has been applied successfully for the detection of the same landslides, determining accurately the extents of the areas affected by the landslides [75]. Additionally, the combination of Sentinel-1A and Sentinel-1B is planned to be studied further, as the temporal baseline of 6 days that the two-satellite constellation offers might provide additional benefits in the detection of landslides.
Despite all the advantages that EO techniques have, certain limitations exist mostly because of the spatial resolution of Sentinel-1 images and the C-band wavelength. Areas covered with vegetation cannot be monitored very efficiently with Sentinel-1 images; thus, L-band and P-band SAR images, which have enhanced capabilities for penetration, should be preferred instead.
It should also be stressed, at this point, that EO-based techniques can only act as an aid to in-situ engineering investigations, proving complimentary information to that obtained from field campaigns and other observation techniques.
The use of open-source software and freely distributed images was also highlighted in the present study. Freely distributed Earth Observation datasets supported by national and international space agencies, such as those of the Copernicus space program, will continue to provide optical and radar images around the world. The need to provide tangible end-products and services to non-experts in remote sensing, which have limited access to accurate spatial datasets, thus addressing societal challenges both on a local and a regional level, is essential to the scientific community.
In the future, new ways of how Copernicus sensors can be exploited further for non-invasive monitoring of critical infrastructure resilience will be investigated. The exploitation of the Copernicus program in its full capacity, including all other sensors, is expected to enhance research in new directions and methodologies. The integration of different SAR interferometry techniques can overcome any obstacles by combining their respective advantages.