Displacement Monitoring in Airport Runways by Persistent Scatterers SAR Interferometry

Deformations monitoring in airport runways and the surrounding areas is crucial, especially in cases of low-bearing capacity subgrades, such as the clayey subgrade soils. An effective monitoring of the infrastructure asset allows to secure the highest necessary standards in terms of the operational and safety requirements. Amongst the emerging remote sensing techniques for transport infrastructures monitoring, the Persistent Scatterers Interferometry (PSI) technique has proven effective for the evaluation of the ground deformations. However, its use for certain demanding applications, such as the assessment of millimetric differential deformations in airport runways, is still considered as an open issue for future developments. In this study, a time-series analysis of COSMO–SkyMed satellite images acquired from January 2015 to April 2019 is carried out by employing the PSI technique. The aim is to retrieve the mean deformation velocity and time series of the surface deformations occurring in airport runways. The technique is applied to Runway 3 at the “Leonardo da Vinci” International Airport in Rome, Italy. The proposed PSI technique is then validated by way of comparison with the deformation outcomes obtained on the runway by traditional topographic levelling over the same time span. The results of this study clearly demonstrate the efficiency and the accuracy of the applied PSI technique for the assessment of deformations in airport runways.


Introduction
Runways are central elements in airport infrastructures, as they are mostly dedicated to the two fundamental and most critical manoeuvres of taking-off and landing of the aircrafts. In light of this, runways must comply with very strict requirements in terms of the construction [1] and the maintenance standards [2]. In this regard, it is worthy of mention that a continuously regular surface must be ensured in runways over their entire life cycle, and no superficial (e.g., cracking or rutting) or deep (e.g., subgrade subsidence) damage is acceptable to compromise safe manoeuvring of the aircrafts.
Although the provision of proper design methods for runways is crucial to minimise future maintenance and rehabilitation of the infrastructure, this can be compromised by other critical factors, such as construction site related issues (e.g., geotechnical instability of the subgrade) or the non-linear action of different heavy loads over the infrastructure life cycle. The occurrence of these events requires

Aim and Objectives
The main aim of the investigation reported in this paper is to verify the effectiveness of the satellite remote sensing technology in gaining vital information about the functionality of airport runways for inclusion in airport asset management systems. To achieve this aim, the following objectives are identified: • to assess runway displacements at the millimetre scale and evaluate their trend on a multi-year scale using the PSI monitoring technique; • to compare the results obtained with the PSI technique and the traditional topographic levelling method.

Runway Monitoring Techniques
Relevant methods employed for the monitoring of airport runways are presented in the following sections. The topographic levelling and the LiDAR techniques are here referred to as "established techniques", as opposed to the PSI method, referred to as an "innovative technique".

Topographic Levelling
The process of measuring variations in elevation is a relatively basic operation in topographical surveys, and it is typically referred to as levelling. Various levelling techniques have been developed over the time. In regard to infrastructure surveying, geometric levelling is the most adopted [14]. The height difference is here obtained from readings on levelling staffs where the level's horizontal sightline intersects them. The level is mostly used in the middle between two levelling staffs (differential levelling), with the objective of the survey being the definition of the difference in vertical distance between the two staff positions.
As demonstrated by [6], the topographic levelling procedure allows us to monitor the long-term settlements of an infrastructure by measuring the variation in the vertical position of the target over multiple surveys. This is performed with respect to a single or multiple stable point, typically referred to as topographic benchmarks.
In recent decades, various levelling investigations on civil infrastructures have been reported in the literature [15,16]. In regard to the monitoring of the airport runways, both advantages and drawbacks can be mentioned on the use of the topographic levelling. As for the advantages, the high accuracy of the measurements [17] and the possibility to perform tests independently from indoor or outdoor environments once a reference point is defined [18] is worthy of mention. The drawbacks include a) a limited productivity [19,20], due to the need of measuring each target separately, b) the necessity to close the runway during testing, c) the impossibility to perform the survey in adverse weather conditions [15], and d) the provision of a clear line of sight between consecutive targets of the survey, implying several practical constraints. Lastly, the observed displacements are relative measurements that might be affected by potential settlements levelling benchmarks, which may result in a distortion of the results.

LiDAR Surveys
LiDAR is a surveying method that measures the distance to a target by illuminating the target with a laser light and measuring the reflected light with a sensor [21]. The differences in the laser return, times, and wavelengths can then be used to create digital 3D representations of the target. Since this technology has started to spread in several scientific and professional fields, successful applications in transport infrastructures monitoring have been reported [22] for both the static [23] and the mobile [24] configurations of the equipment.
With regard to airport runway monitoring, static terrestrial laser scanners are being mostly used for reconstructing the geometry of the pavement surface in order to detect defects and decayed areas [8,25]. This is carried out by illuminating the whole runway surface by means of different scans performed at different and distributed survey stations. A certain rate of superposition of the point clouds must also be provided for matching purposes [26]. The georeferencing operation between different scans is ensured by the use of common ground targets working as reference objects. However, by coupling a GNSS receiver to this ground targets, it is possible to skip from relative to global coordinates and compare successive surveys, in order to monitor the evolution of differential settlements.
Among the advantages of this technique, it is worth mentioning the rapidity and the high resolution of the surveys, as millions of point clouds can be collected in limited time. Furthermore, the use of global georeferenced ground targets allows us to prevent the measurements from relative errors due to the instability of the reference. In turn, LiDAR surveys require the runaway to be closed to traffic during the tests. In addition, the application of this technique, as well as the levelling, requires the presence of operators on site. This represents a considerable safety concern as airport runways and aprons are considered as high-risk environments.

Persistent Scatterers Interferometry (PSI)
The InSAR technique, or SAR interferometry, relies on the measurement of the signal phase variation between images acquired by a satellite orbiting over the same area [9]. Indeed, once the phase contribution related to atmospheric conditions and both temporal and spatial decorrelations are adequately accounted for, it is possible to detect a time sequence of the sensor-target distance along the sight direction of the satellite. This can be related to the surface deformations, e.g., the subsidence. [11].
For this purpose, various processing techniques have been proposed over time. Among these, PSI is one of the most acknowledged [12,13]. PSI is based on the statistical analysis of the signals back-scattered from a network of phase-coherent targets, namely the Persistent Scatterers (PS), which are defined as the points on the ground returning stable signals to the satellite sensor. Indeed, the constant scattering properties of the PS over time and the reflection dominance within a pixel cell are effective in reducing the temporal and geometric decorrelations. In addition, the atmospheric contribute can be estimated and removed using the series of images acquired at different times.
SAR sensors operate at different bands of the microwave domain, namely X, C, and L, corresponding to wavelengths (λ) ranging from 2.4 to 30 cm. Each wavelength is associated to a differently sized resolution cell on the ground, whose displacement trend is described by a single PS.
By means of the PSI technique, satellite remote sensing can provide a continuous monitoring of the overall stability of structures and infrastructures, as well as the surrounding environment. This is obtained by the analysis of multiple SAR images collected at different time stages. Subject-related literature reports the SAR techniques as viable tools for the monitoring of ground deformations, landslides, subsidence, and tectonic motions [27][28][29]. Similarly, in the last few years, various successful applications of the PSI technique have been presented, proving the feasibility of this technique for the assessment of transport infrastructures and surveillance areas, such as highways [30][31][32], bridges [33][34][35][36][37], subways and tunnels [31,38,39], railways [40][41][42][43], and airport runways [44][45][46]. This evidence confirms the wide applicability of these techniques within these specific areas of endeavour.
The use of the InSAR technique in transport infrastructure monitoring holds several advantages [9]. InSAR data can be collected regardless of the atmospheric and lighting conditions. In addition, a single InSAR survey permits the analysis of extended areas, due to the wide footprint of the sensor. The continuous motion of the satellites ensures the availability of regularly spaced images, which allows us to perform very dense analyses as opposed to on-site and low-frequency inspections. The acquisition and processing of SAR images do not require on-site operations, thereby preventing both the closure of the runway to air traffic and the presence of operators on the site, with related economic and safety benefits.
On the other hand, the applicability of the PSI analysis is by definition limited to areas where an adequate number of PSs can be observed. This implies that, in case of varying surface conditions (e.g., frequent repaving or accelerated degradation), the stability of the scattering properties of the target may be compromised, with a reduced number of PSs being detected. Furthermore, reliable InSAR assessments require the processing of various SAR images in order to detect statistically stable PSs. This occurrence involves potential computational-related issues, due to the size of the database required at each survey. Lastly, according to the frequency of the adopted sensor, uncertainties may arise in the scattering source, as it is impossible to recognise the scattering object within the resolution cell.

Site Description
In the present study, a case study is presented where the PSI technique is applied Runway 3 at the "Leonardo da Vinci" International Airport in Rome, Italy.
The airport is located in the area of Fiumicino, about 30 km on the west of Rome, and carries most of the intercontinental air traffic from and towards Italy, that ranks it as one of the major airports in Europe. With an excess of 43 million passengers in 2018 and over 199,000 tons of traffic, the airport is the largest in Italy by number of passengers/year and ranks second by number of cargo flights/year. The airport has three terminals reserved for domestic, international, and intercontinental flights and three runways. Initially, the layout of the airport included two runways (Runway 1 and Runway 2) only. However, due to an increasing traffic demand, it was expanded with a new runway located in the north-east area of the airport, which develops in the north-south direction ( Figure 1) The area where Runway 3 was realized is known to be affected by the presence of clayey and peaty soils. As a result, the runway was affected by long-term differential subsidence effects, which required accurate monitoring and had to be brought to a full-depth pavement rehabilitation in 2015. This was aimed at increasing the bearing capacity of the subgrade in the norther section of the runway, where severe differential settlements had been recorded in the previous years. Nevertheless, Runway 3 remains a critical asset that is levelled every year in order to be able to monitor the deformation trends and plan maintenance activities.

Levelling Data
In regard to the on-site topographic surveys conducted on the survey area, a classical geometric levelling survey of the runway was conducted by means of the DNA03 Digital Level system, manufactured by Leica. The elevation measurements were collected by means of a 2 m high levelling rod made of invar, i.e., a metal with a limited thermal expansion coefficient. The main features of the employed levelling system are summarised in Table 1. In particular, in situ levelling data were collected by multiple closed-loop levelling nets covering the entire Runway 3. The collected nets were automatically compensated, with an average squared root mean error of 0.94 mm. The starting and ending point of the measure was a levelling benchmark located in the north-west corner of the runway area, which was verified to be stable in elevation. This point was connected to the high precision levelling net developed by the Istituto Geografico Militare (IGM), through a levelling line having an average accuracy of 1.0 mm/km.
Tests were performed every year from 2015 to 2019 and covered five sections along the runway with a transversal spacing of 15 metres and a length equal to the entire runway longitudinal development ( Figure 2). Furthermore, the output of these data allows us to detect and exactly quantify the displacements and the average velocity in the investigated time period. Therefore, this information is significant for a validation of the displacements detected by the InSAR technique. The area where Runway 3 was realized is known to be affected by the presence of clayey and peaty soils. As a result, the runway was affected by long-term differential subsidence effects, which required accurate monitoring and had to be brought to a full-depth pavement rehabilitation in 2015. This was aimed at increasing the bearing capacity of the subgrade in the norther section of the runway, where severe differential settlements had been recorded in the previous years. Nevertheless, Runway 3 remains a critical asset that is levelled every year in order to be able to monitor the deformation trends and plan maintenance activities.

Levelling Data
In regard to the on-site topographic surveys conducted on the survey area, a classical geometric levelling survey of the runway was conducted by means of the DNA03 Digital Level system, manufactured by Leica. The elevation measurements were collected by means of a 2 m high levelling rod made of invar, i.e., a metal with a limited thermal expansion coefficient. The main features of the employed levelling system are summarised in Table 1. In particular, in situ levelling data were collected by multiple closed-loop levelling nets covering the entire Runway 3. The collected nets were automatically compensated, with an average squared root mean error of 0.94 mm. The starting and ending point of the measure was a levelling benchmark located in the north-west corner of the runway area, which was verified to be stable in elevation. This

SAR Imagery
In regard to the application of the InSAR technique, a multi-temporal interferogram analysis of SAR images, namely the Persistent Scatterers Interferometry technique (PSI), was applied. To this effect, two different data stacks were acquired in ascending and descending geometries using highresolution SAR imagery acquired in X-Band, which allows the detection of displacements with a millimetre accuracy.
In more detail, a dataset of 72 Stripmap images collected in ascending and descending geometries from the COSMO-SkyMed mission (COSMO-SkyMed Product-©ASI: Italian Space Agency, 2015-2019, All Rights Reserved) were processed. The system operates in the X-band corresponding to a wavelength of 3.1 cm, with a 3 m ground-resolution cell. The radar antenna is a phased array that is 1.4 m wide x 5.7 m long. The system is capable of both single-and dualpolarisation data collection. The central frequency is 9.6 GHz with a maximum radar bandwidth of 400 MHz. The main features of the SAR dataset are reported in Table 2.

Data Processing
These products have been acquired and processed using the PS technique of SARscape Interferometric Stacking Module [47], integrated in the Envi software, within the framework of the project "MOBI: MOnitoring Bridges and Infrastructures networks" (proposal ID 46829), approved by the European Space Agency (ESA).
The processing algorithm includes the following steps [12,13,48,49]:  Furthermore, the output of these data allows us to detect and exactly quantify the displacements and the average velocity in the investigated time period. Therefore, this information is significant for a validation of the displacements detected by the InSAR technique.

SAR Imagery
In regard to the application of the InSAR technique, a multi-temporal interferogram analysis of SAR images, namely the Persistent Scatterers Interferometry technique (PSI), was applied. To this effect, two different data stacks were acquired in ascending and descending geometries using high-resolution SAR imagery acquired in X-Band, which allows the detection of displacements with a millimetre accuracy.
In more detail, a dataset of 72 Stripmap images collected in ascending and descending geometries from the COSMO-SkyMed mission (COSMO-SkyMed Product-©ASI: Italian Space Agency, 2015-2019, All Rights Reserved) were processed. The system operates in the X-band corresponding to a wavelength of 3.1 cm, with a 3 m ground-resolution cell. The radar antenna is a phased array that is 1.4 m wide x 5.7 m long. The system is capable of both single-and dual-polarisation data collection. The central frequency is 9.6 GHz with a maximum radar bandwidth of 400 MHz. The main features of the SAR dataset are reported in Table 2.

Data Processing
These products have been acquired and processed using the PS technique of SARscape Interferometric Stacking Module [47], integrated in the Envi software, within the framework of the project "MOBI: MOnitoring Bridges and Infrastructures networks" (proposal ID 46829), approved by the European Space Agency (ESA).
The processing algorithm includes the following steps [12,13,48,49]: • Generation of differential interferograms out of the stack of SAR images; • Implementation of High definition Digital Elevation Models (DEM) for topographic phase-term removal; In regard to the implementation of the high-definition DEM models, an SRTM v3 "Shuttle Radar Topography Mission" DEM was collected and implemented in the interferometric process [50,51] in order to identify and subtract phase-related parameters linked to the topography of the investigated area. The DEM, with a pixel resolution of 3 arc-second (90 m × 90 m), is made available by NASA in partnership with the United States Geological Survey (USGS) [50].
The outputs of the PSI processing algorithm were exported into a GIS environment and the PSs were displayed as a function of the average annual-motion velocities.
A specific procedure was adopted in order to compare each ground-truth levelled data with a single satellite-derived displacement measure, as follows: In regard to the implementation of the high-definition DEM models, an SRTM v3 "Shuttle Radar Topography Mission" DEM was collected and implemented in the interferometric process [50][51] in order to identify and subtract phase-related parameters linked to the topography of the investigated area. The DEM, with a pixel resolution of 3 arc-second (90 m × 90 m), is made available by NASA in partnership with the United States Geological Survey (USGS) [50].
The outputs of the PSI processing algorithm were exported into a GIS environment and the PSs were displayed as a function of the average annual-motion velocities. A specific procedure was adopted in order to compare each ground-truth levelled data with a single satellite-derived displacement measure, as follows:  Finally, in order to obtain continuous average-velocity maps representative of the runway condition, a geo-statistical gridding method was implemented both to the PS datasets and the levelling measurements.
The geo-statistical ordinary Kriging method was used for this purpose [52,53]. Kriging is applied due to the flexibility and the high accuracy in the gridding methods and the provision of representative maps applied to different types of datasets. Moreover, it can compensate for clustered data by weighting less than the cluster in the overall prediction. Each grid node value is based on the known data points neighbouring the node. Each data point is weighted by its distance away from the Finally, in order to obtain continuous average-velocity maps representative of the runway condition, a geo-statistical gridding method was implemented both to the PS datasets and the levelling measurements.
The geo-statistical ordinary Kriging method was used for this purpose [52,53]. Kriging is applied due to the flexibility and the high accuracy in the gridding methods and the provision of representative maps applied to different types of datasets. Moreover, it can compensate for clustered data by weighting less than the cluster in the overall prediction. Each grid node value is based on the known data points neighbouring the node. Each data point is weighted by its distance away from the node, and consequently, points that are further from the node will be weighted less in the estimation of the node. To compute theẐ(x 0 ) value at a randomly given grid node at position x 0 , the following equation is used [52,53]:Ẑ whereẐ(x 0 ) is the estimated value of the grid node, n is the number of neighbouring data values used in the estimation, Z(x i ) is the observed value at the ith location weighting w i (x 0 ) with i ranging between 1 and n. The values of the weights add up to 1 in order to ensure that no bias occurs towards clustered data points. The weights are intended to summarise two important procedures in a spatial inference process, i.e., (i) to reflect the structural proximity of the samples to the estimation location; (ii) in parallel, they should not have a separation effect in order to avoid bias effects caused by potential sample clusters.

Results
As a qualitative assessment of the reliability of the PSI technique, Figure 4 shows the velocity maps obtained from the levelling and the PSI techniques by interpolation of the displacement velocity data. The similarity between the two maps in low (Figure 4a where Z ( 0) is the estimated value of the grid node, n is the number of neighbouring data values used in the estimation, Z(xi) is the observed value at the ith location weighting wi(x0) with i ranging between 1 and n. The values of the weights add up to 1 in order to ensure that no bias occurs towards clustered data points. The weights are intended to summarise two important procedures in a spatial inference process, i.e., i) to reflect the structural proximity of the samples to the estimation location; ii) in parallel, they should not have a separation effect in order to avoid bias effects caused by potential sample clusters.

Results
As a qualitative assessment of the reliability of the PSI technique, Figure 4 shows the velocity maps obtained from the levelling and the PSI techniques by interpolation of the displacement velocity data. The similarity between the two maps in low (Figure 4a   The correlation and the loss performance of the prediction are reported in Table 3. This includes Pearson's coefficients (r) and root squared mean error (RSME) values for the five survey profiles. It is worthy to note that the survey profile L3 returns slightly less accurate results, as shown in both Figure 5 and Table 3. Especially for the highest values of velocity (which are mainly included within L3), the PSI method seems in fact to underestimate the ground-truth levelled displacement trend. The nature of this coherent error is most likely related to the atmospheric noise contribution, which seems not to be completely filtered out through the applied processing [53,54]. However, the extent of the error is found to be quite limited (around 0.5 cm/year, approximatively). Accordingly, it is observed that the potential of the method in reconstructing the deformation behaviour of the runway is not significantly affected. This allows the provision of very useful information to airport managers for scheduling flight-lists and planning strategic on-site monitoring operations.
As a further confirmation of the above observations, Figure 6 shows the displacement velocity of the survey profiles against the space (WGS84 N Coordinate). Still, it is possible to note that as the displacement trend reaches values higher than 15 mm/year (Figure 6 b,c, and e), the PSI method turns out to slightly underestimate the deformation rate. However, besides this specific observation, Figure  6 succeeds in demonstrating the effectiveness of PSI in reconstructing the actual deformation pattern over the inspected infrastructure. The correlation and the loss performance of the prediction are reported in Table 3. This includes Pearson's coefficients (r) and root squared mean error (RSME) values for the five survey profiles. It is worthy to note that the survey profile L3 returns slightly less accurate results, as shown in both Figure 5 and Table 3. Especially for the highest values of velocity (which are mainly included within L3), the PSI method seems in fact to underestimate the ground-truth levelled displacement trend. The nature of this coherent error is most likely related to the atmospheric noise contribution, which seems not to be completely filtered out through the applied processing [53,54]. However, the extent of the error is found to be quite limited (around 0.5 cm/year, approximatively). Accordingly, it is observed that the potential of the method in reconstructing the deformation behaviour of the runway is not significantly affected. This allows the provision of very useful information to airport managers for scheduling flight-lists and planning strategic on-site monitoring operations.
As a further confirmation of the above observations, Figure 6 shows the displacement velocity of the survey profiles against the space (WGS84 N Coordinate). Still, it is possible to note that as the displacement trend reaches values higher than 15 mm/year (Figure 6b,c,e), the PSI method turns out to slightly underestimate the deformation rate. However, besides this specific observation, Figure 6 succeeds in demonstrating the effectiveness of PSI in reconstructing the actual deformation pattern over the inspected infrastructure.  A quantitative analysis of the comparison between the two methods is instead reported in Table  4. Specifically, for every levelled point, the displacement rate is reported for both the survey methodology, with the reliability of the InSAR measurement being expressed by the root squared mean error.  A quantitative analysis of the comparison between the two methods is instead reported in Table  4. Specifically, for every levelled point, the displacement rate is reported for both the survey methodology, with the reliability of the InSAR measurement being expressed by the root squared mean error. A quantitative analysis of the comparison between the two methods is instead reported in Table 4. Specifically, for every levelled point, the displacement rate is reported for both the survey methodology, with the reliability of the InSAR measurement being expressed by the root squared mean error.

Conclusions
This work demonstrates the applicability of the high-resolution X-band COSMO-SkyMed mission data and the Persistent Scatters Interferometry (PSI) technique in monitoring the deformations occurring on airport runways.
To this purpose, a time-series analysis of a set of satellite images acquired from January 2015 to April 2019 is carried out by employing the PSI technique. To assess the actual reliability of the method in reconstructing the exact deformation pattern of the monitored infrastructure, the outcomes from the PSI technique were compared to those obtained by the traditional topographic levelling technique, which has been applied on the runway at the same time range as the satellite analysis.
The outcome of the analysis clearly indicates that the use of the PSI technique is reliable and accurate for deformation assessment purposes.
More specifically, the application of the PSI technique has proven to be effective in reconstructing the average trend of the annual deformation velocity (see the coloured maps in Figure 4 and the profiles reported in Figure 6) and the monthly deformation time series of each levelled point (see the examples included in Figure 7).
As opposed to the advantages related to the application of this method, it is fair to comment that in the case of the highest rates of deformation observed in the investigated runway, the PSI technique was found to slightly underestimate the subsidence velocity. Such a coherent bias was related to a limited processing of the atmospheric noise contribution, and as a consequence, it is expected to be accounted for by means of an advanced processing phase.
In general, the results presented in this research demonstrate that the PSI technique is worthy of implementation in Airport Pavement Management Systems (APMS), which may profit by a significant increase in the efficiency in the scheduling of maintenance operations.
Author Contributions: Conceptualization, L.B.C. and A.C.; methodology, A.C. and F.D.; investigation, V.G. and C.F.; data curation, L.B.C. and F.D.; writing-original draft preparation, L.B.C. and F.T.; writing-review and editing, F.D., A.C., and F.T. All authors have read and agreed to the published version of the manuscript.