Time-Lapse Landform Monitoring in the Pisciarelli (Campi Flegrei-Italy) Fumarole Field Using UAV Photogrammetry

The utility of new imaging technologies to better understand hazardous geological environments cannot be overstated. The combined use of unmanned aerial vehicles (UAV) and digital photogrammetry (DP) represents a rapidly evolving technique that permits geoscientists to obtain detailed spatial data. This can aid in rapid mapping and analyses of dynamic processes that are modifying contemporary landscapes, particularly through the creation of a time series of digital data to help monitor the geomorphological evolution of volcanic structures. Our study comprises a short-term (in geological terms) monitoring program of the dynamic and diffuse Pisciarelli degassing structure caused by the interplay between intensive rainfall and hydrothermal activity. This area, an unstable fumarole field located several hundred meters east of the Solfatara Crater of the Campi Flegrei caldera (southern Italy), is characterized by consistent soil degassing, fluid emission from ephemeral vents, and hot mud pools. This degassing activity is episodically accompanied by seismic swarms and macroscopic morphology changes such as the appearance of vigorously degassing vents, collapsing landslides, and bubbling mud. In late-2019 and 2020, we performed repeated photogrammetric UAV surveys using the Structure from Motion (SfM) technique. This approach allowed us to create dense 3D point clouds and digital orthophotos spanning one year of surveys. The results highlight the benefits of photogrammetry data using UAV for the accurate remote monitoring and mapping of active volcanoes and craters in harsh and dangerous environments.


Introduction
In this paper, we present the first application of time-lapse unmanned aerial vehicle (UAV) photogrammetry for landform monitoring in the Pisciarelli fumarole field. This active volcanic context, one of the largest and most dangerous in the world, is located inside the Campi Flegrei caldera ( Figure 1) near Naples, Italy. Digital photogrammetry using UAV image acquisition has become an important and widely used technique for investigating and monitoring areas prone to natural hazards. It is exceptionally useful in the geosciences when combined with the Structure from Motion (SfM) algorithm. The main purpose of this paper is to detail the application of well-developed technologies (UAV and SfM) to an under-explored field of geoscience: the landform change monitoring of highly dynamic environments such as volcanic and hydrothermal fields. The integration of time-lapse UAV survey data and time series of available geophysical/geochemical and meteorological parameters helped us to characterize geomorphological changes and investigate their causes. The main advantage of UAV photogrammetry is its ability to reconstruct 3D structures from 2D images; in particular, it permits the evaluation of the shape, position, and size of target elements or areas with great precision and in a shorter time frame than field investigations. This can be attributed to the improved technical performance of drones [1,2] and the development of upgraded artificial vision algorithms [3], or information technologies which significantly accelerate the processing times and the quality of 3D surface reconstructions [4,5]. The Structure from Motion (SfM) algorithm used in photogrammetry opens new application possibilities in the geosciences [6][7][8][9] because it allows for the creation of high-resolution digital elevation models (DEMs). It can facilitate the identification of centimetric-scale features of natural and artificial surfaces quickly at a low cost and in safer conditions for researchers. The number of applications of SfMand UAV-based photogrammetry is growing rapidly in geomorphological research. This methodology is particularly suitable for repeated topographic surveys in remote or poorly accessible areas due to the combined accuracy and reproducibility reached by topographic data based on UAV imagery [10,11]. Some of the most striking uses of the technology include geomorphological mapping [12] and studies of landslides [13], erosion [14,15], and sedimentation [16].
In active volcanic contexts, it is essential to monitor the dynamic topography frequently. Detailed DEM analyses and multitemporal comparisons can aid in this endeavor. Photogrammetric techniques, too, can be successfully applied to study volcanic and geothermal hazards. The applications of these techniques range from change analyses in the morphology and structures of lava domes [17][18][19][20][21] and dykes [22] to volcanic landform monitoring [23], ancient lava flow field surveys [24], volcano-tectonic surveys [25], geothermal research [26,27], and volcanic gas and plumes [28][29][30][31]. These studies can additionally include the use of optical and thermal imaging capabilities [32] to provide frequent and detailed detection and monitoring while simultaneously permitting geoscientists to detect changes in volcanic activity without putting themselves in hazardous conditions. The time-lapse technique-which consists of performing identical photogrammetric surveys multiple times in the same location over defined time intervals-is highly useful for identifying topographic changes and yet has been rarely used in active volcanic areas for multitemporal monitoring [33].
In the context of Italian volcanism, this approach has only been used to monitor active lava flows in Mt. Etna (Sicily) [8,34] and coastal landslide evolution and geomorphic changes along tuffaceous coastal cliffs located around the Campi Flegrei caldera (Naples) [35,36]. The Pisciarelli fumarole field, a locus of both ancient and ongoing morphological changes, is localized near the Solfatara crater at the border between the towns of Pozzuoli and Naples. Historically, a vast amount of geochemical, geophysical, and topographic data is collected-both in discrete and continuous form-from this site, and many of its fumarolic vents have been subject to routine gas sampling and real-time geochemical measurements. In recent years, however, the Pisciarelli fumarole field has become less accessible due to the increasingly frequent openings of small fumarolic vents and the presence of boiling mud, hence the importance of remote sensing. In order to monitor and quantify the fast changes in the morphology of this area, UAV time-lapse surveys and SfM-based photogrammetric processing were performed starting in April 2019. These allowed us to highlight variations in the volume of the mud pools, acidic hot springs, and other topographic changes.
The geomorphological dynamics of this area have been characterized by both vertical and horizontal ground movements at gradients that range from centimeters to meters per year [49,50]. Two recent periods of unrest (in 1969-1972 and 1982-1984) caused a cumulative maximum uplift of over 3.5 m in the central sector of the caldera. A new episode of intense ground deformation started in 2011 and is characterized by a current cumulative maximum ground deformation of about 0.70 m [51].
The Pisciarelli area falls within the central sector of the Campi Flegrei caldera between the Solfatara, Astroni, and Agnano vents ( Figure 2). The Astroni volcanic edifice was built through several eruptions at approximately 3.8 kyr [47]. The large Agnano-Monte Spina eruption occurred at about 4.1 kyr and was followed by the volcano-tectonic collapse of the Agnano plain [52]. The erupted magmas varied in composition from trachybasalt to latite, trachyte, and alkali-trachyte.  [37]. The used coordinate system is WGS84/UTM zone 33N.
Solfatara is a maar-diatreme structure located at the caldera center and crossed by NW-SE and NE-SW faults [53]. It was formed 4.2 kyr ago by phreatic and phreatomagmatic eruptions and overlays several tephra deposits and the Mount Olibano lava dome [53,54]. The Solfatara hydrothermal system lies atop a hydrothermal plume created by large quantities of hot fluids (mainly H 2 O and CO 2 ) derived from a magmatic body located at a depth of about 8 km [55,56]. After mixing with meteoric components in a hydrothermal reservoir located at a 3-4 km depth, these fluids are released through diffuse and direct degassing at the surface of the Solfatara crater [57]. The main surface hydrothermal features are represented by the Fangaia mud pool, the Bocca Grande and Bocca Nuova fumaroles within the Solfatara crater, and the fumarolic area at Pisciarelli located east of the Solfatara crater [56] along the intensely urbanized western margin of the Agnano crater ( Figure 2).
The Pisciarelli fumaroles are located at about 75 m a.s.l. at the base of the slope declining from Solfatara vent volcanic rim (175 m a.s.l.). This slope is characterized by valleys with ephemeral streams separated by hilly divides (Figure 3). A hilly divide is located at 170 m a.s.l. and delimits the west and south sides of the fumarole field. The eastfacing slope of the hill can be impacted by low-magnitude (1 to 10 m 3 in volume) rockfalls, deposits of which can reach the underlying pools ( Figure 3). Multiple small fumaroles, bubbling mud pools, and fractures with diffuse soil degassing and fluid emission surround the main fumarolic vent, which is locally named "Soffione" (gas temperature~114 • C). This vent is comprised of a semi-circular shaped bubbling mud pool (water temperaturẽ 90 • C) [58]. It was conjecturally formed in 2009, accompanied by a small explosion [59]. The bulk CO 2 flux from the Pisciarelli and Solfatara sites reaches 50,000 g/m 2 d over an area of ∼1 km 2 [60]. Fumarole discharges are mainly composed of H 2 O (>70%), followed by CO 2 and H 2 S, and have surface temperatures between 40 and 60 • C. [61][62][63]. This degassing activity is episodically accompanied by local seismic activity at a shallow depth (1-2 km) and by some macroscopic changes in ground topography. The degassing activity is mainly controlled by both meteoric water and shallow aquifers [64]. The Pisciarelli fumarolic field evidenced an increase in hydrothermal activity from 2005 to 2013 [59,65]; this is marked by a nearly linear trend in peak temperatures, from about 97 up to 112 • C, and by increased local seismic activity [66][67][68]. Several ground changes have occurred in recent years [69]. An episode of mud emission, prompted by a weak phreatic explosion, occurred in  [69]. On 1 December 2018, a significant enlargement of the emission area occurred [70]. The activity has remained very high in the last few years and has been monitored through our surveys during the 2019-2020 period.

Materials and Methods
Our fieldwork involved acquisition, reconstruction, and analysis sessions. In order to obtain multitemporal 3D point clouds and high-resolution orthophotographs, we performed repeated UAV surveys over the active Pisciarelli fumarolic field from October 2019 to November 2020. In total, we made ten surveys with interruptions due to the COVID-19 pandemic (March-April 2020) and drone maintenance (July and August 2020).
The UAV imagery was acquired with a DJI Phantom 3 Standard imaging platform with a Sony EXMOR 1/2.3" with effective pixels/photosites of 12.4 megapixels. The sensor is 6.16 mm wide and 4.62 mm high; the maximum created still image size is 4000 pixels wide by 3000 tall. The field of view is 94 • with a 3.61 mm focal length, and the platform has an internal Global Positioning System (GPS) positioning system. The hover accuracy range is ±0.5 m (vertical) and ±1.5 m (horizontal).
The flight path was pre-configured by the flight-controlling commercial software Pix4Dmapper ® [71] and was performed at an altitude of 50 m above the ground surface. An average of 100 photos with an 80% overlap and sidelap was acquired in each flight [72,73]. The flight height of 50 m above the ground was selected by taking into consideration different requirements, such as narrowness, the morphological instability, and the harsh conditions (highly corrosive atmospheric environment conditions with high temperature and humidity) of the investigated area, and the presence of numerous vertical artifacts such as poles or nets. These conditions were considered to select the lowest possible Ground Sampling Distance (GSD) and consequently obtain very detailed images. GSD, generally expressed in cm/pixel, represents the distance between the centers of two consecutive pixels in the image on the ground [74]. For a given flying height H (m), the GSD is given by: where P is the pixel size of the sensor (micro) and f is the focal length (in our case, 3.61 mm).
In this study, the obtained GSD is about 2.13 cm/pixel. The use of orthophotos with a high degree of definition allows the reconstruction of a more detailed three-dimensional model which is able to describe the ground surface in a much more detailed way. The total area covered by the reconstruction of the ground surface was 0.9 km 2 . Two different waypoints were used in each survey with the aim to obtain the best 3D image configuration. The flight time to complete each waypoint by drone was roughly 6 min. Figure 4 shows an example of a rectangular waypoint grid acquisition at a flight height of 50 m with the camera lens in a nadir orientation and an elliptical waypoint acquisition grid at a flight height of 50 m with oblique camera for the best definition of the vertical surface. We repeated the same waypoint for each survey. The geometric accuracy of 3D models and derived DEMs also depends on the acquisition of Ground Control Points (GCPs). These GCPs are essential for geo-referencing the images and their homogenous distribution on both sides of the study area, which will influence the accuracy of the obtained point cloud [75]. In our study, nine fixed GCPs were distributed unevenly over the area for geo-referencing the reconstructed models. The surveyed areas are mostly covered by natural features that do not provide well-defined points to be used as GCPs, such as rocks and vegetation. Therefore, it was decided to put artificial points on the surface. These were mainly plastic bands with crosses in the middle and in contrasting colors. The position of these markers was measured in the WGS84/UTM (Universal Transverse Mercator)-Zone: 33N-EPSG (European Petroleum Survey Group): 25,833, with a static survey using Trimble Geo 7X GNSS (Global Navigation Satellite System), which is characterized by a horizontal accuracy of 3 mm + 0.5 ppm RMS and a vertical accuracy of 3.5 mm + 0.5 ppm RMS.
The Agisoft ® Metashape software package version 1.6.4 [76] was used to process the photographs using SfM photogrammetric methods. The applied SfM image processing consists of extracting tie points from the photographs at a lower resolution, performing GNSSassisted aerial triangulation and bundle adjustment, and generating three-dimensional point clouds and a 3D model ( Figure 5). The SfM-based photogrammetric processing is designed to derive 2D orthomosaic and 3D DEMs data from a set of overlapping images acquired in several directions in order to cover the area of interest [77,78]. The parameters of the flight, imagery, and terrain reconstruction used for each survey are summarized in Table 1. The comprehensive processing required to obtain high-quality results during work sessions required 1 to 3 h with a PC-based workstation featuring a Core i9 9900 k processor, 32 GB DDR4 RAM, and Nvidia GeForce RTX 2080ti 11gb. The image orientation was performed with normal full automation setup defaults by choosing the high-resolution option to have the maximum average number of points per photo. The obtained final value of root mean square error (RMSE), expressed in pixel units, was below 1 pixel RMSE for all the performed projections. It is important to stress that the high resolution of the orthophotos and the accuracy of the entire post photogrammetric processing allowed us to obtain accurate measurements of the fallen rock sizes and pool geometry within an accuracy range of ±2 cm. In particular, their volumes, areas, and perimeters are measured by means of an automatic tool of Agisoft Metashape software.
The software Cloud Compare [79], part of an open source project, was used to align the point cloud to correct for the incorrect georeferencing of the GCPs caused by the loss of some markers. Through the alignment toolbox, we obtained a cloud-cloud distance computation for each pairing of point clouds. The alignment procedure involves selecting the point clouds that will remain static and those which will move. After having chosen at least three constraint points (Figure 6), the Cloud Compare software makes a calculation based on a standard deviation of the displacement. After this procedure, all the raster outputs were aligned with each other. In each survey, we used nine fixed GCPs as control points and three as check points (CPs) for georeferencing the model to calculate the exterior orientation parameters of the point clouds. During the photogrammetric process, we performed manual GCPs positioning on images and then alignment optimization using the selected GCPs. Although the process of identifying GCPs is manual, it is relatively fast because their positions are predicted in all images and the user only needs to correct them. The bundle adjustment was repeated in a fixed reference system with camera self-calibration. The time-lapse analysis of DEMs and orthophotos allowed us to detect changes in the geomorphology and topographic features in the study area. Multitemporal quantitative geomorphic changes in the fumarolic area were identified by comparing the SfM-derived point clouds in month pairs. In this way, we obtained time series of morphometric parameters of the main fumarolic vent. To aid in the interpretation of the detected geomorphic changes, we compared these with the datasets of several parameters. Meteorological data were derived from the monthly rainfall measured at the Pozzuoli rain gauge, which is part of the real-time hydrological monitoring network of the Campania Regional Agency of Civil Protection [80]. The rain gauge is located along the northern sector of the Solfatara crater rim, about 700 m away from the Pisciarelli area. Geochemical data relating to the CO 2 soil fluxes were recorded since 2005 via the FLXOV3, an automatic station located in the Pisciarelli area. The CO 2 flux is measured every 2 h using the accumulation chamber method [81]. Seismological data recorded by the Istituto Nazionale di Geofisica e Vulcanologia-Osservatorio Vesuviano (INGV-OV) seismic network at Campi Flegrei [82] were also used for correlation with the morphological data.

Results
During surveys beginning in October 2019, we performed repeated UAV surveys over the active Pisciarelli fumarolic field that allowed us to assess morphological changes. Two main types of changes occurred: rockfall events and shape variation in the main hot pool. A small rockfall occurred between November and December 2019 on the southwest slope of the fumarole field, as shown in Figure 7. This episode could have been caused by severe rainfall or (more likely) could be a result of the seismic swarms that occurred on 6 December 2019 with a main 3.1 M seismic event localized in the Pisciarelli area at 2.3 km in depth [82]. Another rockfall episode occurred between December 2019 and May 2020. In Figure 8, we show the presence of rock debris at a depth of 30-40 cm thick on the edge of the bubbling mud pool. It is not possible to exclude the possibility that a low-energy phreatic events proximal to the main fumarole vent may have caused the fall of rocky boulders from the slope surrounding the pool.  An additional aspect that we determined from the time-lapse images is the changing shape of the main bubbling mud pool over time ( Table 2). In autumn-and, in general, during more intense rainfall periods-water accumulated in the shallowest topography and the bubbling mud pool appeared larger than in dry months. One of the most relevant changes was observed in the month of September 2020, when the water level of the bubbling mud pool dropped dramatically and the pool appeared completely dry. It is possible that this episode was associated with an unusual seasonal effect, such as very low rainfall (i.e., 3-10 mm during the previous six months). In Figure 9, the graphic differences in the bubbling mud pool in December 2019 versus September 2020 are shown, indicative of the wide range of morphological shapes the pool took during that year. Due to the acquisition of ortho-photos, it was possible to measure the size of the two areas ( Figure 10) located in the center of the area previously occupied with water and bubbling mud.   Multitemporal quantitative geomorphic changes in the fumarolic area were identified by comparing the SfM-derived point clouds in the following pairs: (2019) October-November, November-December, December-January (2020), January-February, February-May, May-June, June-September, and September-October ( Figure 11A-I). Figure 11. Time-lapse of the relative heights (i.e., differences in height) obtained by comparing the SfM-derived point clouds in pairs. The "Soffione" main vent is indicated with a black small triangular symbol. The blue box indicates the area characterized by strong changes.
Overall, the image sequences confirm multiple episodes of rockfall and the accumulation of mud, particularly on the western side of the fumarole field. In Figure 11A, we quantify a 10-12 cm disparity visible on the western slope. This change is related to a maximum volume of 0.010 to 0.0025 m 3 and is associated with the presence of rocks on the ground. In Figure 11C,D,F, the differences in height (13-20 cm) affecting the north-south west side of Pisciarelli are probably due to the accumulation of debris as a consequence of intense rainfall, which is visible in Figure 11H. In Figure 11G, the void left by the drying up of the bubbling mud pool and the accumulation of mud are clearly evident.
In other cases ( Figure 11E,H), changes in volume are concentrated around the main fumarole (black triangular symbol) probably as a result of minor explosion events during intense hydrothermal activity. Through the analysis of the topographic changes detected with the time-lapse surveys compared with the obtained high-resolution orthoimages and 3D model, we identified two main landslide areas impacted by rockfall episodes at the western side of the Pisciarelli fumarole field. The limits of these two areas-respectively, 202.3 m 2 and 119.7 m 2 wide in November 2020-were automatically drawn with the segmentation tool of Agisoft Metashape (Figure 12). The topographic differences between the October 2019 and November 2020 surveys allowed us to evaluate the volumes of fallen materials, estimated at 78 and 22 m 3 respectively ( Figure 12). The draining of the bubbling mud pool, which occurred between July and August 2020, is clearly visible in Figure 13. The June and September 2020 point clouds show a void in the central part of the fumarole field and an accumulation of muddy material around the pool (10-12 cm high and a total volume of~3 m 3 ). Some days later, a spatter of ash deposit was observed close the eastern side of the bubbling mud pool.

Discussion and Conclusions
This is the first study using UAV-based time-lapse photogrammetry to analyze the morphological changes in Neapolitan volcanoes, specifically in the dramatically dynamic and challenging study area of the Campi Flegrei caldera and the Pisciarelli fumarole field. Analyses of the orthoimages and DEM-based volume differences allowed us to identify the morphological and topographic changes that occurred during a year of repeated, systematic survey. This study confirms the presence of instabilities on the west slope of the fumarole field, as demonstrated by the occurrence of rockfall episodes during intense rainfall or seismic swarms. This is crucial in terms of in situ monitoring, particularly for the safety of geoscientists.
Time-lapse photogrammetric data by UAV surveys also allowed us to follow the temporal evolution of the geometry of the main bubbling mud pool of Pisciarelli. The behavior of the bubbling mud pool appears in part to be correlated with rainfall rates. During warmer and drier months, we observed a decrease in the water level (and a minor size of the pool) and vice versa. This suggests that the geometry of the bubbling mud pool may be influenced by seasonal effects that are probably driven by the influx of rainwater into the system.
The macroscopic changes in the Pisciarelli area, including the overall expansion of the bubbling mud pool, can also be attributed to local hydrothermal volcanic activity. An increase in geochemical and geophysical signals, not observed prior to the onset of the unrest in 2005 and with an overall uplift of~57 cm and~448 seismic events, was recorded in 2018-2019 [58]. Additionally, the high pressure of fluids (CO 2 and water) mobilized deeper underground is capable of altering the host rocks mechanically and/or chemically. To consider such connections, analyses of the morphology of the bubbling mud pool and other geophysical data were attempted. The parameters considered here are the CO 2 fluxes, the rainfall rates, and the local seismicity. The CO 2 fluxes were recorded in the Pisciarelli main vent during the period 2007-2020 (starting in 2012, the Campi Flegrei is at a "yellow" alert level in accordance with National Civil Protection, and for this reason the rate sampling of the CO 2 soil flux is now weekly), in correspondence with the FLXOV3 measuring station. The monthly rainfall rates, relative to those measured during 2019 and 2020 at the meteorological station of Pozzuoli [79], are also reported in Figure 15.
The comparison of the rainfall rate and CO 2 soil flux, reported in Figure 15, shows a correlation with the geometry of the bubbling mud pool. A strong seasonal effect seems to influence the behavior of the CO 2 fluxes and the expansion of the bubbling mud pool. Such a correlation could be due to a buffering and/or cumulation of the gas flowing from the geothermal system to the atmosphere in those zones. In fact, the water accumulated in the soil during the wet season could explain the low values of CO 2 fluxes during the hot season. This seasonal effect is also confirmed by the geophysical study, based on Electrical Resistiviy Tomography (ERT) time-lapse surveys applied in the same study area. Troiano and Di Giuseppe [33] demostrate the usefulness of time-lapse ERT to discern the areas dominated by seasonal effects from the ones more influenced by other, most significant endogenous contributions, linked to gas injected into the system from deep locations which cumulates in the subsoil.
The comparison between local earthquakes with the perimeter and area of the bubbling mud pool, reported in Figure 16, suggests a link between the parameters. The contribution of the seismic events seems to have an effect on the expansion of the bubbling mud pool, particularly when a seismic swarm occurs. At least in principle, the contribution of earthquakes to the observed changes in the geometry of the bubbling mud pool, such as fractures or rock alterations as a consequence of the very low depth seismic events (0.5-1 Km), cannot be excluded. Unfortunately, our data collection was interrupted by the COVID pandemic from mid-November 2020, and hence we lack data to confirm a the further expansion of the bubbling mud pool, as happened in 2019. The results of the time-lapse analyses suggest that the observed rock mass deposits were likely caused by small, local seismic swarms in conjunction with rainfall-induced surface erosion. In this scenario, the monitoring technique described here provides useful data for hazard assessment in addition to a valid opportunity to monitor the harsh, inaccessible environment that an instable fumarole field represents. The capability to remotely monitor topographic changes, landslides, and surface mud pools is extremely important in an active volcanic area, especially during periods of unrest, which could be preludes to eruption and which are also time periods in which geoscientists are unable to collect data at close proximity. This technique, clearly, can also be of great use during and after dramatically transformative eruptions to rapidly update the DEMs of volcanic areas. This technique could additionally, however, be used for rapid, preliminary inspections during and after myriad other disasters, including earthquakes, landslides, floods, and fires.
In summary, our research confirms the usefulness of unmanned aerial vehicles for image collection and of SfM photogrammetry process for the long-/short-term investigation of single geo-hazards. These techniques can greatly reduce the bodily risks associated with fieldwork and provide exceptionally valuable, accurate, and high-resolution data to support volcanic hazard monitoring in a densely inhabited volcanic area.

Data Availability Statement:
The data presented in this study are available in this article.