Sequential Earthquake Damage Assessment Incorporating Optimized sUAV Remote Sensing at Pescara del Tronto

A sequence of large earthquakes in central Italy ranging in moment magnitudes (Mw) from 4.2 to 6.5 caused significant damage to many small towns in the area. After each earthquake in 2016 (24 August and 26 October), automated small unmanned aerial vehicles (sUAV) acquired valuable imagery data for post-hazard reconnaissance in the mountain village of Pescara del Tronto, and were applied to 3D reconstruction using Structure-from-Motion (SfM). In July 2018, the site was again monitored to obtain additional imagery data capturing changes since the last visit following the 30 October 2016 Earthquake. A genetic-based mission-planning algorithm that delivers optimal viewpoints and path planning was field tested and reduced the required photos for 3D reconstruction by 9.1%. The optimized 3D model provides a better understanding of the current conditions of the village, when compared to the nadir models, by containing fewer holes on angled surfaces, including an additional 17% surface area, and with a comparable ground-sampling distance (GSD) of ≈2.4 cm/px (≈1.5 cm/px when adjusted for camera pixel density). The resulting three time-lapse models provide valuable metrics for ground motion, progression of damage, resilience of the village, and the recovery progress over a span of two years.


Introduction
A sequence of large earthquakes caused considerable damage to multiple cities surrounding the faults in the Central Italian Apennine Mountain Range [1][2][3][4][5][6][7]. Analysis of high-magnitude earthquake sites assists post-disaster recovery and preventative efforts; an example in Iran is shown by Karimzadeh et al. [8]. Images of the affected area were collected by small unmanned aerial vehicles (sUAV) and 3D digital models of the small mountain village of Pescara del Tronto were made using Structure-from-Motion (SfM). Three models of Pescara del Tronto represent the post-disaster conditions after the August 2016 earthquake, after the October 2016 earthquake, and in July 2018, two years after the series of earthquakes; these models narrate the progression of damage and recovery. Volumetric measurements, movement detection, post-hazard reconnaissance, exterior safety inspections of buildings, documented recovery developments, and more can be gleaned from the sequential 3D reconstructed models. Pescara del Tronto presents a case study for comparing sequential changes of the city, comparing automated grid and optimized camera location methods, and an introduction to additional analysis (statistical and geotechnical).
The 2016 Italy earthquakes devastated communities along the Apennine Mountain Range along a series of normal faults with rake angles of −80 • to −100 • . The August 2016 earthquake caused significant damage in the following villages: Arquata del Tronto, Accumoli, Amatrice, and Pescara del Tronto [4,5]. The October events caused significant damage primarily in the villages of Visso, Ussita, and Norcia [2,6], and exacerbated the damage in Pescara del Tronto. Other analyses of the 2016 earthquake sequence include geohazard forensics of Norcia in Barone and Di Maggio [9] as well as numerical modeling of ground deformation in Valerio et al. [10], the mechanisms of the earthquakes in Zhong et al. [11], and Amatrice building damage assessment using synthetic aperture radar (SAR) in Karimzadeh and Matsuoka [12]. The earthquake events appear in a gap between two recent M6.1 events-the 1997 Umbria-Marche earthquake to the north-west and the 2009 L'Aquila earthquake to the south-east. Figure 1 from Stewart et al. [2] depicts the three major earthquakes in relation to the Umbria-Marche and L'Aquila earthquake events. This gap was determined to be an area of elevated risk prior to the 2016 earthquake events [4,5]. In total, a sequence of 17 major earthquake and aftershocks with moment magnitudes greater than 4.2 Mw were recorded in central Italy between August and November 2016. These events were initially recorded by the Italian National Seismic Network (Rete Sismica Nazionale, RSN) owned by the Italian Institute of Geophysics (Instituto Nazionale di Geofisica e Vulcanologia, INGV) near the villages of Amatrice, Accumoli, Norcia, and Visso [13,14]. Later, the number of recorded events increased to 49 by the United States Geological Survey (USGS) as recorded on the USGS web site (https://earthquake.usgs.gov/earthquakes/). Three main events identify the major earthquakes: M6.2 (24 August 2016), M6.1 (26 October 2016), and M6.6 (30 October 2016) [1][2][3][4][5] with the remainder classified as aftershocks.
Due to the temporal and geographical proximity of the series of earthquakes, the city of Pescara del Tronto was devastated by these events. The city lies within approximately 29 km of the all three main earthquake events and 8 km from the M6. 2 24 August earthquake event (see Figure 2). The peak ground acceleration (PGA) in city for the M6.2 24 August, M6.1 26 October, and M6. 6 30 October events were 0.48 g, 0.10 g, and 0.34 g, respectively [4,5]. The closest ground motion station to the events' epicenters was located at 42.962 • N 13.050 • E at a distance of 29.6 km from Pescara del Tronto, as recorded on the USGS web site for each event (https://earthquake.usgs.gov/earthquakes/). The M6. 2   To the authors' knowledge, the referenced Geotechnical Extreme Events Reconnaissance (GEER) reports follow international guidelines for post-disaster reconnaissance, but the scope of this paper is neither Rapid Mapping (hours to days post-disaster) nor Risk and Recovery Mapping (weeks to a couple months after the event) as described by Copernicus Emergency Management (https: //emergency.copernicus.eu/) [1][2][3][4][5]15]. As the sequential data is as far out as 2 years post-disaster, this paper is an additional ad hoc representation of the earthquake recovery at Pescara del Tronto.

Method of Preliminary sUAV Reconnaissance and Modeling in 2016
The use of sUAVs in disaster reconnaissance prevents researchers and surveyors from having to enter many high-risk areas and reduces the time spent collecting surveyed data at every desired point. A recent example of high resolution sUAV data acquisition is Pavelka et al. [16] of the Pista Geoglyph in Peru, and expertise with sUAVs (and LiDAR) as shown in Pellicani et al. [17] demonstrates that sUAVs support and compliment geomorphic interpretations of disaster sites. Creating useful and accurate models with SfM software requires a series of photos at various normal and oblique angles that eventually tie into ground control points (GCPs) in a 3D reconstructed model. Three methods were used to collect the required imagery set of Pescara del Tronto during the two 2016 reconnaissance missions: basic nadir grids planned using the eMotion software to control the eBee system, multiple manual flights, and autonomous flights through the GSPro app.
During the 2016 reconnaissance trips, multiple sUAVs collected the images for reconstructing 3D models of Pescara del Tronto. Due to widespread use and an easily programmable software developer kit (SDK), the main method of collecting the images at Pescara del Tronto was with DJI products. During the 2016 reconnaissance trips, the platform used for collecting images was a DJI Phantom 4 Pro TM quadrotor sUAV, equipped with a 4 K video camera with a 1/2.3" CMOS sensor, 94 • field of view, 12.4 MP images, and a focal length set to infinity [3,18] was used. The Phantom 4 Pro system weighs 1.388 kg, has a maximum flight time of 28 min. In addition to multiple manual flights, GSPro was used as the autopilot to fly the drone and capture the desired nadir grid images.
The group from the Politecnico di Torino Geomatics provided a SenseFly TM eBee. The eBee is a small fixed-wing sUAV manufactured by SenseFly at about 0.6 kg of platform weight, and can collect up to 8 km 2 of surface area and 40 min of flight time depending on payload for the reconnaissance mission. The equipped payload included a digital camera Canon Power Shot S100, 1/1.7" Canon CMOS sensor, focal length of 5.2 mm, and 12 MP images [6]. While the platform can be navigated manually by a ground pilot, a piece of mission-planning software called eMotion [19] pre-programmed the system with GPS aerial waypoints. In addition to the models developed from the DJI Phantom 4 Pro, the eBee's high altitude images provide an expanded view of the site and the surrounding terrain.
To georeference the models, GCPs were placed throughout the city and matched to multiple photos within the SfM software. Pescara del Tronto contained a total of 23 points with five 50 cm 2 fabricated markers and 18 natural GCPs. The markers were surveyed by a Real Time Kinematic (RTK) configuration using the GNSS Geomax Zenith 35 TM (www.geomax.com) [3]. Images showing the layout of the GCP markers and layout in Pescara del Tronto are shown below (see Figures 3 and 4 from the GEER report) [1].  After the images and ground control points were obtained during the 2016 reconnaissance projects, 3D models developed and recorded current conditions and detected deformations that took place between the two events [6] (see Figures 5 and 6 for the models). The models were made through Bentley TM ContextCapture modeling software which determines 3D point locations within a scene using GPS and camera metadata. The data is stored within the photos and used to calculate interspatial intersections-molding/stitching the photos together by identifying matching figures in the overlapping photos as in Koenderink and van Doorn [20]. SfM assigns point x, y, or z positions to similar features in each picture and spatially stitches each point together into a cohesive object [21]. However, the quality of the final product is generally determined by the skill and experience of the pilot constructing a grid overlap and obtaining additional images through manual flights.

Explanation of the Algorithms Implemented in 2018
Due to the mountainous terrain surrounding Pescara del Tronto, standard grid pattern photo gathering is inefficient and makes the collection of the required images very difficult [22]. Many flight planning applications have assisted autopilot software; however, flight path algorithms rarely consider occlusions, topographical features outside of elevation change, and do not change or optimize the six degrees of freedom data for camera positions in addition to the desired flight path. In post-hazard environments, prompt assessments and responses are essential. The time required to create a model is dependent upon the number of photos evaluated and number of tie points that are matched in a set of photos. Grid overlap flight paths not only take redundant photos, but require a pilot to estimate occluded areas, oblique positions, and take manual photos before generating an optimal model [20,23,24]. Moreover, locations similar to Pescara del Tronto have large hills, steep slopes, and cliffs that leave missing model sections with nadir grid flight paths. The nadir flight paths may go through a high unforeseen hill or may not collect the required images that are normal to a cliff or complex 3D surfaces. The optimized view and path planning algorithm developed at Brigham Young University (BYU) chooses optimal camera locations and oblique angles to capture the required images at various locations and orientations, thus, reducing the redundant photos and required time for 3D reconstruction. Genetic principles are incorporated into view-planning and field tested as part of the sequential earthquake reconnaissance case study. The basics of the discussed algorithms (assuming that the reader already has a general knowledge of combinatorial optimization) are followed by an explanation of the sUAV SfM work by Martin et al. [23] and is implemented at Pescara del Tronto.

Optimized View and Path Planning Algorithms
To overcome the downsides of the nadir grid approach and obtain better orthogonal views of Pescara del Tronto, a view-planning algorithm delineates optimized photo locations and angles through greedy heuristics [25]. View-planning identifies the best sensor locations for observing an object or site by planning oblique images for optimized SfM [26,27]. As shown in Martin et al. [23], when comparing models that use optimized view-planning with ones that use standard grid patterns, models using additional view-planning yield a greater accuracy of up to 43%. These results led to the use of algorithmic view-planning (specifically model-based view-planning) at Pescara del Tronto.
The sUAV SfM used in Italy builds off previous work by Palmer et al. [28]. Various computational algorithms optimize what the sUAV needs to "see" as well as the best route to do so [23,[29][30][31][32][33]. The path planning algorithm, follows a weighted Traveling Salesman Problem (TSP). Combining greedy and TSP algorithms effectively addresses non-deterministic polynomial-time hard (NP-hard) problems [34]. These computational tools organize the data points in a practical manner that produces photogrammetric results that look down vertically as well as at various oblique angles. The optimized approach is safer to people and equipment by taking into account the 3D terrain and produces more complete photogrammetric models. View and path planning algorithmic execution smooths out the sUAV mission and results (see Figure 7 for a basic graphic of the algorithms). The view-planning algorithm follows the basic structure of a greedy algorithm [35]. While relatively simplistic, it quickly decides from a histogram where an sUAV can take a picture of the environment at the few locations that are worth taking a picture. The greedy algorithm follows the structure shown in the left of Figure 7 demonstrating a sub-optimal (highlighted by the red circles) but still effective solution by taking the smaller number in each binary option. While quick and generally effective, the downside to greedy heuristics is that no attempt to remember the previous choices is made, so a sub-optimal scenario could be a chosen result.
The path planning algorithm solves a weighted TSP (see the right portion of Figure 7) [36]. The Christofides algorithm is used to solve the TSP to a certain guaranteed degree of optimization as explained in Hoppe [37] and is used despite further research into computational solutions to the TSP as in Bożejko et al. [38]. The TSP path planning replaces the nadir grid approach in a manner that optimizes resources.

Evolutionary View-Planning
The optimized method applied at Pescara del Tronto deals with evolutionary model-based view-planning and closely follows Martin et al. [23] using a priori terrain data, a fitness function, inheritance, diversity and elitism, and convergence.

A Priori Terrain Data
A priori data is obtained from iterative flights and/or rough elevation models from open sources such as Google Earth and the USGS [39,40]. From the data, camera angles (azimuth and image elevation angles while neglecting roll (analogous to pitch, yaw, and roll with 1 less degree of freedom)), and X, Y, and Z coordinates are determined from an outlined polygon of the area of interest. The optimized approach directs the sUAV and obtains oblique angles by using the camera gimbal [41]. Figure 8 from Martin et al. [23] graphically shows a tiny data set (histogram) of possible camera locations and angles. At Pescara del Tronto, the similar histogram contains thousands of potential views with more noticeable camera angles.

Fitness Function
The histogram data is organized following fitness function f (see Equation (1)) as an objective where C is coverage (total visible points) and Fc is functional coverage (a percentage that describes what the sUAV needs to "see" on its mission) [23]. For coverage C (see Equation (2)), v is visible points from the images, and T is the points describing the total sUAV mission area. For functional coverage Fc (see Equation (3)), P is points that will be "seen" by the sUAV, and T is the same as in Equation (2). In other words: P is the points that are captured by at least 3 camera shots; Fc is a normalized subset of C as P is divided by the points representing the entire terrain space T. The objective seeks to choose a minimum subset of the histogram data set, from the reduced histogram we generate sUAV camera locations and poses to maximize model quality (in accordance with similar procedures [42]).

Inheritance
To achieve the optimal sUAV camera locations and orientations, the genetic algorithm uses inheritance as in Martin et al. [23]. From the histogram, two random image solutions become mother and father, undertake a blend crossover with 0.7 probability, and the weighted average becomes the new image solution. Equations (4) and (5) further describe image crossover for each image variable ( 1 5 where each variable has 0.5 crossover probability) and r is randomly chosen between (0, 1). The inheritance continues until there are two new image solutions sets.

Diversity and Elitism
To ensure diversity in the population (image sets) simple random selection prevents premature solutions (that occur from a seeded tournament selection). A 5% chance to mutate occurs during crossover on each variable (akin to a chromosome mutating [43]). A dynamic mutation parameter of α = 1 keeps initial values in variable bounds. An intermediate variable a accompanies the dynamic mutation (see Equation (6)) where g is current generation number and G is total generation number. Using the same example as in Martin et al. [23] using Equation (7), E is current image easting, E new is new image easting, E max is maximum allowed northing, E min is minimum allowed northing, and r e is a random number between (E min , E max ).
The diversity of the image solution set is culled by elitism towards optimal results. Traits of better parents are weighted higher than lesser parents such that the better 25% of solutions replace the lesser 25% of solutions in each child solution set (e.g., an angle of a camera pose that views more points is more likely to propagate across generations). Preserving better potential solutions leads to convergence.

Convergence
Every generation of the algorithm is tested for convergence using Equation (8) from Martin et al. [23] (with appropriate parameter tuning). In the convergence equation, a c averages the current generation score and m is the mean of the average scores from the previous ten generations.
Upon convergence, the solution removes obsolete images that cannot "see" enough terrain points and grooms the image set by removing image locations that pass through or are under terrain. Cleaning the solution image set increases safety and avoids potential collisions that may have been generated. Once optimal, the matrix of remaining solutions is ready for the TSP as in Section 3.1, and, once completed, the data set is ready for a sUAV mission.

Why Evolutionary View-Planning Matters
A priori terrain data, fitness function, inheritance, diversity and elitism, and convergence produce a superior image set from genetic principles. The main downside of the optimized algorithms is that the use of evolutionary view-planning in the field is a bit slower than a manual approach (but is faster than a comparable grid approach). This is because of the required time at each camera location is longer due to stopping to stabilize while positioning the gimbal each time. However, even though the optimized approach requires fewer photos, the sUAV "sees" more than the grid approach that has the camera pointed in a single direction at all times. The optimized method captures pictures in areas that need extra detail rather than taking pictures at every location (which fetches excess data).

Application of the Algorithm to Pescara del Tronto in 2018
Pescara del Tronto serves as a prime location for testing the optimized path and view-planning algorithms. The 2018 visit to Pescara del Tronto yielded three models: one representing the data collected using an industry standard automated nadir flight planning application (a lawnmower approach with the camera pointing down); another model using the view-planning algorithm [44], and the final model created by combining all collected images, including those from manual flights, into a single model-the third model is not included in further quantitative discussion as the resolution is a skewed average of the other two models that use many more resources, but the combined model is used in Sections 5 and 6 for qualitative comparisons. A sequential comparison of the 2016 and 2018 models is discussed in later sections of this paper.
During the July 2018 trip, an additional sUAV was used to collect images for reconstructing 3D models of Pescara del Tronto. A DJI Inspire 2, equipped with a Zenmuse X4S 4K video camera with a 1" CMOS sensor, 84 • field of view , 20 MP images, and a focal length of 24 mm to collect the images for the nadir grid model. The Inspire 2 system weighs 3.44 kg, and has a maximum flight time of 27 min with the gimbal and Zenmuse X4S payload. Originally, the plan was for the Inspire 2 to fly the optimized view and flight path due to the Zenmuse X4S being significantly better than the Phantom's camera. However, due to unforeseen events, the Optimized View and Flight Path algorithms were used on the Phantom 4 Pro with the same payload as the 2016 flights (all specifications from the cameras are kept the same except for pixel density).
Multiple ground control points were surveyed over the span of a week at Pescara del Tronto. During the first few days, images collected by using the nadir approach and manual flights were obtained in addition to 22 GCPs. Later in the week, the site was flown with the optimized view and path planning algorithms and was surveyed with 12 new GCPs (see Figure 9). The different number of GCPs was due to different groups performing the surveys because the GCPs were removed before the optimized view and path planning algorithms were used in the later part of the week.
The surface used to generate the cameras positions for Pescara del Tronto combines a Google Earth keyhole markup language (.kml) file and a point cloud (.ply) file of the area. The .ply file was obtained from the previous models of the area in 2016. By using the point cloud in addition to the .kml file, the view-planning algorithm identifies better camera locations for the flight compared to generated optimal locations based solely on the .kml file.
In 2018, the Phantom 4 Pro carried out the optimized mission and the DJI Inspire 2 flew the nadir grid mission. The Inspire 2 had a Zenmuse X4S camera attached that rates up to 4K like the Phantom 4 Pro camera but with 20 MP instead of 12 MP [45]. Figures 10 and 11 demonstrate some of the camera locations for the nadir grid and optimized mission plans, respectively. Once the drone captured all the images at the desired locations, each image with the georeferenced metadata and additional GPS data was processed by Bentley's ContextCapture 3D modeling software [46]. Tables 1 and 2 summarize quantitative results using the GSPro autopilot (a typical industry nadir grid) and the optimized approach (targeted automated flight).    Table 1 reviews GCP, photo count, Surface Area, and ground-sampling distance (GSD). GSD is a measurement of model resolution derived from pixel density and camera distance from the target. The model created by the optimized view-planning algorithm, uses 30% fewer photos per square kilometer while providing a more detailed model. The 3D model made from the optimized view locations had significantly fewer holes and a greater surface area by 17.1%. Even though the nadir 3D model's GSD is slightly better than the optimized equivalent, it is because of camera quality. This makes it difficult to compare the two models from the GSD. As the Inspire 2 is ≈ 167% pixels compared to the Phantom 4 Pro, the adjusted GSD (as if both flights were done with the Inspire 2) based on MP differences would be 1.5 cm/px for the optimized approach and 2.0 cm/px for the nadir grid mission. The differences in camera qualities explain the discrepancy of the nadir-based results, and if projected would have higher average GSD than the optimized flight. However, even though both models have similar GSD, an analysis that does not include the GSD follows (see Table 2). Table 1. Data from Nadir Grid-Based and Optimized Models. The quantitative comparison in Table 2 combines a four by five grid over Pescara del Tronto. At intersections of the grid, the photo count is found. Areas further up the slope have fewer photos in the nadir model, and areas of overlapping flights paths over the steeper slopes had a maximum photo count of 67. The analysis reveals that the model built from the optimized approach delivers a lower standard deviation of 6.78 when compared to the nadir approach that yields a standard deviation of 17.47 (258% the spread of the optimized approach). On average, each point delivers a mean photo count of 23.2 when optimized and a photo count of 28.13 from the nadir approach. These results correspond to the total photo counts of 716 and 912 respectively; this confirms that sloped regions of high interest and flux are better represented by the optimized approach.

Model Observations of Temporal Changes Two Years after the Earthquake
With the 2018 reconnaissance mission, the new techniques involved in collecting the images at desired locations with viewpoints greatly reduce the holes in the models and provide a more complete model of the city with its many slope failures near heavily forested regions (see Figure 12). The textured models, derived from the dense point-clouds, aid in post-reconnaissance visual inspections. The sequential models create "snapshots" in time that are revisited and compared with other "snapshots" to highlight patterns, differences, and continual new discoveries. The 2016 and 2018 models of Pescara del Tronto allow changes in the environment over time to be observed. The optimized flight path for the sUAVs improves observations and makes further study of the area possible. The optimized approach analyzes drop offs, rockslides, and structures at high spatial resolution and is qualitatively clearer than the grid approach. Going to the site to collect data, fly the optimized path, and later create a model economizes time.
The models demonstrate the structural rebuilding process, geotechnical changes, and environmental regrowth that allows changes to be quantified (see major regrowth and environmental changes at Pescara del Tronto in Figure 12). Temporal changes have damaged infrastructure within the city. Reconstruction efforts focus mainly on stabilizing slopes and cleaning the rubble off roads that provide access to bridges. The following images focus on particular points of interest: specifically, the retaining wall and lower landslide (see Figure 13).
As recorded in the models, areas throughout the city are being reclaimed by vegetation. Areas where this is most evident are the upper slopes and lower retaining walls (see Figure 12). New trees and vegetation have begun to reclaim the damaged areas and provide natural protection and stabilization to the failed slopes (see Figure 13). Moreover, while rubble covers the majority of the 2016 models, the 2018 models record vast numbers of trees and vegetation that continue to reclaim most of the damaged areas. Because the soil is loose near the retaining wall, it is steep and thick vegetation has yet to develop. These small details are difficult to study when on site but are clearly recorded in the models-smaller details are clearer from optimized than grid reconstructions. For further research and questions, each location can easily be "revisited" as well as other places of interest that may have been missed while on site through means of the sequential models. The central part of the city against a steep slope and retaining wall noticeably changed from the 2016 to 2018 models (see Figure 13). Over the span of two years, the steeper sloped area of the lower landslide (Section A) did not progress in its failure mode. Rather, in 2018, the slope was overgrown with foliage. The foliage is beneficial to the slope stability because the slope is reinforced by vines, bushes, and young trees. The main section of the slope failure (Section B), begins to be reclaimed by the foliage. However, with the consecutive quakes and further erosion, the slope continues to fail and slide onto the lower road. The third area of note is the large retaining wall failure (Section C). After the retaining wall fell during the October quakes, wind and rain continue to erode the shear sections of the cliff faster than foliage can reclaim and strengthen the slope. At the toe of the failure, foliage begins to grow on the rubble and likely will grow up the slopes in coming years (the optimized results reconstruct vegetation better than grid results). See http://prismweb.groups.et.byu.net/gallery2/ for access to the 3D models themselves [44].

Results of Observed Damage over Time
The value of these models for earthquake studies is greatly enhanced because of their sequential nature. From the three models of Pescara del Tronto, temporal "snapshots" are available between the August and October 2016 earthquakes, just after the October 2016 earthquake, and then 21 months after that same October earthquake (July 2018). This sequential data is rare, and allows for the opportunity to study changes over time in structures, slope failures, and other earthquake induced dangers and phenomena. The results of observed damage are similar to post-disaster situation maps and grading as used by Copernicus Emergency Management, but due to a longer scope, is an ad hoc addition for sequential earthquake recovery analysis at Pescara del Tronto [15]. The two older grid-based models provide needed data, but the third optimized model increases clarity at points of interest. Nadir-based models are used for the 2016 models and the combined nadir/optimized-based model is used for the 2018 model for qualitative comparisons. Due to the prohibitive size of the digital surface model (DSM) point-clouds without decimation, the primary focus of this section is qualitative visual analysis of the sequential damage over time (see Section 4 for quantitative comments. The study captures moments in time that give key visual information about the behavior of infrastructure developments and recovery after earthquakes.

Damage Observations between the August and October 2016 Events and the July 2018 Follow-Up Visit
Observations and data from this section are largely referenced from the GEER 2017 report of these events [1]. Significant damage was sustained in Pescara del Tronto after each of these earthquake events. The August events (MaxPGA = 0.48 g) left nearly half of the structures collapsed. After the October event (MaxPGA = 0.34 g) structures with damage levels of D4 or D5 increased from 65% to 70%.
Between the October event in 2016 and the 2018 July visit, only about 20% of the shown building/structural areas (see Figures 14 and 15) have been completely cleared away with about 9% partially cleared. Given that a span of over 18 months had passed between the earthquake and the return visit, the time lapse shows how slowly disaster areas can recover.  For analysis and detailed inspection of the damage, Pescara del Tronto has been divided up the into seven focus areas (see Figure 16). The priority areas are grouped by similar damage characteristics. The most damaged areas lie within the dense city center; area one studies the structural failure of densely spaced houses. Areas two and three comprise the geotechnical failure of a retaining wall and failed slopes. Areas four and five are slope stability failures of nearby residences. Area six records the exposed pipeline near a gravel pit on the north side of the city. Area seven shows signs of a slope stability failure next to a bridge abutment. Location 1 (see Figure 16) is characterized by talus deposits with dense structural damage and rubble [1]. After the first earthquake in August 2016, this area still had many damaged or compromised structures standing and was scattered with rubble. After the October 2016 earthquake, almost no structures were left standing and the whole area was covered in dense rubble (see Figure 17). The rubble is often dense enough that individual building locations become hard to distinguish. Specific damage to structures is shown in Figure 14 and cleanup is graphically analyzed in Figure 15.
In the 2018 visit, much of the area was left as found following the October 2016 event; however, a portion was cleared out towards the north-west. Data from 2018 shows just under 70% of Pescara del Tronto has yet to be cleared and recovered (see Figure 15). As previously stated, ≈20% of the area has been completely cleared and ≈9% has partially been cleared of debris. The process of clearing out the rubble and improvements for future reconstruction are documented two years since the cataclysmic series of earthquakes in 2016.  Figure 13). The angle of repose of the native talus behind the retaining wall was measured as 53 • . Again, note that the sharp angles and vegetation of the optimized model are cleaner than previous grid-based models.   [3] and the third (c) is from a combined nadir/optimized model generated July 2018. After the first earthquake, the gully wall depressed slightly; however, after the second earthquake in 2016, the 9 m thick and 20 m wide slope catastrophically fell leaving the 12 m scarp at an angle of nearly 52 • , and the residential structure after August 2016 has crumbled into the gully as of November 2016 [3]. As of July 2018, the slope continues to erode and fail, but there have been some repairs to the road and slope. The stairs in the upper right of the shown portion for all three models, are more accurate in the optimized model than in the older models. Towards the northern end of Pescara del Tronto rests a shallow landslide in the desolated cityscape addressed in Stewart et al. [3] and is paraphrased and expounded upon (see Figure 22, the first two  As discussed in Stewart et al. [3] and shown in Figure 23 the pipelines off the haul road and near the gravel pit were affected by the earthquakes and over time ((a) and (b) are referenced in the GEER report). Additionally, Figure 23c shows the same location as of July 2018. From August 2016 to November 2016 about 8 m of the larger pipeline became exposed by weathering; the lower and smaller pipeline became visible in the same time frame. From November 2016 to July 2018 below the haul road, a hole in-between the larger pipeline and the haul road developed. The lower pipeline became more prominent with continued erosion. The pipelines sag and are more prominently exposed with the passage of time; the exposure is especially evident because the optimized approach forms a more detailed reconstruction of the pipes than the grid approach. The haul road has weathered overtime and has not been repaired since the slope of the hill degraded after the earthquake events until July 2018.

Conclusions and Future Work
The sequence of earthquakes caused considerable damage throughout Pescara del Tronto. The widespread progression of damage in focus areas becomes difficult to track and manage. However, by the aid of sUAVs, digital imagery sensors collected data and provided crucial details after the August and October 2016 earthquakes. In addition, the July 2018 post-recovery reconnaissance sUAV mission, two years later, assessed the progression of damage. These images generate dense point-clouds and textured mesh models that are visually and analytically compared to obtain information about the construction of the city, geotechnical background of the city's foundation, and the surrounding countryside.
During the 2018 reconnaissance trip, optimized view and path planning algorithms, increase the completeness of the 3D models by providing needed oblique angles (despite increased vegetation growth). When compared to the other model created from data obtained by a typical flight path developed by GSPro, optimized results (using genetic principles) cover 17% more surface area, better represent objects, contain fewer holes on the slopes and sites of interest, requires 30% fewer photos, has a higher resolution adjusted GSD of 1.5 cm/px, and yields an overall GSD of 2.4 cm/px, and is a case study of the applied algorithms for post-earthquake sUAV remote sensing. Optimized sUAV missions as a resource in post-hazard reconnaissance strengthens the ability to delineate and rate areas of focus affected by earthquakes. The constructed models provide a detailed reference source for the tracking of recovery efforts. Seven areas of interest were identified (see Figure 16). The focus areas depict two main causes of damage: geotechnical slope failures; and structural failure of old masonry buildings. Sequential changes to each area of interest are compared using models from August 2016, October 2016, and July 2018. A specifically promising avenue of future work is to compare complete DSMs of sequential changes with the quantitative tools available from CloudCompare. The progression of damage and recovery are assessed and show many areas being reclaimed by vegetation. While complete event details are difficult to analyze in person, event specifics and recovery are recorded in the sequential models. Further study of each location can easily be "revisited" (as well as other places of interest) and research of focus areas, throughout the models, promotes better understanding of the recovery that takes place after a series of earthquakes devastates a city.