ADAtools: Automatic Detection and Classification of Active Deformation Areas from PSI Displacement Maps

Abstract

This work describes the set of tools developed, tested, and put into production in the context of the H2020 project Multi-scale Observation and Monitoring of Railway Infrastructure Threats (MOMIT). This project, which ended in 2019, aimed to show how the use of various remote sensing techniques could help to improve the monitoring of railway infrastructures, such as tracks or bridges, and thus, consequently, improve the detection of ground instabilities and facilitate their management. Several lines of work were opened by MOMIT, but the authors of this work concentrated their efforts in the design of tools to help the detection and identification of ground movements using synthetic aperture radar interferometry (InSAR) data. The main output of this activity was a set of tools able to detect the areas labelled active deformation areas (ADA), with the highest deformation rates and to connect them to a geological or anthropogenic process. ADAtools is the name given to the aforementioned set of tools. The description of these tools includes the definition of their targets, inputs, and outputs, as well as details on how the correctness of the applications was checked and on the benchmarks showing their performance. The ADAtools include the following applications: ADAfinder, los2hv, ADAclassifier, and THEXfinder. The toolset is targeted at the analysis and interpretation of InSAR results. Ancillary information supports the semi-automatic interpretation and classification process. Two real use-cases illustrating this statement are included at the end of this paper to show the kind of results that may be obtained with the ADAtools.

1. Introduction

The Multi-scale Observation and Monitoring of Railway Infrastructure Threats (MOMIT) project (see [1] for details) aims to develop and demonstrate a new use of remote sensing techniques for railway infrastructure monitoring. MOMIT solutions are targeted at supporting the maintenance and prevention processes within the infrastructure management lifecycle. The overall concept underpinning MOMIT is the demonstration of the benefits brought by Earth Observation and Remote Sensing data to the monitoring of railways networks both in terms of the infrastructure and the surrounding environment, where activities and phenomena impacting the infrastructure could be present. MOMIT leverages state-of-the-art techniques in the fields of space-based remote sensing and remotely-piloted aircraft systems (RPASs) to perform different kinds of analysis thanks to the wide variety of sensors they may be equipped with.

To achieve its goals, six demonstrators showing how these data and techniques may contribute to such objectives have being built, namely:

• Ground movements nearby the infrastructure.
• Hydraulic activities nearby the tracks.
• Global supervision for natural hazards.
• Electrical system monitoring.
• Civil engineering structures monitoring.
• Safety monitoring.

The Division of Geomatics of the Centre Tecnològic de Telecomunicacions de Catalunya (CTTC) took care of building some of the components integrating the first demonstrator, whose objectives are detailed in [1]. The University of Alicante (UA) mainly defined the methodology used in the first demonstrator, supported CTTC in the tuning of the software, and applied the tools developed in several case studies. The goal most relevant to the work presented here is to introduce the ADAtools, a set of software components targeted at the detection and interpretation of active deformation areas (ADA) using displacement maps created by means of persistent scatterer interferometry (PSI) method.

Several PSI approaches have been developed in the last twenty years; a review of them is provided in [2]. The basic concept behind the PSI techniques is to collect a stack of several synthetic aperture radar (SAR) images for the same area, with very similar acquisition angles; these images are later processed to detect the so-called persistent scatterers (PS), which are artificial and natural structures that show interferometric coherence over the time; through advanced data processing, it is possible to determine the displacement time series (TS) of the detected PS. PSI results have been applied for several applications at different working scales, ranging from wide-area processing at national level [3,4] to single infrastructures monitoring [5,6,7]. PSI is widely accepted as a reliable tool for the precise measure of a variety of geohazards, including landslides [8,9,10], natural and anthropogenic subsidence [11,12,13], sinkholes [14,15,16], earthquakes [17,18,19], and volcanoes [20,21,22]. Nowadays, the increasing number of PSI applications and the tendency to process wide areas with millions of measurement points require the definition of reliable semi-automatic and automatic tools to ease the analysis and interpretation of the PSI results [23,24].

This paper describes in depth a set of tools for PSI data analysis and interpretation developed at CTTC, the so called ADAtools. The package is composed of 4 modules, namely, ADAfinder, ADAclassifier, THEXfinder, and los2hv [25]. The main goal of such applications is a semi-automatic extraction and preliminary interpretation of the areas affected by deformation detected by the PSI technique. The goal is to update and assess the activity of phenomena related to geohazards (volcanic activity, landslides, or ground subsidence, among other phenomena) or human activities of a given area. The first one, ADAfinder, is a tool dedicated to the detection of the ADA, extracted from a PSI-derived displacement map. ADAclassifier and THEXfinder go a step beyond ADAfinder, trying to classify the kind of deformation process undergone by the ADA, and trying to answer the question: is the detected deformation process a landslide, a sinkhole, or something else? Finally, los2hv computes the East-West horizontal and vertical components of the movement measured along the satellite’s line of sight (LOS). The East-West horizontal component of motion is also one of the inputs for the ADAclassifier.

Finally, but no less important, it is worth noting that although the development of these tools was motivated by the needs of the MOMIT project, the applications described in the next sections are not limited to this purpose. On the contrary, these may be used to track deformation processes occurring anywhere, as is shown by the real use-cases proposed in Section 7.

2. ADAtools

The four tools presented in this paper rely in preexisting methodologies, i.e., these applications have automated a set of procedures that already existed and were executed manually, step by step, using Geographic Information System (GIS) tools.

In the case of ADAfinder, the methodology to perform the identification and assessment of ADA was presented in [26,27] back in 2017 and 2019, respectively. In these works, the authors explain in detail the procedure to identify the active deformation areas and to assess their quality. The input is a set of persistent scatterers(PS; or “the points”) covering the area to analyze. Each ADA is then defined on the basis of the location and density of PS, depending on some thresholds as the minimum number of points making an ADA or the area of influence of each PS. A quality index describing the noise level and the consistency of the displacement time series of the PS forming each ADA is calculated as well. The concept of ADA has been exploited to estimate geohazard-related risk (mainly landslides) in different European environments (see [28,29,30]). It is worth underlying that the ADA extraction does not overcome the intrinsic limitations of the PSI technique: ADAfinder extracts what PSI technique can detect. As an example, the absence of ADA does not necessarily mean no-deformation; it could be absence of deformation as well as a not-detectable movement due to an unfavorable geometry [31] or absence of information for low coherence (noisy area).

ADAclassifier and THEXfinder rely on a methodology developed by the UA that is under constant research; although it is mature enough to be automated, which eases the experimentation process and helps to improve the methodology itself. The task of ADAclassifier and THEXfinder is to identify the kind of geological or anthropogenic process motivating the presence of ADA. Up to six different phenomena (see Section 2.2) are probed to obtain an estimation of the causes motivating the ground deformation; these are landslides, sinkholes, subsidence, constructive settlements, expansive soils, and thermal phenomena. ADAclassifier takes care of the first four, while THEXfinder is responsible for the last two. For more details about the methodology on which these two tools rely, see [32].

Finally, los2hv computes the East-West horizontal and vertical components of the ground displacement measured with PSI along the satellite’s LOS (the horizontal component is one of the inputs needed by ADAclassifier). los2hv performs a tessellation of the area of interest to obtain an averaged result for each of the resulting tiles. Such an approach means that overcoming the strict requirement stating that data from ascending and descending orbits are necessary to perform such decomposition of the movement; treating the PS in the same tile as a single coarser point makes possible to obtain an averaged but still useful value for the purposes of ADAclassifier. See Section 2.3 and [32] for a detailed description of this tool. As stated above, the aforesaid procedures usually rely on the heavy use of GIS tools and the expertise of the operator; the number of steps these procedures consist of makes them time-consuming, error-prone processes, which entails the need for qualified personnel.

The goals of the ADAtools package are: (1) to automate the respective procedures and to limit unnecessary human errors; (2) to reduce the time needed to identify and pre-classify the ADA, thus opening the door to more frequent updates and analyses; and (3) to reduce the expertise required to obtain such results, making possible to integrate the process in a semi-automated production workflow, if necessary.

2.1. ADAfinder

This application, based on the methodology explained in [26], is used to identify the main areas where a displacement has been measured by the PSI processing, squeezing the information contained in the input deformation map and assessing the quality of the time series information (i.e., spatial-temporal noise) of each ADA. Additionally, ADAfinder includes the option to filter the input displacement map from the isolated or potential outlier PS. We refer to [26] for more details on the quality assessment and the outlier filtering methods. This tool allows the amount of data to be managed in terms of both hardware space and computational time to be reduced.

As inputs, ADAfinder requires:

• An ESRI shapefile containing the PS that will be used by the detection algorithm. Besides their coordinates, ADAfinder needs some attributes defining the PS; their average velocity expressed in mm/yr, and the deformation time series measuring the movements undergone by these.
• Optionally, the user can upload a polygon to resize the area of interest. All PS in the input shapefile are considered when such a polygon is not provided.
• A set of parameters defined by the methodology described in [26] (Figure 1).

The outputs are two ESRI shapefiles, the first containing the polygons defining the boundaries of the ADA; the second including the set of filtered PS. For what concerns the PS output, ADAfinder allows choosing between two different options, depending on the user necessities. The output can include all the filtered PS of the area of interest or only the ones included in the extracted ADA.

Figure 1 depicts the graphical user interface (GUI) version of the ADAfinder tool. Note that some of the values to be provided by the user are already set; the reason is that the tool loads a default options file (see Section 3.3 for details on options files) that serves two purposes: providing some hints to the user about the values of some parameters and saving time spent in typing many of them when working in the same project. The default options file may be changed by the user, to adjust the defaults to their preferences at any time. Finally, the user may also save a set of options differing from those included in the default in a separate options file that later may be loaded as many times as needed.

It is worth remarking that the output ADA shapefile includes among its attributes the results of the quality assessment of the PS displacement TS in terms of both spatial and temporal noise, based on simple statistical analysis of the TS within each ADA [3]. Such assessment is represented by a four-level classification of TS reliability of the extracted ADA, where 1 means “very reliable”, 2 means “reliable” but an analysis of the TS is suggested, 3 means “not so reliable”, i.e., a deeper analysis of the TS is necessary, and 4 means “not reliable”. Figure 2 shows some ADA; their colors map the four levels just mentioned.

For a detailed description of the ADAfinder tool, please refer to its user guide [33].

2.2. ADAclassifier and THEXfinder

These tools try to determine the kind(s) of deformation process(es) undergone by the terrain.

Up to six different kinds of deformation phenomena are detected. The ones analyzed by ADAclassifier are landslides, sinkholes, land subsidence, and constructive or consolidation settlements. THEXfinder takes care of the identification of deformation processes due to expansive soils and temperature effects. A different algorithm (or sub-methodology) has been devised for each of these. Figure 3 depicts the workflow of the algorithms implemented by ADAclassifier to detect the four aforementioned phenomena.

ADAclassifier and THEXfinder need a substantial number of inputs. For ADAclassifier these are:

• The ADA and PS files created by ADAfinder (see Section 2.1). Strictly speaking, these files do not need to be created by ADAfinder; any other tool or manual process identifying ADA may be used instead. However, the set of attributes included in the attribute table of the shapefiles must match those required by ADAclassifier –attributes that ADAfinder does include in its output.
• A digital terrain model (DTM), to compute slopes.
• A series of polygon vector maps (inventories from now on,) in the form of ESRI shapefiles, to check whether an ADA has already been catalogued as belonging to any of the four aforesaid deformation processes. The required inventories are those for landslides, sinkholes, land subsidence, and infrastructures. A geologic map (another polygon vector map) is also needed. In this last case, a read-map file defining how the inventory is structured is also needed to point to the attributes stating the kind of soil covered by each polygon in the inventory. See Section 3.4 for a detailed description about the so-called read-map files.
• An ESRI (polygon) shapefile storing the horizontal component of the movement for the study area. This is the output of los2hv (see Section 2.3 for details).
• The set of parameters—typically thresholds—needed by the different algorithms in charge of the classification processes must be supplied. Examples of such parameters are slopes, determination coefficients to state whether some statistical check is positive, or the minimum percentage of overlap of an ADA and the polygons in some inventory to consider that they do intersect. These thresholds appear as Th1 to Th11 in Figure 3.

The output of ADAclassifier is another file with ADA, where the attribute table has been extended to include four additional fields. Each of them states the probability that the ADA belongs to the corresponding deformation process. This is so because all of the detection algorithms are applied to each ADA. Consequently, and although it could be considered incongruous, some ADA might be classified as positives in more than one deformation process. This is due to the fact that the ADAclassifier is making a most probable cause analysis.

Four levels of certainty are defined by the methodology: “it is an X”, “it could be an X”, “it is not an X” and “unable to check whether it is an X”, where “X” stands for any of the four deformation detection processes executed by ADAclassifier; for instance, when talking about landslides, “it is a landslide” would be one of the values of the corresponding attribute.

Note that one of the four values is “unable to check whether it is an X”. This is because ADAclassifier may decide not to apply one or more of the four detection processes because of the lack of data. As stated above, a noticeable number of inputs is required. Just the inventory files already amount to five. Taking also into account the DTM and the horizontal components of the movement, it is easy to realize that in many cases the full set of files will not be available.

To limit this common problem, ADAclassifier makes almost every input file optional (this includes the DTM, the inventories, and the horizontal components). Each time the application is run, it analyzes the dependencies of each sub-algorithm and decides which of these may be executed as a function of the inputs provided by the user. Consequently, it is necessary to add the “unable to check whether it is an X” value as one of the possible outputs of each classification process.

Note that this behavior makes the tool much more flexible: while concentrating the detection of four deformation processes in a single application, it may be used to check just only one of these, just providing the available set of data for the target deformation process.

Figure 4 and Figure 5 depict the GUI of the ADAclassifier tool; the first one shows the options tab where the set of thresholds controlling the application are set; in the second one the list of (mostly optional) files may be appreciated.

As it happens with ADAfinder (and with los2hv too), a default options file may be set by the user to load a predefined set of threshold values. It is also possible to save/load option files created by the user. The ADAclassifier tool is fully described in its user guide [34].

THEXfinder, the second classification tool, does not work in the same way that ADAclassifier does. ADAclassifier takes the output of ADAfinder as input, i.e., it deals with ADA that have already been identified; starting here, it checks whether these ADA correspond to one or more deformation processes.

On the contrary, THEXfinder performs both tasks, that is: (a) identifying the ADA themselves (thus ignoring ADAfinder); and later (b) running the checks targeted at classifying such ADA either as expansive soils or thermal effects.

In fact, the first version of the toolset did not include THEXfinder as a separate tool. Instead, ADAclassifier was responsible for detecting not the current four, but the whole set of six deformation phenomena, including the two for which THEXfinder is now responsible. The availability of the first version of ADAclassifier made possible to attempt to classify ADA much faster, thus obtaining results more easily. Because of this extra abundance of output data, it was discovered that the mechanism to identify ADA used by ADAfinder was not appropriate when the phenomena to track were expansive soils or thermal effects. Therefore, the detection of these phenomena was removed from ADAclassifier and THEXfinder was created, implementing an appropriate detection mechanism for this kind of ADA together with the classification algorithms themselves. This means that THEXfinder identifies and classifies ADA in a single process.

The inputs for THEXfinder are the following:

• The original set of PS files as well as the read-map file defining the structure of the PSI file (see Section 3.4). Note that, in this case, the tool starts from the original PSI data set, not from the ADA. See the discussion above.
• An optional polygon defining the area of interest (shapefile).
• Optional ESRI shapefiles representing the infrastructures (buildings, bridges, etc.) and geologic inventories. In the case of the geologic map, an extra read-map file is also required.
• The parameters (thresholds) controlling the behavior of the application.

There are up to four output files, two for the expansive soils’ analysis and two more for the thermal expansion case. Each set contains one file with the ADA and another one with the points inside these ADA.

From the structural standpoint, the output files with ADA are almost identical to those created by ADAclassifier; practically, the whole set of attributes is identical. The only difference is that ADAclassifier includes a set of four fields to state the probability of an ADA matching each of the deformation processes tested, while the files created by THEXfinder contains only one of these fields per file. This is so because the ADA files created by THEXfinder correspond to just a single deformation process.

In the same line that ADAclassifier, THEXfinder will only try to check whether an ADA corresponds to some deformation process when all the input files required to do it exist. This is the reason why most of the input (and output) files are optional, making the tool more flexible.

Figure 6 and Figure 7 show the interface of the THEXfinder application.

2.3. Los2hv

The los2hv tool is targeted at the computation of the separate East-West horizontal and vertical components of the ground displacement measured with PSI along the satellite’s LOS. These components are inputs required by ADAclassifier to execute some of the classification algorithms it implements. Both ascending and descending datasets are required.

The tool accepts (input) and produces (output) ESRI shapefiles for compatibility reasons. On output, los2hv creates two files, to store, respectively, the East-West horizontal and vertical components of the movement as observed along the LOS. It should be clarified that North–South displacements are not calculated since, due to the nearly North–South orbit direction of SAR satellites, InSAR is insensitive to this displacement component [35]. It is noteworthy that this limitation is associated with the InSAR information itself instead of the los2hv tool.

Figure 8 shows the GUI of los2hv. As usual, the values shown are obtained from a default options file; also, the user may save/load their own options very easily.

los2hv performs a tessellation of the whole area covered by the PS included in the two input files. The size of the tesserae (grid spacing) is decided by the user. Consequently, each PS belongs just to one of the resulting tesserae.

There may be tesserae where: (1) there are no PS; (2) there are only PS from the ascending input file; (3) there are only PS from the descending input file; and finally (4) there are PS from both input files (see Figure 9). For those tesserae of type (4), the ground movement is averaged using the values of all the PS included in the tile. The resulting amount, measured along the LOS, is then converted to East-West horizontal and vertical components according to the formulae described in [36]. The value of the East-West horizontal and vertical components of the ground movement corresponds now to the whole area covered by the tile where the points involved in the computation where located.

These values are saved to the corresponding output files. Note that the user may select to store points (centroids of the tesserae with data) or polygons (more precisely, squares, the boundaries of the tesserae). The second kind of output (squares) is the one used by ADAclassifier. The centroid output is provided for whatever other purposes.

Please refer to the complete user guide of the los2hv tool [37] for more information.

3. Implementation and Integration

3.1. The Language of Choice

The toolset has been implemented in C++ to boost performance. The ADAfinder tool, for instance, may need a sensible amount of resources when working with large datasets (see Section 5 for details on performance of the whole toolset). Therefore, the use of a compiled (not interpreted) language was of capital importance. Other popular languages, as Python, have been avoided precisely for this reason, despite their popularity.

Although developed using the C++ compiler included in Microsoft’s Visual Studio, special precautions have been taken to make the source code portable, particularly for the most popular C++ compiler used in the Linux operating system, gcc. The use of Qt (see below) is also a factor contributing to the portability of the code.

Several open source libraries were used to implement the toolset. These are:

• Qt (see [38]). Although it has been used with several purposes in mind, the main target was to guarantee portability. Since the applications have a GUI, it was very important that such GUI was built using a portable library to avoid the need to write different code for each of the platforms which these tools are targeted at (at least Windows and Linux). Qt is a framework that guarantees such portability; in fact, developing cross-platform applications is its motto.
• Shapelib. This library is a very convenient tool to read and write ESRI shapefiles. See [39] for further details.
• Clipper. A library available for the Delphi, C, C+++ and Python, used for clipping and offsetting lines and polygons. For a complete description of this library, please refer to [40].
• Dlib [41] and Eigen [42] to implement some mathematical algorithms required by the detection processed included ADAclassifier and THEXfinder, such as curve or plane fitting.

3.2. The Three Incarnations

The tools described in this paper may be used in several quite different work environments. For instance, a user could use the ADAfinder application repeatedly to play with the parameters controlling the algorithm and then decide what would be the best strategy to identify the ADA in some area. This would imply the use of an ergonomic tool, where changing such parameters should be extremely easy and proof safe. On the contrary, once such parameters have been found, ADAfinder could be used routinely, with no human intervention, to detect ADA as one more step in an automated batch process. In this case, the parameters would be mostly the same for the whole dataset, and there would need to change only a few of them. There are other situations where the detection of ADA could be seen as a part of a much bigger process, being thus interesting to be able to embed the logic implemented in ADAfinder (or in any of the other tools) in another application.

These are the reasons why each of the applications described in Section 2 are available in three different flavors or incarnations:

• As a C++ class (one for each application) in a library. Third party (C++) software willing to embed the logic of ADAfinder, ADAclassifier, THEXfinder, or los2hv as a black box only needs to instantiate the corresponding class. Thus, embedding the necessary logic to be able to identify or classify ADA or to compute the horizontal components of the movement is just one procedure call away. Only software components developed in C++ will be able to integrate the logic in the library, since no bindings for other languages have been developed.
• As a command-line utility. This makes possible to integrate these tools in batch workflows, since no human intervention is required to run them. See Section 3.3 for details on options files, the mechanism used to obtain the information controlling the behavior of the applications.
• As an application featuring a GUI. This flavor is the best one for experimenting because of its ease of use. GUI-based applications, however, cannot be integrated in batch workflows.

There also exist pluggable versions of the three applications for the Quantum GIS (QGIS) version 2 tool. These plugins just call the GUI versions of the toolset thus avoiding the need to quit the GIS environment when running any of its applications. The tools are connected by means of some glue code written in this case in Python; since this is the unique choice when working with QGIS. Note that the plugins have not been yet migrated to QGIS 3 due to the changes affecting how this kind of software must be built for the latest version of this tool.

It is worth remarking that the command-line and GUI versions of the applications are just interfaces calling the classes in the library that actually implement the logic of the tools. This approach allows for a simplified maintenance process; the logic is concentrated in just one place, no matter how this logic is used (library, console or GUI application). Changing the classes implies an immediate update in the three flavors of each application.

3.3. Option Files

Both the command-line and GUI incarnations of the applications in the toolset rely on option files to retrieve the information defining how to behave; this includes input or output files and thresholds controlling some condition, among other data. Note that this is so for the GUI-based tools too; in fact, the GUI is just a mechanism to fill the gaps in an option file template, the so-called defaults files. This simplifies the design of the classes implementing the logic, since only one interface (the option file) needs to be taken care of. The command-line incarnations of the three tools have a single parameter: the name of the options file with the program’s parameters. The option files used by the toolset are uncomplicated, plain text files including couples of labels and values. Comments to clarify the purpose of each of these couples may be included just by adding a leading “#” character preceding the descriptive text. Figure 10 depicts a real, quite self-descriptive options file used to control the los2hv application.

3.4. Real-Life Shapefiles: Read-Map Files

ESRI shapefiles, although standardized, may include variable sets of attributes. Even when a shapefile includes the full set of attributes needed by a tool (for instance, the x and y coordinates as well as the velocity and the deformation time series in the case of ADAfinder), they may appear in different columns of the attribute (.dbf) file. This usually depends on the provider of the files.

This variability might become a serious problem, since the input module should be adapted for each kind of shapefile to process.

The solution to avoid this problem are the “read-map files” defining how the relevant attributes in a shapefile are organized. Of course, the attributes needed by the tools to work properly must always be present in the files, but thanks to these read-map files it is possible to deal with changes in the positions where these appear.

Thus, a read-map file is just an extra options file (see Section 3.3) where, by means of couples of labels and values, the positions of the attributes on which each tool rely are specified. Figure 11 is an example of an ADAfinder read-map file for the input PS.

The meaning of the pairs or labels and values in Figure 11 is the following:

• the x-coordinate of the PS must be read from column 5 in the .dbf file,
• the column to read to obtain the y-coordinate is the sixth one,
• the velocity may be found in column number 9 and, finally,
• the set of values making the time series start at column number 11 and there is a total of 50 of these values.

Obviously, if a shapefile with input PS is organized in a different manner, the values in the example above must be changed to match the actual positions where the fields of interest lie.

Note how such a simple mechanism makes the toolset much more flexible and able to deal with real-life datasets.

4. Quality Assurance

Prior to the delivery of the tools to the MOMIT consortium, these went through a rigorous testing process to guarantee that they worked correctly. Note that both a test plan and a test report concerning the whole toolset are described, respectively, in [43,44]. There, the full details on how testing was devised and took place are provided. However, and for the sake of completeness, the next sections briefly present the most relevant steps taken to guarantee the quality of the code.

Note, however, that the tests described in this section do not include THEXfinder as the independent application that it is now. This is so because the testing process took place before ADAclassifier was split in two (see Section 2.2). Consequently, only ADAclassifier, and not THEXfinder, was tested.

This apparently poses a question about the validity of the tests related to ADAclassifier, since it took care of classifying two extra phenomena (expansivity, thermal effects) that should not be implemented by this tool. The problem is caused by the mechanism used by ADAfinder to identify the ADA. It has been said (see Section 2.2) that such an algorithm is not appropriate to find ADA that are the result of expansivity or thermal phenomena. Therefore, if ADA affected by these phenomena cannot be properly identified by ADAfinder, it will be impossible for ADAclassifier to catalogue these properly.

Although this is essentially true, the tests applied to ADAclassifier where synthetic ones (see the description of the tests for ADAclassifier later in this section). Data were prepared in such a way that it was possible to predict the results of testing the six kinds of phenomena, assuming that ADA had been properly identified before. Consequently, the problem related to the appropriate identification of ADA did not invalidate the tests for ADAclassifier. This makes the tests for the specific algorithms for the identification of expansivity and thermal effects valid, since these were executed using correct (although synthetic) input data. This means that the results that were considered correct at that moment are still valid, including the two specific cases for which now THEXfinder is responsible.

THEXfinder has not been fully checked yet, however. From the discussion above it should be clear that the specific algorithms to diagnose expansivity or thermal effects are correct, since these were checked as components of ADAclassifier; nevertheless, the procedure used by THEXfinder to detect (not classify) ADA has not been rigorously tested yet. This is also the reason why performance results for this tool are not given in Section 5.

Going back to the tests themselves, in the case of ADAfinder, a manual methodology had been used for some time in GIS environments when this tool was developed, so datasets including both inputs and outputs were available. Therefore, the tests consisted essentially of comparing the results of the manual procedure with those created by the tool.

A mechanism to quickly compare the results produced by the manual and automated solutions was devised. From the numerical standpoint, it basically consisted of exporting the values of the attributes to check for both outputs (manual and automated), sorting these to easily match the attributes in each file and then computing the differences of their values, which, in all cases, were under the threshold set by the precision of the typical 8-byte IEEE 754 double (around the 15th decimal position). This means that the results (for instance, the coordinates of the output points) were equivalent.

To check the correctness of the algorithm, two kinds of tests were performed. Firstly, the values of the attributes standing the level of certainty were checked for absolute equality, since these were represented by integer (non-floating point) magnitudes. Then, the number of ADA and their contents were also checked. This means that both the manual and automated algorithm had to identify the same number of ADA and that the set of points included in each of them had to be the same. Both numerical and correctness tests were passed satisfactorily.

The situation for ADAclassifier and los2hv was different; no previous results existed, so no reliable dataset to compare their outputs was at hand. The approach, therefore, was to create synthetic datasets for both applications. The rationale behind the way these datasets were created was to organize the information in easily identifiable geometric patterns, so, when combined, the area(s) where positive results were obtained was (were) also predictable.

For instance, and to check the landslides algorithm in the ADAclassifier application, all ADA were created as identical square-shaped polygons containing a regular grid of 5×5 PS. These ADA were distributed in a 24×24 checkerboard pattern (see Figure 12a); the horizontal components of the movement were distributed in 4 adjacent vertical stripes covering 6×24 ADA each. A known value for each of them was set (respectively, from left to right, 0, 10, 20, and 30 mm/year). This is shown in Figure 12b. Other easily identifiable patterns were used for the rest of input files; such as the values of the slopes in the DTM.

Then, for example, if the landslide algorithm imposes that the value of the horizontal component of the movement must be higher than 25 mm/year, then only the rightmost column made of 6 × 24 ADA will meet the requirement. This may be seen in Figure 12c where the ADA and the horizontal components of the movement (Figure 12a,b) have been overlapped.

A similar setup was used to check the los2hv application. In this case, the area to test was tessellated using a 7 × 7 grid. There, two groups of PS belonging to the ascending and descending datasets were arranged in alternate rows and columns so only in a well-known subset or the tiles defined by the checkerboard of PS coming from both datasets would coincide. See Figure 13 for a graphical depiction of this setup. There, the tiles surrounded by a green square are the only ones containing both ascending and descending PS and therefore the unique places where the horizontal and vertical components of the movement may be computed by los2hv.

Known values were set for the velocity (v_a and v_d for ascending and descending datasets respectively, see again Figure 13) and the (ascending/descending) deformation time series (def_a and def_d) for each PS were assigned. This made possible to compute manually (and easily) the results that should be expected from los2hv and thus validate its correctness.

5. Performance Evaluation

Only one of the applications developed and tested, ADAfinder, relies on a procedure for which previous results already existed. Therefore, it was the only one for which a performance reference was available for comparisons. Since such a procedure was executed manually (a series of steps performed by an operator using the tools offered by a GIS), a noticeable improvement of performance was expected due to automation.

There were no performance references for los2hv and ADAclassifier (Section 4 explains that when the tests took place, ADAclassifier took care of the detection of expansivity and thermal effects. Now, these two algorithms have been moved to THEXfinder. This implies that the performance results given for ADAclassifier in [44] are, actually, pessimistic, since these include the time needed to check two extra phenomena that now are not the responsibility of this tool.), so no specific expectations about performance improvements existed. All tests took place using a computer with the following characteristics: Windows 10 64-bit, Intel Core i5-5300U @ 2.3 GHz, 2 cores, 4 threads, 8 Gb RAM, 500 Gb magnetic (non-SSD) hard disk.

Table 1. Performance of the toolset.

Tool.	Dataset	Time (s)
ADAfinder	20,351 PS	2
ADAfinder	926,916 PS	179
los2hv	2 (ascending, descending) × 135 PS. Grid: 7 × 7 tesserae	55
ADAclassifier	144 ADA, 3600 PS Between 4 - 8 polygons per inventory DTM with 14411441 z values	125

Note that there are no performance results for THEXfinder yet; this tool has not been fully tested at the moment this work was published. Note as well that in the case of ADAclassifier, the process includes the identification of the whole set of deformation phenomena (that is, the test covers all the possible classification processes available). Testing for fewer options will reduce the time needed to process data.

shows the performance of three of these tools.

6. Availability

Although the authors are not yet ready to offer the set of tools described herein as an open source project, it is possible to obtain, strictly on demand, a free, executable version (Windows) of ADAfinder, ADAclassifier, THEXfinder, and los2hv directly from the authors, for both commercial and non-commercial purposes. Please, contact the first author of this work for more details.

7. Real Test Cases

In this section two case studies are shown to illustrate the tools presented above. Unfortunately, THEXfinder has not been applied to these areas, since this tool has not been fully tested yet (see Section 5). In the near future, a new paper will deal with THEXfinder applied to some area of interest.

The first case, which has been processed and studied in the framework of the RISKCOAST project (SOE3/P4/E0868), is an area of around 100 km² located in the coastal area of the province of Granada, in southern Spain (Andalusia), including Salobreña and Motril towns. For this test site a limited number of auxiliary data were available for the classification analysis. The second area, for which most of the optional information for ADAclassifier was available, is located in southeastern Italy, between Tropea and Zaccanopoli, in the Calabrian coast of the province of Vibo Valentia, and extends for 12 km².

Figure 14 depicts the location of these two test areas.

7.1. Southern Spain

7.1.1. Input Data

The InSAR-derived displacement map, which is the main input of the ADAfinder, has been generated applying the PSIG approach using the processing tools developed from the Geomatic Division of CTTC. The processing allowed the annual velocity and the displacement time series to be estimated; refer to [45] for more details on the PSI processing. A set of 138 images acquired by C-band Sentinel-1 satellites (developed by ESA for the Copernicus Programme), in ascending geometry and Wide Swath acquisition mode, has been processed at full resolution covering the period from March 2015 to September 2018. The resolution of Sentinel-1 data is approximately 4 × 14 m² and the temporal sampling is 6 days. The auxiliary data used for the ADAclassifier are: (1) the digital terrain model of the project PNOA-LIDAR (from the National Center of Geographic Information; CNIG) with a resolution cell of 5 m, (2) the geologic map from the Instituto de Estadística y Cartografía de Andalucía, at 1:400.000 scale; (3) the cadastral data from the Spanish General Directorate for Cadaster to select the urbanized and construction areas; and (4) the Corine Land Cover 2018 to complete the cadastral inventory information about structure and infrastructure areas.

7.1.2. Results

The input shapefile (PSI displacement map) for the ADAfinder includes 61,802 measured points. From this map a total number of 82 ADA have been extracted, where 32 have QI = 1, 20 have QI = 2, 13 have QI = 3 and 17 QI = 4 (Figure 2). The QI class allows the user to easily understand the noise level (i.e., the reliability) of the displacement time series (TS) within each extracted ADA. Based on the QI classification we have decided to be more restrictive by selecting only the ADA with higher reliability of TS information (QI equal to 1 or 2, for a total of 53 ADA) for the ADAclassifier, where the decision of whether or not the ADA is a potential settlement is based on the mean TS trend. For the output classification (Figure 15 below and

Table 2. Summary of ADA extracted with QI < 3 and classified in the test areas.

Dataset	# ADA	# Active PS			Area ADA (m2)			Landslide			Subsidence			Settlement			Sinkhole
		Avg	Min	Max	Avg	Min	Max	C	P	N	C	P	N	C	P	N	C	P	N
Southern Spain
ASC	53	19	5	255	13,216	3620	134,852	0	41	12	0	17	36	6	1	46	0	0	35
Southeast Italy
ASC	38	8	5	23	92.2	52.0	312.3	7	6	25	0	14	24	2	0	36	0	0	38
DESC	133	11	5	97	124.8	51.4	704.8	66	20	47	0	28	105	13	0	120	0	0	133
							TOTAL	73	26	72	0	42	129	15	0	156	0	0	171

on page 22), we want to underline that very few auxiliary data have been selected as input. For example, no landslide or subsidence inventories have been used, nor the horizontal-vertical decomposition was available (we have worked with only the ascending geometry of acquisition). For this reason, for both landslide and subsidence phenomena, only results tagged as “potential” will be possible at most. The minimum slope angle to be a potential landslide has been set to 5 degrees, while the maximum slope angle to be a potential subsidence has been set to 10 degrees. We remind that the considered slope is the mean slope value within the ADA. The different slope thresholds allow a double classification of the ADA that are at the limit between the 2 potential phenomena. Moreover, an ADA is classified as potential subsidence only if it lays on a quaternary lithology. For what concerns the consolidation settlement due to new constructions, 6 ADA have been classified as settlement due to the clear inverse exponential trend of the mean TS and the intersection with polygons related to infrastructures and urban areas. One ADA has been classified as potential settlement because it presents an inverse exponential trend, but it does not intersect with any polygon of the cadastral inventory. Among the settlement ADA, 5 lay on the recently built A-7 railway and one on a building close to the port of Motril; both structures were built between 2014 and 2015. For the ADA classified as potential landslides, some of them are already known slope instabilities affecting the coast of Granada [46,47,48], in this test area movements affecting the urbanizations of Los Almendros and Alfamar are included.

We consider that a multiple classification of the same ADA lets the final user know that the detected movement is a complex case. For example, the ADA classified as both potential settlement and landslide needs to be further analyzed, as it could be a landslide movement that has been stabilized (causing the inverse exponential trend), or a settlement due to a construction not yet inventoried.

7.2. Southeastern Italy

7.2.1. Input data

In this case, the InSAR-derived deformation maps were obtained using the Persistent Scatter Pairs approach (PSP, [49]). The processing involved three COSMO-SkyMed Stripmap frames, one ascending and two descending; X-band COSMO-SkyMed Stripmap images have a ground resolution of 3 × 3 m. The ascending frame is composed of 71 images, acquired between March 2013 and October 2018 with an average incidence angle of 29°. The first descending frame is composed of 29 images acquired between March 2017 and March 2019 with an average incidence angle of 26°; the second descending frame is composed of 94 images acquired between October 2011 and October 2018 with an average incidence angle of 29°.

Complementarily, a Digital Elevation Model (DTM), a landslide inventory map, a geological map and a land use map have been integrated into the classification process performed by ADAclassifier. The DTM consists in a 5 m cell resolution of the Calabria Region. A geological map –scale 1:25.000– and the Italian Landslide Inventory map—IFFI project—of the study area have been used also used. Finally, the Corine Land Cover, Level IV has been used to map urban areas.

7.2.2. Results

A total of 38 and 133 active deformation areas (ADA), all them exhibiting a QI = 1, were detected, respectively, for ascending and descending InSAR datasets. The sizes of the ADA vary from 51.4 to 704.8 m² and they are mainly located in the urban areas and the reliefs located at the southeastern part of the analyzed area.

The classification of the ADA identified using ascending and descending datasets by means of ADAclassifier shows that the ADA are classified as: (1) landslides (42.7%) and potential landslides (15.2%); (2) potential land subsidence (24.6%); and (3) consolidation settlement (8.8%). No potential sinkholes were identified by ADAclassifier. Table 2 summarizes the obtained results, depicted in Figure 16.

The landslides are mainly located in the areas with a certain slope in which the inventory map indicates the presence of gravitational processes. Most of the landslides (66 and 7, for descending and ascending frames, respectively) have been confirmed since they are known phenomena included in the IFFI database. Complementarily, 26 ADA (20 and 6, for descending and ascending frames, respectively) are potential landslides since they fall within areas with slopes higher than 10° and/or exhibit not-negligible horizontal displacements. A literature review confirms that this area is particularly prone to weathering processes, lateral spreading, and landslides [50].

In the urban areas of Tropea and Zaccanopoli, 42 potential subsidence areas (28 and 14, for descending and ascending frames, respectively) have been identified. These ADA are located over very flat areas geologically composed of Quaternary sediments. Some of these ADA classified as subsidence areas can be related to the overexploitation of underground water resources, typical for the coastal plains of the Calabria Region [51,52]. The ADA classified as consolidation settlements are mainly placed in the harbor and in some specific locations of the urban area. These areas fit an inverse exponential function with a determination coefficient (r²) higher than 0.8 and thus are probably associated with post-constructive consolidation settlements of the dikes of the harbor and of some recent buildings.

8. Discussion & Conclusions

This work presents the set of tools integrated by ADAfinder, losh2v, ADAclassifier and, to a lesser degree, THEXfinder, since this tool could not be fully tested at the time this work was published. These applications are targeted at automating the identification and classification of ADA. All of them rely on methodologies that had been thoroughly used in real (production) use cases prior to their materialization as software modules. Their implementation is based on well-proven techniques. Although the tools—not the underlying methodologies—were originally conceived to track ground deformation phenomena close to railway infrastructures, they may be used in any other context, as presented by other authors [28,29,30].

The applications have been sufficiently described in enough detail, considering aspects such as the algorithms implemented, inputs, outputs, or the way these may be used (GUI or console incarnations, or classes in an embeddable library). The procedures to guarantee correctness, as well as the performances of ADAfinder, los2hv, and ADAclassifier, have also been discussed.

Note that, especially in the case of the so-called classification tools (ADAclassifier and THEXfinder), the ability of these applications to work with a big number of optional inputs (i.e., inventories of several kinds) leads to a high degree of flexibility, making the tools useful in a wide variety of situations that are usual in real-life production environments. In other words, ADAclassifier can deliver results for just one deformation phenomenon or go up to four of these, depending on how many optional inputs are provided. The flexibility of the ADAclassifier is displayed in the two use cases (see Section 7.1 and Section 7.2), where the availability of input data for the classification process was completely different. The situation in the Spanish use case (Section 7.1), where only a few ancillary data are available, is unfortunately more usual than that of the Italian case (Section 7.2), for which a noticeable number of data sources were available. The fact to note here is that ADAclassifier (and THEXfinder) is (are) able to deal with these two extreme situations.

It is worth noting that the results provided by these tools are affected by the drawbacks inherent to the InSAR technique used as input, such as the loss of coherence and the low sensitivity to measure North–South directions displacements.

No less important is the fact that automating preexisting methodologies by means of reliable tools (Section 4), provides good response times (Section 5) which may be completely controlled by the user via a whole set of editable parameters, paving the way to better analysis of the areas of interest. Every area is different, and adjusting the parameters controlling the behavior of the process to obtain reliable results is therefore a task that depends on where such areas are located. Avoiding error-prone manual processes thanks to these automated tools, and having results in a few seconds or minutes, make it possible to assess whether or not the starting set of parameters is or is not the most appropriate, it thus being possible—and fast—to repeat test after test at almost no cost, thus increasing the quality of the results. This is a technique that is also applied when using manual procedures, for instance, relying on GIS environments; therefore, being able to do this at a much lower cost and risk allows to obtain better results in less time.

Concerning the real use cases, the results obtained illustrate the applicability of the tools developed to automatically and successfully identify and classify the geological and geotechnical processes affecting wide areas of the territory. Landslides, potential landslides, consolidation settlements, and potential land subsidence have been identified in both test sites using inputs with different characteristics (i.e., InSAR and DTM resolutions, inventory maps, and land use maps). Thus, the possibility of using different input data reveal a high versatility of the proposed tools. Furthermore, this methodology has the advantage that results can be periodically updated by incorporating new InSAR datasets and updated ancillary information.

The outputs provided by ADAtools consist of a set of maps that can be of high interest for geological hazard management to be incorporated into land management and planning.

Moreover, ADAtools can be considered a first step towards satisfying the need for post-processing tools due to the increasing use of PSI data at a regional, national, and European level [3]. The ADAtools will allow different stakeholders to have a fast selection and preliminary interpretation of deformation maps composed of millions/billions of points. This will reduce the amount of data to be managed in terms of hardware, space, and time of analysis, thus increasing the operational use of PSI displacement maps over wide areas. To finish, it is worth noting that the ADAtools may be obtained on demand, completely free of charge, for both commercial and non-commercial purposes (refer to Section 6 for details).