UAVIMALS: THE "OPEN" REMOTE SENSING SYSTEM FOR SURFACE ARCHAEOLOGICAL INVESTIGATIONS

: Today, there is an increasing use of airborne sensors in archaeology, especially to investigate the surface of more or less vast territories quickly and accurately. In particular, airborne laser scanning technologies from small remotely piloted aircraft are rapidly developing towards increasingly high-performance solutions for the investigation of archaeological evidence hidden by vegetation or more or less substantial soil deposits. The proposed contribution intends to insert itself within this field of archaeological research by presenting "UAVIMALS" (Unmanned Aerial Vehicle Integrated with Micro Airborne Laser Scanner), an aerial remote sensing system of "soil marks", designed for surface archaeological investigations and the result of an Early Career Grant from the National Geographic Society. The system, consisting of a customised drone based on an open architecture and software for vehicle control and data processing, integrates a solid-state laser sensor, commonly engineered for obstacle avoidance, but here exploited to process accurate DTM (Digital Terrain Model) of small land surfaces with a significant reduction in cost and acquisition time. The system, whose engineering was contributed by the BioRobotics Institute of the S. Anna University of Pisa, was tested within the archaeological context of Leopoli - Cencelle (Tarquinia, Italy). A mediaeval city that has been researched for about 25 years by the Chair of Christian and Medieval Archaeology at the 'Sapienza' University of Rome. Experimentation missions carried out on the site, which is still only partially investigated, have been successful in bringing to light some urban areas that had not yet been investigated.


INTRODUCTION
The UAVIMALS system 1 is the result of interdisciplinary research between archaeology and biorobotics, carried out between the Sapienza University of Rome and the Scuola Superiore Sant'Anna of Pisa, which led to the creation of a small size aerial laser scanner prototype, useful for light archaeological investigations. The project was funded by a National Geographic Society grant awarded in September 2018, with an Early Career  by Federica Vacatello (Principal Investigator and Post Doc Researcher in Archaeology and Post-Classical Antiquities at Sapienza University of Rome) who coordinated an interdisciplinary research team made up of archaeologists, engineers and biorobotics technicians. The main goal of the project was to create a simple and fully Open Source system focus on cushioning the cost and size of instruments used in surface archaeological investigations. By experimentation with an engineered LiDAR sensor for autonomous vehicle guidance, it was possible to create a specific aerial system for the detection of archaeological traces from 'micro-relief'. This type of anomaly is usually detectable only under given light conditions or through the use of specific aerial instrumentation (Piccarreta -Ceraudo 2000). The ambition of the UAVIMALS project was therefore not to create a lowcost, lower-performance airborne LIDAR than those already on the market, but rather to create an instrument that is easy to transport, less expensive and equally accurate. We believe that the solution can provide an advancement of research in the field of airborne laser scanner technology. The acquisition of three-dimensional images with very high morphometric resolution has proved to be a practice of fundamental importance for the study of various contexts on our planet, but in the field of archaeology, in particular, remote sensing by drone represents a practice of extreme * 1 Corresponding author importance for the investigation of ancient structures, sometimes still unexplored, that cannot otherwise be investigated by other means, such as excavation and reconnaissance activities, due to uncomfortable geomorphological conditions, places of difficult access and traces invisible to the human eye especially in particular climatic and vegetative conditions (Pfeifer -Gorte -Elberink 2004;Optiz -Cowley 2013;Optiz -Herrmman 2018). Despite this, most commercial instruments are still prohibitively expensive for archaeological research, as well as having unfavourable dimensions for meeting transport needs in inaccessible locations without means of travel. Experimental tests in the archaeological context of the mediaeval city of Leopoli -Cencelle (Tarquinia, Italy) were successful in identifying a number of as yet unknown urban structures not yet intercepted by excavation surveys. In association with the hardware system, a software application is also being developed which is useful not only for controlling the vehicle during flight, but also for monitoring the data acquired with a view to initial graphical processing. At this moment, it is only possible to process the point clouds of LiDAR data by means of dedicated software CloudCompare (CloudCompare, Development Team, 2023); 3DF Zephyr (3DF Zephyr, Development Team, 2023); QGIS (Development Team, 2022) etc. which, not being connected with the drone, do not allow a real time visualisation of what is seen by the sensor, preventing a preliminary monitoring of any archaeological presence hidden in the flyover area. The DTM, meshes and point clouds obtained from the sensor can later be uploaded into geospatial software such as QGIS, thus enabling spatial, territorial and geomorphological analyses of the data acquired through the use of specific tools. If for other application contexts this activity may be useless, in the field of archaeology, a system such as the one developed may represent a concrete possibility of extending archaeological investigations, which would thus be speeded up by an observation tool as well as facilitated by a cost widely accessible to university research funds.

MATERIALS AND METHODS
It is now fairly well known in the literature that the investigation of the human past requires an interdisciplinary approach in which remote sensing data represent an irreplaceable investigative tool in the search for traces of the past, sometimes explaining the signs of the present. The anomalies, the colour, the patterns visible from a remote sensing scene are all the result of the interaction between natural phenomena and human activities whose ultimate outcome is the landscape (Küçükdemirci et al., 2021). In this regard, it is now widely accepted that certain types of buried archaeological deposits can be identified due to their inherent ability to produce different proxy indicators visible on the ground, even though physical and micro topographical changes that can be intercepted from an aerial view (Crawford 1929;Dassie 1978;Wilson 1982;Optiz and Cowley 2013). In particular, the so-called 'soil marks' (Fig. 1) represent one of the nodal aspects of the 'UAVIMALS' project, which envisaged the realisation of an aerial instrument designed precisely for the identification of the latter in the different territorial contexts with the presence of sparse and low-trunk vegetation. Soil marks are nothing more than micro-relief anomalies caused by archaeological deposits still buried under layers of earth of moderate thickness (Masini -Soldovieri 2017). Light archaeological investigations, today representing a large part of archaeological research, are increasingly being conducted with the use of high-performance, airborne instruments which, however, while yielding a high degree of metric accuracy, are still too limited in their use due to high costs, inconvenient dimensions for inaccessible locations and the limitations of flight licences. Moreover, all aerial scanning systems currently on the market are designed only for data acquisition and not for data processing, a possibility that could lead to a significant advancement in research, especially in the archaeological field, as well as a reduction in investigation times by allowing an initial targeted analysis of any 'soil mark' already in the flight phase through the use of specific software. A workflow capable of automating the elaboration of point clouds and the construction of meshes depicting all surface characteristics of the area being flown over.

System engineering
The 'UAVIMALS' system was engineered by the BioRobotics Laboratory of the Scuola Superiore S. Anna in Pisa, which specialises in the design and development of customised drones capable of flying autonomously to complete specific tasks through the use of Artificial Intelligence algorithms. For the UAVIMALS project, the Institute integrated a LIDAR sensor on one of its best drones (Cascarano et al.,2021) and developed a code to process the acquired data to calculate point clouds and reconstruct the 3D soil's geometry. The drone is designed as a highly customisable modular structure in both hardware and software. Hardware-wise, the flying platform is a quadhelicopter Tarot 650 Sport mainly equipped with an on-board computer (NVIDIA Jetson TX2) both connected to a solidstate LIDAR sensor for scanning sites of interest and a Flight Control Board (FCB -Pixhawk 1) including its own flight sensors and a precision laser altimeter. (Fig.8) The on-board computer integrates the navigation system with the scanning sensor, collecting both flight and survey data and synchronising them directly on board in real time during the flight. From the software point of view, the modular structure is fully adaptable to different hardware configurations and makes it possible to separate the navigation task, assigned to the FCB, from the automatic mission management and data acquisition, provided by two dedicated Python codes running on the on-board computer. Finally, the flight mission, in terms of GPS waypoints over the site, is sent to the drone from an integrated open-source Ground Control Station (Mission Planner) via an easy-to-use operator interface running on a laptop computer. The data collected during the flights consist of the position and attitude of the drone and the relative distances acquired by the scanning sensor. First, a Python GUI (Graphical User Interface) allows the user to prefilter the dataset and select the data corresponding to a specific region of interest within the digitised area. The resulting data are the input for a MATLAB® (MATLAB, 2018) code that reconstructs and displays a raw point cloud.
The cloud was subsequently exported to CloudCompare and here cleaned of points corresponding with tall vegetation, most of them obscuring the morphological data of the underlying terrain. (Fig.9).

Data acquisition
Adapting a sensor not designed for remote sensing but rather for obstacle detection and automatic vehicle braking (ADAS) was a challenge that involved finding specific software solutions. Specifically, the chosen sensor returns only distance measurements along 16 cones, as you can see in Figure 2 A, the 16 cones are aligned on a straight-line Figure  2B and thus producing a one-dimensional GF ground footprint, the one-dimensional nature of the GF prevented the direct use of previous and subsequent acquisitions for software corrections e.g. Iterative closest point algorithms, which require a common area of 2 acquisitions to correct the point cloud.
The amplitude A_gf of GF i.e., the distance between bin 1 and bin 16 measured in a plane, is directly proportional to the distance of the measured surface and can be calculated as: The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLVIII-4/W7-2023 FOSS4G (Free and Open Source Software for Geospatial) 2023 -Academic Track, 26 June-2 July 2023, Prizren, Kosovo The A_gf is extremely important because it determines the number of passes required to cover a given area and thus the duration of the flight for point cloud creation and its resolution on the ground. So first variable we had to resolve was the height of the flight H to be kept during the acquisition. Suppose for example that we want to fly at a speed = 2 ⁄ and we want to examine a field of 100m*50m ( F_l, F_h ). Flying with a height H=5m the A_gf is circa 1.67m with a distance between the center of bins of approximately 10cm along the y-axis in the drone's coordinate system (Fig.4). An A_gf = 1.67 involves a number of passes P above the field calculated as = ⌈ _ℎ _ ⁄ ⌉ = 30. The drone will have to approximately travel a distance = ⌈ * _ ⌉ = 3000m which, at speed S corresponds to a flight time FT=25 minutes. By repeating the same calculations with a H = 6m we will obtain instead: A_gf =2.01 center bin distance 12.6cm → P = 25 → D=2500m → FT = 20 minute and 50 seconds. So, a single meter difference leads to a much shorter flight time then after numerous tests carried out at the training camp, we decided that an H between 10m and 12m was the best compromise between flight duration and resolution on the ground.
The distinctive feature of the sensor chosen, compared to others on the market of the same cost, which influenced its choice was the ability to return up to 3 measurements in the same bin associating them with a reliability index based on the quality of the measurement performed. As you can see from Figure 3 which schematizes an example of acquisition, in the single bin represented with an orange cone, 3 different distances are acquired: the distance from the tree in blue, the distance from the tufts of grass in green and the distance from the soil. By acting on the software parameters that the sensor provides, we were able to make an initial screening of the measurements that did not reach a certain quality. And subsequently, through the cloud points creation program, we went to choose only the measurements concerning the ground, limiting the overlying vegetation when possible, simply choosing the largest among the 3 with the goodness above a certain threshold. This element has allowed the use of the instrument for measuring the distance from the ground, thus eliminating the measurements received from the overlying vegetation (Fig.3).
Therefore, as explained, the sensor returns a set of lengths, from 16 to 48. In order to transform these lengths into points in a three-dimensional space it is necessary to have fixed external points so that they can be used as a reference and a common coordinate system CCS. The first fixed point for each flight is the take-off and landing point of the drone (ground origin HOME CCS→ HCS home coordinate system). The others are instead the set of points which form the perpendicular plane to the straight line passing from the HOME and its zenith, distant from HOME the chosen flight altitude. Figure 3. Example of data acquired within a single bin, the orange cone is the area scanned by the bin under examination, in (blue) the distance between the drone and the tree, in (green) the distance from the tufts of grass and in (red) the distance from the ground.
On this plane it was necessary to control the movement of the drone in order to minimize its vertical deviation.
We have therefore specialized a state of our state machine which is responsible for the various phases of flight: take-off, The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLVIII-4/W7-2023 FOSS4G (Free and Open Source Software for Geospatial) 2023 -Academic Track, 26 June-2 July 2023, Prizren, Kosovo landing, reaching waypoints, safety maneuvers in case of emergency, etc.
The status for reaching waypoints thanks to the use of information received from the GNSS system, from the IMUs, from the laser altimeter and from the sensor, analyzes the deviation from this plane and progressively corrects the altitude. For the HCS we instead opted for a metric coordinate system ENU(https://en.wikipedia.org/wiki/Axes_conventions#Groun d_reference_frames:_ENU_and_NED) which has the origin of the coordinate system on the starting point of the vehicle, both because this system is available as the drone's coordinate system via mavros, and because works in a threedimensional space metric facilitated testing during the early stages of development.
Going from a set of lengths to a georeferenced point cloud in HCS needs some trigonometry and some multiplications of inverted roto translation matrices.
For each acquired length we can generate a point in the 3dimensional space having as origin the sensor itself using the formulas schematized in figure 6: = 0 = * sin = − * cos l : misured lengh, α : angle on the y,z plane, enclosed between the z axis and the straight line passing through the center of the bin cone which generated the length itself. This transformation alone does not return an overall result since each measurement has a different origin O_sensor_1, O_sensor_2, ..., O_sensor_n as the sensor moves in space during the acquisition. In order to merge these spaces into one common to all(HCS) we need the information generated by the FCU of the drone so we need to move each space O_sensor{k} to a new space having as its origin the center of the FCU O_drone{k}. Since the sensor is solidly mounted to the body of the drone, this transformation is calculated as a translation of the origins of x_diff, y_diff, z_diff given by the position of the sensor and the FCU. Finally each point of an O_drone{k} space can be roto translated in HCS using the information of latitude, longitude, altitude, roll angle ɸ, pitch angle θ and yaw angle φ acquired simultaneously with the length l by the FCU. The transformations just explained are schematized in figure 6. The need to also consider the ɸ, θ and psi angles in the last change of coordinates is due to the fact that the sensor is solidly mounted to the drone body. The desire to contain The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLVIII-4/W7-2023 FOSS4G (Free and Open Source Software for Geospatial) 2023 -Academic Track, 26 June-2 July 2023, Prizren, Kosovo costs and weight of the sensor prevented us from engineering it with a gimbal. The complete system would have required a much larger drone frame which would have limited its portability in archaeological areas. Having the sensor fixed has therefore introduced all those errors due to a flight not parallel to the acquisition plane. Figure 7 shows, for example, a roll error due to a crosswind. In fact, in both examples the acquired measurements are identical L1=L3 and L2=L4, but in A the measurements are different L1 > L2 because the ground under the drone is inclined. In B instead L3>L4 because the drone itself is inclined while the ground is flat. The inclination of the drone was therefore considered during the transformation of the Odrone_{k} spaces into HCS.

Figure 7.
Simulation of two different measurement scenarios. The measurements made in situations A and B are the same L1=L3 and L2=L4. However, while in image A L1> L2 depends on the inclined course of the ground flown over, in B the ground is perfectly flat and the difference between L3 and L4 is due to the roll of the drone produced by the wind.

Data analysis
The system, capable of acquiring a cloud of points with centimetric precision, made it possible to extract a DTM reproducing such a detailed morphology of the terrain of the areas overflown that an accurate slope analysis was possible using QGIS open source software (QGIS Development Team, 2022). The cloud was first rasterized using CloudCompare's 'rasterize' tool in order to create a digital terrain model with 20 cm ground resolution cells. The reference parameter set in the " Step Grid" may change in function of the final objective for which the DTM was created, however, in the case of archaeological analyses, the lower the ground resolution of the cells, the more detailed will be the mapping of the possible archaeological indicators present in the area under survey (Masini et al., 2011, pp. 263-290).
In fact, most archaeological deposits, both earthy and builtup, have thicknesses greater than 20 cm, so the decision to create such a DTM depended on the need to be able to visualise as many variations as possible. In addition, during the rasterization phase, the density of the starting point cloud allowed the empty cells to be filled in by a fast interpolation process.
The DTM obtained was loaded into the QGIS 3.30 software to modify the visualisation using different colour scales based on the Z value of each cell. By varying the different scales and degrees of brightness and contrast, all those more superficial archaeological indicators present within the overflown area were first highlighted (Fig.10).
Finally, in order to increase the visualisation potential, an even more accurate model for highlighting terrain anomalies was graphically elaborated using the QGIS slope analysis tool. Slope Analysis is an algorithm capable of highlighting the steepest points on the analysed territory (Brogiolo and Citter 2018, p. 601). The analysis carried out on 1cm DTMs with Z factor = 1.000000, showed quite accurately all the major height difference points, sometimes overlapping with the densest vegetation, but also, with indicators of an archaeological nature.

RESULTS
The engineered system saw its first experimental application within the archaeological context of Leopoli-Cencelle, a medieval town founded by Pope Leo IV in 854 A.D. and for several years the subject of systematic archaeological surveys by the chair of Christian and Medieval Archaeology at the University of Rome 'Sapienza' (Stasolla 2012). Even though numerous town spaces have been brought to light over the several years of research, the overall urban extension is still not entirely clear, except for the presence of a few outcropping walls that suggest the perimeters of some town units.
The investigations with the UAVIMALS system were concentrated on the top of the hill, where the civil pole with the public palace and the religious pole with the municipalage church dedicated to St. Peter were found. The survey area was chosen for the presence of thin burying layers and partially outcropping structures covered by shrub-type vegetation, representing one of the best test fields to verify the functioning of the system created. Over a period of three working days, 30 flight missions were carried out to survey both the structures already revealed in the religious and civil sectors, as well as to investigate any archaeological deposits not yet covered by the excavation investigations (Fig.3). The different flight missions were distributed within three large areas, depending on the degree of disruption each building presented in terms of soil and vegetation cover. In fact, starting from the easternmost area, the civil and religious poles almost entirely excavated, it proceeded westwards until it reached the area behind the city walls only partially investigated through excavation activities (Fig.9). This made it possible to establish the degree of accuracy of the relief of each archaeological evidence visible in the DTM produced by the instrument according to the degree of coverage present on each deposit. The results presented here, however, refer only to the westernmost area of the hill (Fig.9, yellow area), which yielded the largest number of unknown traces. The DTM obtained from the LiDAR point cloud was the starting base for the elaboration of a Slope Analysis, one of the GIS tools that most easily allows one to highlight the steepest points of a surface, in this case coinciding with possible indicators of buried archaeological deposits. The method adopted, as in some Trieste contexts (Forlin 2012;Bernardini et al., 2018) made it possible to clearly distinguish the presence of structures, linked to a series of environments, already partially identified in the excavation years preceding the remote sensing survey and in part still completely unknown. Slope Analysis revealed anomalies of a possible archaeological nature, homogeneously distributed throughout the overflown area ( Fig. 11-12). The most evident marks, those highlighted in blue (Fig.12), correspond to the wall septa of the rooms found during the archaeological excavations carried out between 2016 and 2019, near the western walls of the site (Annoscia et al., 2020).
However, the reading of the slopes also made it possible to distinguish with greater graphic clarity the vegetation indicators from the rest of the visible marks, thus highlighting additional traces (traces in red in Fig.12) that, due to their type of shape and course, can be identified as archaeological indicators of further wall fragments of buildings belonging to the urban quarters in the westernmost area of the site that is only partially known today.

CONCLUSION
UAVIMALS was a particular application for this type of sensor, never before proposed, but our experimentation showed the suitability of an instrument not designed for remote sensing but rather for obstacle avoidance, allowing us to conclude that the system is capable of providing meaningful information without requiring a large budget. The open architecture was key to further reducing costs because flight data from the drone's on-board IMU was integrated with sensor measurements and calculations were performed to produce the point cloud. The use of a dedicated IMU did not prove necessary and the accuracy of the system was sufficient for the planned application. The archaeological markers that emerged from the surface investigations were still completely covered and unknown at the time of the remote sensing in May 2020. To date, the continuation of the archaeological excavation has confirmed that many of the marks found did indeed correspond to archaeological features discovered during the recent campaigns of 2021 and 2022. In particular, the wall septa of medium and large sized rooms of the quarters of the western urban pole and some portions, not yet fully known, of the same city walls, which are difficult to access from the west precisely because of the extremely steep morphology of the plateau on which the site stands, were found (De Lellis 2015).