High resolution orthophotos and a digital surface model of the Roman city of Pollentia (Mallorca, Spain) using RPAS imagery, aerial images, and open data archives

This contribution presents an approach to generate high resolution orthophotos and a digital surface model combining RPAS imagery and ground control support derived from aerial and publicly accessible aerial data instead of organizing a specific surveying campaign. The scope of the presented research is the multi-resolution geospatial data management; its goal, to contribute to this topic by assuring geo-referencing consistency between RPAS imagery and aerial and historical orthophotos. The approach has been experimentally tested and validated in the Roman city of Pollentia (Alcúdia, Mallorca, Spain) and might be further applied to any other archeological site where aerial imagery and auxiliary data are available.


I. INTRODUCTION
Latest advances in Remotely Piloted Aircraft System (RPAS) technology, Commercial-of-the-self miniaturized cameras, and photogrammetric software have allowed the archeological community the access to this technology beyond research and experimental projects [1,2,3]. Recent works have shown the potential and utility of the outputs of these technologies (orthophoto, point clouds and Digital Surface Model (DSM)) for the management of archeological sites, both for documenting and a 3D modelling perspective [1,4,5]. Thanks to the capability of flying at low altitudes, dense point clouds and orthophotos with a high level ground sampling distance (GSD) can be generated. With that, the GSD is increased from meter (satellite imagery) or decimeter (aerial imagery) to few centimeters [6]. This improved resolution may allow to digitalize an archeological site at a level where stones may be clearly distinguished [1,6] or provide detailed scale of the structures (Figure 1).
These outputs provide a suitable alternative to total stations by relaxing requirements in terms of metric accuracy. Moreover, these technologies may also help to detect potential buried archaeological remains when multispectral orthophotos are generated. The advantage in this case is the capability to acquire data in the optimal time window for detecting soil or cropmarks anomalies [7,8].
From an archeological perspective it is important to have historical data, acquired during many years, properly geoferenced and co-registered. It is essential to have high resolution orthophoto and DSM, but they might be meaningless if such data, obtained in different time frames, are not properly registered or aligned with aerial orthophotos. An example of these archives are those provided by the Spanish Plan Nacional Ortofotografía Aérea (PNOA) program. The PNOA is a national program led by the Spanish National Geographic Institute. The goal of this program is to generate and make publicly available high resolution orthophotos and Digital Terrain and Surface Models. These products are generated for the entire country and updated every 2-3 years. The program also provides access to raw imagery and auxiliary data used to generate the aforementioned products. The auxiliary data includes the position and attitude of the raw imagery as well as the camera calibration.
The standard workflow for orthophoto and DSM generation includes the following steps: aerial triangulation, point cloud and digital surface model generation and finally orthophoto generation. The aerial triangulation step involves the estimation of exterior orientation (position and orientation) parameters as well as camera calibration parameters. These parameters are usually estimated using homologous, exterior orientation observations provided by on board Global Navigation Satellite System (GNSS) receivers and Ground Control Points (GCP). These points are usually surveyed with differential or Real Time Kinematics (RTK) GNSS techniques [9]. However, it may happen that resulting orthophoto and DSM are non-proper registered with the aerial cartographic archives. Thus, an additional step may be required to register multi-temporal and multi-resolution orthophotos using an affine or projective model. Registering the DSM may be more difficult or even impossible due to the difficulty of identifying common points.
Alternatively, several solutions have been proposed to deal with the co-registration of multi-temporal datasets during the orientation step, prior to orthophotos and DSM generation [10,11]. Both solutions focus on the automatic detection of common points between RPAS datasets [10] or between RPAS and aerial datasets [11], but differ in the way in which the ground control information is generated. The first solution generates the ground control information by performing initially the image orientation step for one of the dataset [10]. Then some images in such dataset are used as anchor images to constrain the orientation step of remaining dataset (or datasets) without using ground coordinates of common points. [11] uses planimetric coordinates of common points extracted from available orthophotos and elevations rom Digital Surface or Digital Terrain Models (DTM) [11]. The main limitation of the aforementioned solutions is that it relies on the capability to automatically identify common points between datasets. This cannot straightforward for datasets with different GSD, long temporal gaps where few common points can be detected. In this context, the reliability of the previous approaches depends of the number of detected points and their distribution.
In this paper, an approach to generate high resolution orthophotos from RPAS imagery and DSM is presented avoiding the use of GCPs measured with GNSS techniques. Nevertheless, it is still being possible to co-register orthophotos from aerial archives and RPAS imagery. The innovation of the approach is to avoid the use of GNSS based ground control points. They are replaced by GCP derived from aerial imagery and exterior orientation data. Note also that these data were used to generate open and broadly accessible aerial orthophoto archives (PNOA). The key point of the approach is the capability to use reliable ground control points triangulated from aerial images that were already used to generate aerial orthophotos. Thanks to this point, the RPAS orthophotos, as well as the DSM, will be co-registered with the PNOA orthophotos avoiding extra co-registration steps. Additionally, the proposed approach may help to potentially reduce the surveying cost by avoiding the use of coordinates of GCP surveyed with GNSS RTK techniques.

II. PROPOSED APPROACH
The proposed workflow includes several steps ( Figure  2). The first one is the generation of ground points, that is GCP and Checkpoints (CHP). The GCP are used in the triangulation step while the CHP are used only for validation purposes. In this approach, the ground points are observed both in aerial and RPAS imagery. The coordinates of these points are obtained using image triangulation techniques instead of deriving the ground coordinates of these points from the planimetric coordinates of the orthophoto and the height from the DSM. This can be done in this way because the exterior orientation of the aerial images in the archives as well as camera calibration data are already open and available to anyone upon request. The planimetric coordinates of the points derived directly from the orthophotos are very good (25 cm resolution GSD in the PNOA program). However, this approach is used because it provides better estimates of the height component than the ones provided in the DTM. Open and available DTM (also from PNOA) provides height in a grid of 5x5 m although it can be slightly improved by interpolating height from nearby cells.
The selection of the number and type of GCP is done taking into account two considerations: spatial distribution and temporal stability of points. Regarding the first point, GCP should be distributed covering all the area surveyed by the RPAS imagery. The selected points should be easily recognizable in both aerial and RPAS imagery and invariable in time. Examples of potential candidates are points belonging to fixed structures, permanent road markings. Examples to avoid are points belonging to vegetation or temporal road marking. The next step is the aerial triangulation of RPAS imagery involving the following observations: image coordinates of homologous points between images, initial orientation provided by a GNSS receiver and image and ground coordinates from GCP. The exterior orientation of every image as well as the camera calibration parameters are estimated using a bundle adjustment approach. This step is validated using root mean square error (RMSE) between the estimated coordinates of set of checkpoints from aerial triangulation process and the ground coordinates of these points. The final steps, common in many photogrammetric software tools, involve the dense point cloud, DSM and orthophoto generation.
The proposed approach has several advantages: first, the capability to fuse multi-temporal and multi-resolution datasets to obtain a historical perspective of the archeological excavations. This is done thanks to integrating all data available in a common georeferencing frame. The second advantage is that the structure-from-motion software pipelines may be used with no modification by this approach. The third advantage is the capability to have centimetric precision (less than 10 cm), allowing the clear identification of the different elements composing an archeological structure, instead of the general form that can be observed from aerial orthophotos. Another advantage is that with the proposed approach, the georeferencing accuracy of the PNOA orthophotos is kept while the resolution might be increased to few centimetres depending of the camera used and the flight altitude. The main limitation of the approach is the lack of automation on the detection of GCP. It relies on the manual identification of GCP and their image coordinates.
The methodology has been tested and validated in the Roman city of Pollentia with the aim to be extended to other archaeological sites.

III. CASE STUDY: POLLENTIA
Pollentia is located at the city of Alcudia, in the island of Mallorca (Spain). The ancient city was identified in the 19th century and continuous excavations have been carried out since 1923. Excavations have uncovered a residential area (Sa Portella), a theatre, part of the forum, several necropolises, and other remains of the city (Figure 3). Nevertheless, several areas are still covered and potential buried remains might be located also in the nearby fields beyond the current limits of the Roman city [12,13].

A. Datasets
The datasets include PNOA aerial archives and RPAS imagery. Regarding the aerial archives, three different overlapped aerial images, used to generate the 2015 PNOA orthophoto of the area together with the exterior orientation and camera calibration values, were available. The 4-band images acquired with a high-performance metric aerial camera (Vexcel UltracamXP) were taken providing a GSD of 21 cm. The RPAS imagery (5-band images) were acquired with a Micasense RedEdge multispectral camera that was flown in a fixed-wing RPAS over the Roman city of Pollentia and nearby fields. 150 images were collected corresponding to a GSD of 10 cm. No laboratory camera calibration was available beyond the information provided in the EXIF file.

B. Experimental Results
From the aerial archives, 33 points covering Pollentia and its surroundings were selected and their coordinates were triangulated. The triangulation was performed using the image coordinates of these points, the position and attitude of the aerial images and the aerial camera calibration (Figure 4). 11 of them were selected as GCP while the remaining have been selected as checkpoints (CHP) to evaluate the quality of the RPAS imagery aerial triangulation step.
The typology of the selected ground points was diverse. Most of them belonged to edges or clear visible elements from structures of the Roman city or structures from nearby fields. In addition, some points belonging to structures such as edges of swimming pools, road markings were also selected. Points subsceptible to vary such as vegetation were avoided.
The Agisoft Photoscan software was used for performing the aerial triangulation step using the ground and image coordinates of the GCP, image coordinates from homologous points together with initial orientation of the images provided by the on-board GNSS receiver of the camera. The image coordinates of the GCP and CHP were manually identified in the RPAS imagery. The ground coordinates of the GCP were input with a very low standard deviation (0.001 m) to constrain the bundle adjustment. Camera calibration parameters (focal length, principal point, radial and tangential distortions) were also estimated by the adjusment. The final step involves the generation of georeferenced maps. Orthophoto and DSM were also generated with the same photogrammetric software. The generated products have a GSD of 10 cm.
In order to evaluate the quality of the results, three different strategies have been used: relative precision evaluation, accuracy evaluation and visual analysis.
The aim of the first strategy is to evaluate the quality of the aerial triangulation, and so, assure the geometric consistency between the aerial and RPAS datasets. To do so, the analysis of the CHP residuals was the tool to evaluate the relative precision quality. The residuals of these points show a planimetric error better than 1 (aerial) GSD (0.1 m and 0.15 m) and slightly higher than 1 GSD for the height component (0.3 m). The residual for the height is worse than the planimetric components and (slightly) higher than 1 GSD and. The possible explanation is that the ground points were estimated using only three aerial images belonging to the same flight strip. The use of images from different flight strips would improve the precision of the ground points height component. And thus, it would help to decrease the height component residual below 1 GSD.
These results suggest that that aerial triangulation was not only able to estimate the orientation of the images and camera calibration parametes but to provide reliable estimates. In addition, it is also a way to assure that both aerial and RPAS datasets are geometrical consistent.
The proposed methodology relies on the accuracy provided by the Aerial (PNOA) imagery. According to specifications [14], the planimetric accuracy of the PNOA orthophoto should be better than 0.5 m (RMSE) and 1 m for the height component. However, for the seek of completeness, an assesment of the global accuracy was done. To do so, 5 different points were surveyed with differential GNSS techniques in the area of Forum ( Figure 5). The residuals between the triangulated coordinates and the differential GNSS coordinates were computed. The results show a residuals of 0.07 and 0.13 m for the planimetric coordinates. The RMSE for X, Y coordinates 0.45, within the expected accuracy of the PNOA program.   The last procedure for the quality checking of results was a visual testing. This visual checking seeks to assess the coregistration between PNOA and generated RPAS orthophoto. After the RPAS orthophoto generation, some structures belonging to Sa Portella area were manually digitalized into a vectorial layer. The overlaping of this vectorial layer with the aerial ortophot was the tool to check visually the coregistration between both orthophotos. Figures 6 and 7 show the geometric consistency of both datasets.

IV. CONCLUSIONS AND FURTHER RESEARCH
An approach to assess not only the proper geoferencing of high resolution cartography generated with RPAS imagery but also the multi-temporal co-registration using available aerial georeferenced imagery has been presented. The approach has been tested and validated with RPAS and aerial imagery of the Roman city of Pollentia. The aim of this work was not to assess the potential of RPAS imagery, for monitoring and documenting the site, already known by archeological excavations or remains, but to assess the proper integration/fusion with the available historical georeferenced data.
The approach relies on the use of ground control support derived from aerial and publicly accessible data instead of organizing a specific surveying campaign. The preliminary results show that a relative precision, between aerial and RPAS orthophotos, better than 1 GSD is achieved for the Roman city of Pollentia and their surrounding area. The results show also an absolute accuracy better than 0.5 m (RMSE) for the planimetric coordinates.
Further research will take care of a more complete validation to prove the reliability of the approach. This validation will include a comparison of the results obtained with the proposed methodology and the ones obtained using GCP coordinates measured from a specific surveying campaign. Beside this, the same methodology will be applied to additional RPAS imagery from Pollentia as well as other archeological sites.

ACKNOWLEDGMENT
This contribution is part of the results of the project Archaeology, Remote Sensing, and Archaeometry: A multidisciplinary approach to landscape and ceramics from the Roman to the Medieval period in Mallorca (Balearic Islands) (ARCHREMOTELANDS) (HAR2017-83335-P), PI: Miguel Ángel Cau Ontiveros, funded by the Ministerio de Economía, Industria y Competitividad, with contribution from the European Regional Development Fund from the European Commission. This is also part of the activities of the Equip de Recerca Arqueològica i Arqueomètrica de la Universitat de Barcelona (ERAAUB), Consolidated Group (2017 SGR 1043), thanks to the support of the Comissionat per a Universitats i Recerca del DIUE de la Generalitat de Catalunya.
The work of Catalina has been possible thanks to a Beatriu de Pinós postdoctoral fellowship from the AGAUR, with the support from the Secretaria d'Universitats i Recerca of the Departament d'Economia i Coneixement of the Generalitat de Catalunya.