*Article* **Relative Radiometric Calibration of Airborne LiDAR Data for Archaeological Applications**

#### **Christopher Sevara 1,\*, Martin Wieser 2, Michael Doneus 1,3 and Norbert Pfeifer <sup>2</sup>**


Received: 19 March 2019; Accepted: 16 April 2019; Published: 19 April 2019

**Abstract:** Airborne laser scanning (ALS) data can provide more than just a topographic data set for archaeological research. During data collection, laser scanning systems also record radiometric information containing object properties, and thus information about archaeological features. Being aware of the physical model of ALS scanning, the radiometric information can be used to calculate material information of the scanned object. The reflectance of an object or material states the amount of energy it reflects for a specific electromagnetic wavelength. However, the collected radiometric data are affected by several factors that cause dissimilar values to be recorded for the same object. Radiometric calibration of such data minimizes these differences in calculated reflectance values of objects, improving their usability for feature detection and visualization purposes. Previous work dealing with calibration of radiometric data in archaeological research has relied on corresponding in-field measurements to acquire calibration values or has only corrected for a limited number of variables. In this paper, we apply a desk-based approach in which radiometric calibration is conducted through the selection of homogenous areas of interest, without the use of in-field measurements. Together with flight and scan parameters, radiometric calibration allows for the estimation of reflectance values for returns of a single full-waveform ALS data collection flight. The resulting data are then processed into a raster reflectance map that approximates a monochromatic illumination-independent true orthoimage at the wavelength of the laser scanner. We apply this approach to data collected for an archaeological research project in western Sicily and discuss the relative merits of the uses of radiometric data in such locations as well as its wider applicability for present and future archaeological and environmental research. In order to make the approach more accessible, we have developed a freely available tool that allows users to apply the calibration procedure to their own data.

**Keywords:** ALS; amplitude; radiometric calibration; reflectance; archaeology; Sicily

#### **1. Introduction**

Airborne laser scanning (ALS) data have become a common component of the archaeological prospection toolkit. The increasing availability of such data, along with the concomitant development of expertise in their collection, processing, and visualization for archaeological investigation continues to have a significant impact on the detection and interpretation of the present remains of human activity in landscapes around the world. Often, the main sources of such information are visualizations derived from digital terrain (DTM) or digital surface (DSM) models, which are geometric models that are calculated from the range data produced during a laser scan. However, there is other information in laser pulses that, though often neglected, can be of great value to archaeological and paleoenvironmental research. In addition to range data, laser pulses reflected from a given surface also carry information about its physical properties. This information is included in the data in the form of radiometric values (e.g., the amplitude of the received laser pulse signal).

As these values are distinct from those of surface models, they can provide important additional information about the presence and composition of archaeological and environmental features. Radiometric values can be used to derive physical properties of surface objects related to the laser wavelength, independent of local lighting conditions [1] (p. 15). As they are collected, however, these radiometric values are distorted by many variables. These include factors related to the atmosphere, the instrument, its position in relation to the object being scanned, and the environmental context of the target. This makes it difficult to compare and quantify values as collected by the scanner. In order to make full use of the radiometric information, the received amplitude values collected for each laser pulse can be used to derive calibrated reflectance values so that radiometric differences of objects within and between strips are minimized [2] (p. 335), [3,4]. Once calibrated, these values can then be used for purposes such as visualization, interpretation, and detection of archaeological features.

Although a common component of ALS data capture, radiometric data sets have yet to become as integrated into the prospection toolkit as their geometric counterparts. Importantly, calibrated reflectance data in the infrared range have been shown to provide increased information regarding vegetation health, facilitating the appearance of vegetation marks associated with archaeological and paleoenvironmental features, when compared with other data sources [5–9]. Additionally, in cases where optical imagery is not concurrently recorded, reflectance-based imagery can be used as an approximation of simultaneously acquired monochromatic imagery for interpretation purposes [10] (p. 418), [11] (p. 50). While the usefulness of this information has been repeatedly demonstrated by such studies, its application in archaeological research continues to be explored in a very limited manner when compared to the uptake of other ALS-derived data products [7] (p. 161). This may be due to several factors, including the limited usability of amplitude data as acquired during flight, and difficulties involved in data processing.

While recent work has demonstrated the benefit of radiometrically calibrated full waveform (FWF) ALS data regarding increased archaeological feature detection, this work relies on the collection of simultaneous in-field radiometric measurements for calibration [5]. This may not always be a possibility (particularly when working with previously acquired data) or even a necessity for archaeological research that does not utilize and directly compare radiometric values between multiple ALS flights or to physical properties of objects. Often, single-flight ALS data sets are acquired 'readymade' by archaeologists from service providers who may have flown data collection sorties months in advance for other clients and are simply reselling the data. Moreover, whether looking for direct evidence of archaeological features or indirect indicators (e.g., vegetation marks) that may indicate buried archaeological or paleoenvironmental features, many archaeologists are usually only interested in the relative difference between the vegetation producing the marks and the surrounding environment (the contrast). Therefore, in many cases it may be sufficient for archaeologists to work with relatively calibrated reflectance values. However, even relative calibration of reflectance data requires metric information about the scan variables and target properties in order to be successfully applied.

In order to address this issue, the research presented in this paper outlines a desk-based methodology for relative radiometric calibration of FWF ALS data sets through the GIS-based selection and use of spectrally similar areas of interest (AOIs) across multiple ALS strips. This allows for the possibility of calculating reflectance values across strips without the need for simultaneously acquired in-field measurement. To make this approach more accessible to a wide variety of users, we have developed a freely available tool, based on the OPALS system [12], to assist with the calibration process. This approach can therefore be of particular use for single-flight 'off the shelf' FWF ALS data sets that contain radiometric information, may not have originally been acquired for archaeological purposes, and which lack contemporary optical imagery.

In this paper, we first outline the technical specifics of radiometric data as an element of ALS data capture, and its prior applications in archaeological research. Next, the desk-based methodology for relative calibration is presented and applied to FWF ALS data acquired from an archaeological case study area in western Sicily. Using the results, we evaluate the methodology, as well as the usability of such data as a component of archaeological research in western Sicily, which is a Mediterranean environment with considerably different environmental properties to that of previous archaeological studies dealing with the use of radiometric data. Considering this evaluation, we discuss the present and future potential of radiometrically calibrated data as a component of the archaeological prospection toolkit.

#### **2. Radiometric Data from LiDAR: Theoretical Background and History in Archaeological Research**

Radiometric calibration of amplitude data collected via ALS is a well-established process in the discipline of remote sensing and has also been explored in various capacities by the archaeological community. However, as can be seen in literature dealing with the subject (Section 2.3), the term 'intensity' is often used as a blanket term for diverse values, including uncalibrated values as collected by the scanner as well as subsequently calibrated data. Additionally, many authors have used the term 'calibration' for a wide range of attempts to compensate for distortions in uncalibrated data. The broad application of these terms can cause confusion, therefore in this paper uncalibrated values (as collected during a laser scan) are referred to as 'amplitude', and calibrated values as the 'reflectance' of materials or objects as a physical characteristic. The word 'intensity' is avoided, except when quoting the original authors of a study where the term was used. The following sections describe the background, properties, and terminology related to amplitude, reflectance, and calibration of reflectance data as well as its application in archaeological prospection.

#### *2.1. LiDAR and Radiometric Data*

ALS data capture is an optical 3D measurement technique that operates independently of local illumination, using a laser source for acquisition and measurement of a target surface [1] (p. 1). During an airborne laser scan, laser pulses of a certain wavelength are sent toward the ground surface and the backscattered radiation of each individual pulse is recorded. The scanning process results in two kinds of information:


A laser scanner operates using electromagnetic (EM) radiation at a specific, well-defined wavelength (usually either 532 nm, 1064 nm or 1550 nm). As a function of the wavelength, EM radiation interacts in different ways when it hits object surfaces. Amplitude information, as a component of ALS data capture, represents the amount of energy reflected from an object that has been scanned with a focused pulse of laser radiation [1] (p. 15). Depending on the object's material, a certain percentage of energy will be absorbed, transmitted, or reflected. The specific wavelength of the radiation being emitted by the scanner will also affect the transmission of the material being scanned. For example, the backscattered radiance from the 1550 nm short-wavelength infrared (SWIR) band has been shown to produce information about archaeological vegetation marks [5] (p. 170), which may be due to effects such as the greater absorption of SWIR radiance by water within plant materials [13]. Therefore, in theory, amplitude data contains wavelength-specific information about

the physical properties of the reflecting object. In practice, this information is distorted by several other factors that also affect the strength of the backscattered signal, such as the incidence angle or object size.

#### *2.2. Radiometric Calibration*

Material information within the received amplitude values from a laser scan can be very useful for identification of discrete objects, such as certain types of vegetation marks resulting from slightly higher levels of moisture retention or plant stress than their surrounding vegetation. Indications of plant stress may be more apparent in the near-infrared region of the electromagnetic spectrum, appearing earlier and more extensively than in the human visual range [14]. Thus, by extension, they may also appear in the amplitude values of a laser pulse in a similar wavelength. However, as it is collected the amplitude information is highly distorted. This is due to a number of variables, some of which are machine or workflow specific and thus relatively static, and some of which change with each new line or day of flight. Chief among these are the range and the influence of the laser pulse incidence angle on amplitude values, which is caused by changes in the position of the scanner relative to the surface being scanned. Figure 1 shows the influence of the incidence angle between two overlapping strips. On the border of a strip, the incidence angle is larger than in the middle, causing a visible decrease in the amplitude. Topographic variation (different range) and longer flight times (changing atmosphere) can also have a significant effect on amplitude values [3] (p. 416). Normalization of reflectance values of similar objects over an entire ALS data set (consisting of multiple parallel and cross strips) therefore requires calibration procedures that correct influences stated in the following sections [3].

**Figure 1.** Uncalibrated amplitude map of two overlapping strips. Note the extreme differences between the edge of Strip 2 and Strip 1, as indicated by the red arrows. Visualization: 1 m spatial resolution raster using the median amplitude last echo point values for each cell. Coordinate System: WGS84 UTM Zone 33N.

#### 2.2.1. Amplitude

The amplitude of the received signal is the remaining energy that is recorded by the system after the laser pulse is sent out, travels through the medium (e.g., air), is reflected by the object, and travels back (Figure 2). The amplitude, as described in Equation (1) [15,16], is therefore dependent on three types of parameters: (1) Scanner specific characteristics such as the emitted signal strength and the aperture size of the receiving optics, (2) signal loss due to the range, atmosphere and the incidence angle, and (3) target object characteristics (called the backscatter cross section).

Equation (1) Derivation of amplitude [15,16].

$$P\_r = \frac{P\_t D\_r^2}{4\pi R^4 \left\| \begin{matrix} 2 \\ t \end{matrix} \right\|\_t^2} \text{ or } \eta\_{\text{sys}} \text{ } \eta\_{\text{atm}} \tag{1}$$


**Figure 2.** Airborne laser scanning (ALS), basic principle of transmission.

#### 2.2.2. Reflectance

Reflectance is the physical property that denotes how much energy is reflected by an object's material. Depending on the wavelength, EM radiation interacts in different ways when it hits object surfaces. As noted above, a certain percentage of energy will be absorbed, transmitted, or reflected. The backscatter cross section in Equation (2) [15], combines all target parameters, including the object's material and its directionality of the scattering. It is a measure of the surface's interaction with radiation. Within Equation (2), the reflectance is part of the backscatter cross section. For surfaces with uniformly bright, diffuse reflectance behavior (Lambertian surfaces, see [17] (p. 441)), the reflectance is therefore the physical property of how much energy is reflected by a material. For non-Lambertian surfaces, the value can be interpreted as a diffuse reflectance measure at the specific measurement geometry. Thus, this value only delivers the real reflectance for surfaces with Lambertian backscatter characteristics. In practice, these could be surfaces such as asphalt or unfinished wood. In order to calculate an

incidence angle corrected reflectance, Lambert's cosine law [18] is applied, as the reflected energy decreases with the incidence angle.

Equation (2) Derivation of reflectance [15].

$$
\sigma = \frac{4\,\pi}{\Omega} \,\rho \, A\_i \tag{2}
$$


#### 2.2.3. Calibration to Derive Reflectance

The process of deriving material information in the form of the attribute reflectance from the received amplitude data is called radiometric calibration. Reflectance is derived by reformulation of Equations (1) and (2). Thus, the calibration process corrects the influences on the received amplitude and those from the object's geometry described above. Some of the parameters can be derived from laser echo information itself (e.g., range correction) or approximated using the incidence angle, which requires its neighboring points for estimation of the surface normal, scanner trajectory, and point of reflection information. All other parameters that are unknown but can be assumed to be invariant during one ALS campaign can be summarized in one constant, the calibration constant *Ccal* (The full derivation of the reflectance can be found in Appendix A.). Therefore, *Ccal* needs to be estimated during the calibration process. This is done by using areas whose reflectance values are known (either through separate measurement or estimation). Once *Ccal* is obtained, the reflectance data for each laser scan strip can be calculated.

#### 2.2.4. Absolute and Relative Calibration

The derivation of the reflectance depends on the estimation of *Ccal* and therefore on areas with known reflectance values. Depending on the source of the known reflectance values, the calibration either results in absolute or relative calibration and reflectance values [2] (p. 335):


Absolute calibration can be used in order to compare values between multiple ALS data sets flown at different periods or wavelengths, or to relate radiometric data directly to physical properties of objects [2,5]. An absolute calibration approach is clearly necessary when comparing multiple data sets from different flights and/or wavelengths, or when the intent is to connect the physical properties of an object to the measurements (i.e., when calculating plant biomass or leaf-area indices). However, when using data from a single flight that has been previously acquired or purchased from a supplier after acquisition, it is not possible to obtain simultaneous in-field measurements. Furthermore, if the intent is simply to utilize the data to quantitatively or qualitatively look for relative differences in reflectance within a given data set, absolute calibration may be unnecessary. This includes the bulk of archaeologically relevant uses of radiometrically calibrated data, in which the user is only interested in examining the relative contrast between objects and where the real-world values of an object are not necessary.

#### *2.3. Radiometric Data and Calibration in Archaeological Research*

The technical specifics of FWF ALS data capture and its application in archaeological research for topography-based purposes (i.e., as described in point 1 of Section 2.1) have, by now, been widely

reported [19–23]. While the overview of various uses of LiDAR reflectance information presented by Höfle and Pfeifer [3] (p. 417–419) summarizes its use for the improvement of ALS data processing as well as its application in classification of land use/land cover, considerably less attention has been given to such characteristics of laser pulses in archaeological research. Coren [24] visualized 'calibrated intensity' data in the area of Roman Aquieia in Northern Italy, locating a number of "anomalies already seen in the digital orthophoto". However, their archaeological significance was not discussed. Challis has used 'intensity' data in the interpretation of river valleys [6–9]. Challis' work also discussed the wide range of factors contributing to its inherent distortion. Furthermore, attempts were made to compensate for surface variation [6,7] (p. 164). From the resulting normalized image, Challis could identify cropmarks suggesting a Romano-British villa site at Cromwell that was not visible in the DTM [6] (p. 7,11). At the river Trent, similar data revealed cropmarks of paleochannels and hinted at parched crop and "highly saturated ground" [7] (p. 163).

Although both Challis and Coren applied calibration procedures to their data sets, these only accounted for a limited number of parameters. Challis [6] (p. 6) corrected for topography using a range-based formula adapted from Luzum et al. [25]. The impact was, however, minimal [6] (p. 9). This may be due to the fact the calibration process applied in the study changed only one of many distorting factors. Based on this experience, Challis proposed further research into calibration, especially with regard to the use of FWF ALS, and to make further use of the results. This includes ALS acquisition flights during the vegetation season and combining the resulting data with other prospection data [7] (p. 168).

In 2014, Briese et al. [5] presented a workflow for absolute radiometric calibration of FWF ALS data that provides calibrated multi-wavelength reflectance images useful for archaeological prospection. In contrast to the findings of Challis and Howard [7] (p. 167), Briese et al. were able to show significant improvement in usability of calibrated reflectance data acquired from an ALS scan in low-relief areas of fairly uniform vegetation. It should be noted, though, that the scan parameters, scanner type, and flight conditions were very different to that of Challis and Howard's study. The approach used in [5] is based on measurement of in-situ reference targets and recording of associated meteorological data, decomposition of the full-waveform to acquire echo parameters, and subsequent calculation of local variables and correspondence of ALS data with in-situ measurements in order to archive a calibration coefficient, which is applied to all strips [5] (p. 166). While the study showed a significant increase in the detection of vegetation marks indicating archaeological features of interest in the area of (formerly) Roman Carnuntum (Lower Austria) when compared with interpretations based on information from existing orthoimagery, the requirements of the absolute calibration process exceed the limits of data acquired 'off the shelf' or when simultaneous in-field measurements are not possible.

#### **3. Methodology**

Our workflow consists of a desk-based approach to derive the reflectance of a multi-strip ALS data set, ideally collected during a single flight period. Preferably, a geometric strip adjustment should be applied to the data set prior to application of the workflow, as was the case with our data set. The workflow used for the derivation of calibrated reflectance maps from FWF point clouds was implemented in the software package OPALS [12]. Instead of using in-field measurements, which were not available, estimated Lambertian reflectance values for given surfaces were used. However, the process for calibration is the same, regardless of the source of the reflectance values. Therefore, this approach may also be used for absolute calibration if the necessary measurements are available.

In addition to the point cloud (containing amplitude values) and estimated reflectance values for given locations, another key point attribute is the beam vector (Figure 3), which is needed to calculate the range and incidence angle. The beam vector is calculated using data from a trajectory file generated during the laser scan, which contains time, position and orientation information about the laser beam origin as acquired from the GNSS/IMU of the instrument. It is also possible, to some extent, to use approximate estimations of such values. However, this may result in less accurate or unsatisfactory calibration results and more effort in pre-processing.

**Figure 3.** Illustration of the beam vector and normal vector of a laser pulse. The beam vector is calculated using absolute position and orientation information about the laser beam origin.

The calibration process is divided into three parts: (1) Pre-computation of the calibration attributes, (2) derivation of the calibration constant (*Ccal*), and (3) reflectance calculation and generation of a reflectance map from the calibrated data in the point cloud (Figure 4):


To apply the workflow, we began by creating a shapefile to digitize AOIs in a GIS environment, using a raster map created from the uncalibrated amplitude values in the point cloud as a reference. In our case, we used an asphalt road and assumed a standard reflectivity of 0.2, which has been shown to be the mean reflectivity for worn asphalt surfaces at the wavelength of the scanner used in our study [5] (p. 337), [28]. While such surfaces vary due to factors such as age, moisture, and material inclusion, assuming a mean value for calibration was valid enough, as we were only interested in relative values. What was more important was that the AOI areas chosen were similar in terms of wear and material. Each AOI was digitized as a separate shape with its own corresponding attributes in the attribute table. In addition to the unique identifier, the only other necessary field in the shapefile table was a numeric field in which the estimated Lambertian reflectance value for the material within the AOI was input. An initial pre-computation (1) was run to acquire the necessary attributes for

each echo. The OPALS module RadioCal [29] was then used to acquire the calibration constant in conjunction with the data from the AOI shapefile (2). RadioCal was then re-run in order to apply the reflectance calculation to all echoes in the data set (3). The resulting point cloud was used to generate the calibrated reflectance map in the OPALS module Cell [30], using the median reflectance values from the last return points in each cell of the raster grid. A batch file containing the necessary commands was used to automate all stages of the process subsequent to the generation of the AOI shapefile.

**Figure 4.** Workflow for relative calibration of reflectance data.

**Figure 5.** The Prospecting Boundaries project area, western Sicily. Background image © ESRI.

#### **4. Case Study and Data Acquisition**

#### *4.1. The Mazaro River Corridor, Western Sicily*

In order to evaluate the applicability of the relative calibration process for archaeological research, we applied the method to data collected during the Prospecting Boundaries project. The Prospecting Boundaries project [31], in cooperation with the Soprintendenza per i Beni Culturali ed Ambientali di Trapani, seeks to explore the diversity of land use in the region along the Mazaro river corridor using integrated archaeological prospection techniques to investigate the present remains of past activity in the landscape. These techniques include active and passive remote sensing data collection and interpretation, geophysical prospection, surface survey, and geoarchaeological evaluation. Centered on the Mazaro river, the project encompasses an area of roughly 70 km2, stretching inland from the river mouth at modern-day Mazara del Vallo (Figure 5). The project is concerned with human activity during all periods in the region, working backward from the present to document and deconstruct modern and historical land use in order to try to connect the relict fragments of prior human activity to continuity and change in the wider landscape. The hinterlands of Mazara del Vallo are a prime location in which to explore these topics, as they contain a rich matrix of known archaeological sites from the Upper Palaeolithic to the modern era [32–38]. However, while areas in the wider region of western Sicily have been the subjects of more systematic field surveys in recent years (e.g., [39,40]), the area around the Mazaro has received comparatively little systematic attention. This also makes it an ideal location to test the applicability of new integrated archaeological prospection techniques.

Modern land use in our project area has a variable impact on the surface visibility of the archaeological record. In the southern and northern parts of the project area, land use is largely given over to agricultural production. This includes the growing of cereal crops, vineyards, and olive and fruit trees. The central part of the project area with its rocky limestone outcrops, known locally as *sciara* (pl. *sciare*), is mainly used for pasturage and quarrying, which seem to be a common practice in the region from antiquity to modern times. A switch to more intensive agricultural and resource extraction practices since the 1950s, as evidenced by historical aerial imagery, has had a large impact on surface visibility in the region [41,42]. Furthermore, land fragmentation further reduces the visibility of manifesting vegetation marks between irregular fields containing discontinuous vegetation (Figure 6). Thus, the environmental conditions and land use in this study area are heterogeneous and noticeably different to those in the abovementioned studies.

**Figure 6.** (**a**) The sciara di Mazara del Vallo; (**b**) land fragmentation resulting from the subdivision of parcels and agricultural diversity. Photographs: C. Sevara.

Due to the range of human activity and the physical properties of the modern landscape, archaeological sites tend to appear in a variety of direct and indirect ways in our study area. For instance, features of interest manifesting as rock-cut or rock-worn relief objects in ALS-derived digital elevation models (DEMs) generally tend to appear in the southern and central parts of our project area where the land use is more given over to pasturage and quarrying, and the Miocene/Pliocene rocks run closer to

the modern ground surface. In the northern part of the project area, topographic data is more indirectly useful, providing general information about the physiographic environment. This diversity also prompts us to look for new methodological solutions that can support our research goals. Therefore, one of the key goals of this research is to evaluate the applicability of radiometric data sets in extremely heterogeneous environments such as our case study area along the Mazaro.

#### *4.2. Data Sets*

A single flight FWF ALS data set was acquired specifically for the project. Airborne Technologies GmbH (ABT), a commercial provider of ALS surveys with experience in data acquisition for archaeological research, conducted the survey on the morning of 21 February 2016 using a Riegl LMS-Q680i scanner [42] (p. 618), which operates at a wavelength of 1550 nm [43]. Laser data and corresponding orthoimagery were collected in 26 longitudinal strips with an overlap of 20%, and two cross strips, resulting in an average unfiltered point density of 16 points per square meter (Table 1). Initial postprocessing, including strip adjustment and calculation of the 3D point cloud, was carried out in Terrascan [44]. While amplitude data was recorded during the flight, no simultaneous ground acquisition of reflectance values was performed.

**Table 1.** ALS and orthoimage acquisition parameters. Acquisition date: 21 February 2016 0746-124. Source: Airborne Technologies, GmbH/Riegl Laser Measurement Systems.


#### **5. Application and Results**

#### *5.1. Results of the Calibration Process*

The calibration workflow was applied to our entire 60 km<sup>2</sup> dataset along the river Mazaro. A substantial minimization of per strip differences can be noted in the calibrated reflectance image as compared to one using the amplitude values as acquired from the scanner. When calibrated, the values for discrete objects become more uniform (Figure 7). While the differences between strips are not altogether eliminated, they are nevertheless greatly minimized.

This can be observed in the values for a highly reflective feature, such as the modern quarry illustrated in Figure 8. The area shown was used to calculate pixel wise differences of amplitude and reflectance maps from two overlapping strips, and a spot comparison was calculated for the two points illustrated in the figure. Point 1 was calculated from the surface of the modern quarry strip, while Point 2 was calculated from an asphalt road surface similar to that used as an AOI for the calibration process. The amplitude differences deviated significantly in comparison to the calibrated reflectance, as seen in the histogram displayed in Figure 9. The amplitude values were normalized by 600 (i.e., the maximum amplitude within the strips) for comparison purposes. The values for amplitude and calibrated reflectance for each point in three overlapping strips are provided in Table 2.

**Figure 7.** (**a**) Uncalibrated amplitude map of strips; (**b**) calibrated reflectance map of strips. Visualization: 1 m spatial resolution raster calculated using the median value of each cell. Coordinate system: WGS84 UTM Zone 33N.

**Figure 8.** Calibrated reflectance map showing locations of spot comparisons. Point 1: Modern quarry strip. Point 2: Asphalt road. Visualization: 1 m spatial resolution raster calculated using the median value of each cell. Coordinate system: WGS84 UTM Zone 33N.

**Figure 9.** Histogram of differences of two strips. Red: Amplitude scaled by 600. Blue: Reflectance.



#### *5.2. Archaeological Applications*

The calibrated reflectance data in the point cloud have provided useful information about archaeological features within our project area. A 50 cm spatial resolution raster calculated from the reflectance values in the last echoes of the calibrated ALS point cloud was used for the specific purpose of examining the viability of the data for feature detection via visual inspection of the data and comparison with RGB imagery and terrain-derived visualizations. From an archaeological standpoint, we find the reflectance data to be useful in our project area in two main ways: For direct detection of rock cut archaeological features, and for the indirect detection of features through the proxy of vegetation and soil marks.

Rock-cut and rock-worn features are some of the most prominent classes of features in our project area that are detectable in remote sensing data. In addition to the fairly prominent modern quarry sites mentioned above, other subtler rock-cut and rock-worn features could be delimited in the reflectance map. This included older quarries, foundations of habitation structures, and relict evidence of transport in the form of vehicle ruts worn into the exposed bedrock. This information has proven particularly helpful in differentiating between human made rock cut structures, such as the foundations of houses and other buildings, and the geological features visible in the area. The builders of these structures took advantage of the natural geological framework of the area, carving their dwellings into the softer components of the bedrock and following the rectilinear structure of the geological substrate. As the human-made structures often follow the geological pattern and are of similar sizes and shapes to the natural features, they can be difficult to delimit in remote sensing data. Unlike the quarries, many of these multi-room structures are also comparatively shallow, at as little as 20 cm below the surrounding

rock. However, like the modern quarries, the rock cut structures and their fill exhibit more uniform reflectance than their surrounding natural counterparts. This makes them easier to differentiate from the surrounding natural geology, as is the case with the structures in Figure 10, which depicts a group of rock-cut structures in the central part of the image (Figure 10a–c)

**Figure 10.** Rock-cut structure foundations (representative structures are to the left of red dots). (**a**) Contrast stretched calibrated reflectance image of rock cut structures, which can be identified as more homogenous rectangular objects amongst the bedrock outcrop; (**b**) contrast stretched RGB orthoimage; (**c**) multiple visualization (Multi hillshade, Slope, Sky-view factor, positive openness) of a 50 cm spatial resolution digital terrain model (DTM) specifically filtered to preserve archaeological structures; (**d**) ground image of one of the rock-cut structure foundations (position indicated by red arrows). Coordinate system (**a**–**c**): WGS84 UTM zone 33N. Photograph (**d**): C. Sevara.

Additionally, the reflectance values clearly delimit the extensive groups of animal-drawn vehicle ruts that crisscross the exposed bedrock in our project area (Figure 11). These anastomosing networks are relicts of the former landscape infrastructure in the region, and serve as important indicators of trade, transport, and other movement throughout the area. While such networks have been extensively studied in other areas, such as Malta and the Azores [45–47], their history has yet to be explored in western Sicily. Although the route relicts sometimes reach lengths of over 100 m, the individual tracks do not always manifest in terrain-derived visualizations of the ALS data, possibly due to a combination of effects that include their narrow feature dimensions (individual tracks have widths of approx. 50 cm

and depths of up to approx. 60 cm), filtering and interpolation of the elevation models. Some are shallow and appear to have been less frequently used, while others are set much deeper into the bedrock. Furthermore, they are often filled with soil and/or emergent vegetation, both of which could contribute to their low relief manifestation (Figure 11d). However, the terrain-based visualizations have other advantages, such as depicting the outlines of the more topographically prominent of the linear features, which may indicate that particular strands in the network were used more heavily than others. Figure 11 shows the complimentary nature of the terrain, orthoimage, and reflectance-based visualizations. While the terrain visualization indicates the direction, extent, cross-cutting, and depths of the more prominent parts of the network, the reflectance and RGB images detail individual tracks and allow for an estimation of the relative use-wear of tracks based on their color and reflectance. Therefore, we found it advantageous to use the reflectance map in conjunction with other data sets to delimit the properties of these features.

**Figure 11.** Cart tracks made by animal-drawn vehicles. (**a**) 50 cm spatial resolution contrast stretched calibrated reflectance image of track network; (**b**) 50 cm spatial resolution contrast stretched RGB image of track network; (**c**) multiple visualization (Multi hillshade, Simple Local Relief Model, Sky-view factor) of a 50 cm spatial resolution DTM specifically filtered to preserve archaeological structures; (**d**) ground image of one of the more ephemeral cart track pairs. Note that while the DTM visualization depicts the general outline of the deeper features in the network, the individual tracks for all but the deepest of the features are only visible in the reflectance and RGB images. Coordinate system (**a**–**c**): WGS84 UTM Zone 33N. Photograph (**d**): C. Sevara.

In some circumstances, radiometrically calibrated images are more informative than RGB orthoimagery for crop mark visibility [5]. However, this was not the case for the data acquired in our project area. This may be due in part to the time of year the ALS data were collected, but it was also due to the general environmental and agricultural conditions in the region. Although no new archaeological soil or vegetation marks were identified based solely on the reflectance data, we observed comparable marks in the reflectance map to those identified in the RGB imagery. At the site of Guletta, a key research site for our project, data from remote sensing, geophysical prospection, intensive surface survey, and artifact analysis indicate the presence of a multi-ditched archaeological structure occupied between the mid-2nd and 1st millennium BC [41,48]. Based on this information, we interpret the marks in the image as remnants of structures constructed during various phases of the site's occupation. Here, the calibrated reflectance data have helped to delimit soil and vegetation marks indicating the location of the prehistoric diched enclosure, although the reflectance data alone did not provide a significant amount of new information about these features (Figure 12). Nevertheless, it can be seen from this example that, like the RGB image, the reflectance image contains information about the vegetation contrast that is indicative of archaeological features. This is, therefore, a good example of how reflectance data can also be used as a substitute for simultaneously acquired orthoimagery in cases where orthoimagery is not available.

**Figure 12.** Vegetation marks indicating the presence of a ditched settlement structure. (**a**) 50 cm spatial resolution contrast stretched calibrated reflectance image of vegetation marks indicating ditched structure; (**b**) 50 cm spatial resolution contrast stretched RGB image of vegetation marks indicating ditched structure. Red arrows indicate location of vegetation marks. Coordinate system: WGS84 UTM Zone 33N.

#### *5.3. Paleoenvironmental Applications*

Reflectance data have also proven to be a useful source of paleoenvironmental information in the project area. Paleomeanders of the Mazaro are in evidence in the northern part of our project area, where the two spring-fed streams that form the Mazaro converge into the agricultural plain. Here, the calibrated reflectance image clearly indicated former channels of the Mazaro. While some could be attributed to the channelization of the river in the late 1960s and early 1970s as corroborated by historical aerial imagery from the time, others seem to predate this activity. While the presence of such paleochannel indicators in the *sciara* area, where the Mazaro runs deep in its river valley, was not highly evident, the reflectance map indicates further paleodrainage areas that may have once served as aggregation points in an area where access to water has long been a primary concern (Figure 13a,b).

**Figure 13.** Paleoenvironmental features. (**a**) 50 cm spatial resolution contrast stretched calibrated reflectance image of the sciara di Mazara del Vallo, in which a paleodrainage can be observed running through diverse land types; (**b**) 50 cm spatial resolution contrast stretched RGB image of the sciara di Mazara del Vallo depicting the same location as (**a**). Red arrows indicate locations of paleochannel/paleodrainage. Coordinate system: WGS84 UTM Zone 33N.

#### **6. Discussion and Future Research**

#### *6.1. Merits and Limits of the Calibration Approach*

Relative radiometric calibration of amplitude information can increase the usability of LiDAR data by providing normalized reflectance values for similar features. This is useful for the delineation of human-made and natural features within an individual data set. As the resulting reflectance values are relative and based on estimation, they may not be able to be quantitatively compared to other reflectance data sets in the same area or directly to physical properties of the target objects themselves. However, this should not generally pose a problem for many applications. This includes those in which the user wants to qualitatively or visually assess the data in some way, generate true, shadow independent orthoimagery, manually or automatically detect and extract features within a single radiometric data set, perform statistical analyses, or calculate vegetation indices through a combination of the reflectance map with orthoimagery (since vegetation indices are generally relative calculations). Therefore, data generated using relative values can be applicable to a wide range of archaeological and environmental studies.

One thing that the approach we apply does not compensate for is significant change in the reflectance of an object due to changing environmental conditions such as moisture. If, for example, an archaeological feature is partially scanned while the ground surface is very dry, and partially scanned while the ground surface is very wet, there will naturally be a change in the reflective properties of the object. Furthermore, depending on the EM wavelength of the scanner, more radiation may be absorbed by the wet surface, meaning that less will be reflected. While this will generally not be an issue for data collected in a single flight, it may be problematic for data collected over an area when significant weather changes affect ground surface moisture levels during the day, or over multiple days (e.g., if the ground surface is very dry on one day that it is scanned, and very wet the next). This could be an issue for off-the-shelf datasets in which flights covering the project area are performed over multiple days. Scanner wavelength, acquisition times, and dates should be checked, and if no environmental data are present then historical weather records for the area in question should be consulted. In such cases, calibration could be performed separately for each epoch.

Currently, the calibration process requires access to proprietary ALS data processing software. Therefore, we have provided an interactive batch script and instructions for radiometric calibration in OPALS. Users who do not have access to a full version of OPALS can download the software and a key which unlocks the modules needed for radiometric calibration for a limited period. The batch file, test data, and information about how to download OPALS can be found at the University of Vienna Aerial Archive webpage [49,50]. In the future, we plan to release a freely available specialized tool with a graphical user interface that will streamline the calibration process. This application, currently in development, will be available for download in mid-2019 at the abovementioned website [49].

#### *6.2. Archaeological and Environmental Potential*

In addition to providing a workflow for calibration, our research provides further evidence that calibrated reflectance data can be of particular value in archaeological research. One of the chief merits of the approach outlined here is in its ability to unlock unique information about archaeological and environmental features that is stored within certain unexploited parameters of ALS data. In this way, ALS data sets can be seen as archives of not only vegetation and terrain, but also potential reflectance data indicating locations and properties of features that may no longer be extant. Therefore, the 'serendipity effect' [51] can also apply to archival ALS data.

In the context of our study area, archaeological features detectable using remote sensing techniques tend to manifest mainly as rock-cut and rock-worn structures and earthworks. Visualizations derived from the calibrated reflectance data have been very helpful in delimiting such features, providing clear differentiation between the texture of the structures, their fill and the surrounding areas compared to other data sets of similar spatial resolution. This is one way in which such data sets can be useful in environments such as ours, where it can be difficult to delimit visible features within the complex, noisy texture of exposed bedrock areas such as the *sciare*. As laser scanners are themselves a source of radiation, reflectance maps have the additional benefit of providing information similar to that of orthoimagery but without any shadows. This can be particularly advantageous for delimiting features that may be partially occluded by shadowing in orthoimagery. Thus, the calibrated reflectance maps provide complimentary information to terrain-based and orthoimage datasets.

The reflectance data in our study were less productive regarding detection of archaeological soil and vegetation marks in our project area. Although much of this may be due to the type of vegetation and terrain in our study area, it may also partly be due to the time of year in which the ALS data were collected. It would be worthwhile to investigate the usability of data collected at other times of the year and at different wavelengths to isolate the potential reasons for this. Nevertheless, we do see archaeological vegetation and soil mark manifestation on par with that seen in the corresponding orthoimagery, and a particular applicability for paleochannel/paleodrainage detection in both the *sciare* and the inland agricultural environments.

Regardless of the feature type, from a feature detection perspective, the ground sample distance (GSD) of a data set must be relatively high in order to clearly delineate small features and pick up subtle variations in reflectance value. The high density of returns in our ALS data set allowed us to generate a reflectance map at a GSD of 50 cm. This spatial resolution allows for the generation of data sets in which the resolution is high enough to identify fairly discrete features. Such resolution may not be possible for data collected at coarser ground spacing, although this should not limit the usability of reflectance information for the detection of features of larger size, including paleoenvironmental features.

#### *6.3. Reflections on the Future of Radiometric Data in Archaeological Research*

Through the development of a straightforward approach for relative calibration of amplitude values in LiDAR data, we hope to provide a wider range of users with access to an unexploited facet of information within their laser-based data sets. We see this parameter as being essential to a number of applications, yet chronically underused due to either the difficulty in processing it or an ignorance of its potential. While we have so far used our radiometrically calibrated data set for fairly conventional purposes, we see further potential for such information as both direct and indirect aids to archaeological and environmental research. These uses include using reflectance values as parameters for the improvement of filtering routines for ALS data, where stone surfaces (representing surfaces devoid of vegetation containing archaeological features, such as incised cart ruts or stone walls) might be automatically excluded from filtering. Reflectance values can also be useful as parameters for automated feature detection applications [3] (p. 417), [52]. In our project area, we see particular potential for such applications with regard to automated detection and classification of modern and historic land use such as quarrying, which has significantly affected the modern ground surface. We envision this as part of a routine incorporating modern and historic data sets for continual monitoring of land use change, and as a way to inform our ongoing historic landscape characterization of the region.

There may also be potential in the combination of reflectance data with simultaneously acquired optical imagery to produce an approximation of false color infrared imagery [11] (p. 50), which could yield an increase in the detection of archaeological features, particularly vegetation marks. The success of such an approach would be dependent on several factors, including the laser wavelength and time of year of the flight. Image fusion techniques, such as those proposed by [53], may also be useful ways to combine reflectance map data with other remote sensing and geophysical prospection data sets. Additionally, reflectance data could be used for the calculation of vegetation indices [54] (p. 253), mitigating effects related to BRDF (bidirectional reflectance distribution function) [55]. Statistical calculations based on reflectance data could also be used for a number of other purposes, e.g., in order to calculate percentages of coverage of certain types of land cover over a given project area.

Calibrated reflectance maps may also be useful for terrestrial laser scanning applications in archaeology. For instance, as true orthoimages can be produced independent of local lighting conditions, reflectance data may be useful for generation of illumination-independent imagery for documentation and feature detection in low-light environments such as building interiors and caves. Thus, we see many applications for reflectance data, both as a standalone information source and in concert with other prospection information. These and numerous other applications for reflectance data could be very useful not only in the context of archaeological research in the Mediterranean, but for a wide range of disciplines working in a variety of physiographic environments.

#### **7. Conclusions**

In this paper, we have presented a desk-based methodology for relative radiometric calibration of amplitude information in ALS data sets and demonstrated its applicability for archaeological research. From our perspective, the lack of uptake in the use of radiometric data in archaeological research so far may be due to two issues: ignorance of its potential and difficulty in processing. We have addressed both of these issues, and have provided a workflow and freely available tools to help users to apply the approach to their own data sets. This approach provides an additional, accessible way to exploit often unused information already present in many LiDAR scans. Therefore, it is not only applicable to newly collected data, but also to the growing archive of LiDAR data sets that have a radiometric component. Furthermore, we see great potential for relative calibration of reflectance data not only in archaeological research, but in a wide array of disciplines that utilize LiDAR data for environmental analysis.

One of the chief merits of calibrated reflectance data lies in the potential to provide an illumination-independent true orthometric image that approximates an orthophotograph in the EM range of the scanner. We have investigated the use of such information in the context of western Sicily and found that reflectance data clearly provide us with important supplemental information in archaeological and paleoenvironmental contexts. This includes narrow linear features that are difficult to detect in corresponding terrain data, rock-cut habitation, and resource extraction features, paleodrainages, and paleomeanders indicative of the location of former water sources. What makes the reflectance data particularly useful as a comparative data set is that they are based on wholly different data values than visualizations derived from terrain values. While the terrain-derived visualizations highlight relief features in unique ways, they are all based on variations of the same core values. Reflectance values, on the other hand, are based on different physical properties. Therefore, they can not only act as complimentary data sets to terrain-based visualizations, they can provide independent information about features whose reflectance may be uniform even if their relief displacement may be low. Furthermore, even when relatively calibrated, reflectance data have the potential to provide useful information for the improvement of LiDAR filtering, quantitative analyses, and automated feature detection. In light of this, we see reflectance data as being an important addition to the archaeological prospection toolbox in Mediterranean environments and beyond, providing additional ways to take even greater advantage of the information already present in LiDAR data.

**Author Contributions:** Paper Conceptualization C.S., M.D., M.W., N.P.; Writing—Original Draft C.S., M.W., M.D.; Methodology, M.W., C.S. Analysis, C.S, M.W., M.D. Visualization, C.S., M.W.; Writing—Review and Editing C.S., M.W., M.D., N.P.; Software M.W., C.S.; Funding Acquisition, M.D.

**Funding:** This research was funded by a grant from the Austrian Science Fund (Project ID: P-28410).

**Acknowledgments:** The Prospecting Boundaries project is carried out in cooperation with the Soprintendenza per I Beni Culturali ed Ambientali di Trapani. We would like to thank our project partner Sebastiano Tusa, Sicilian Regional Assessor of Cultural Heritage, and Dott.ssa Rossella Giglio, for their continued support of our research project. Programming of the software tool was carried out by the Vienna University of Technology Department of Geodesy and Geoinformation - Research Groups Photogrammetry and Remote Sensing. We would like to thank Dr. Geert Verhoeven, the anonymous reviewers and the special issue editors for their time and constructive comments on earlier drafts of this paper.

**Conflicts of Interest:** The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

#### **Appendix A**

The process of radiometric calibration is based on the following LiDAR adapted formulation of the radar equation, which describes the relation between the transmitted *Pt* and the received power *Pr* of an ALS system:

$$P\_r = \frac{P\_t D\_r^2}{4\ \pi\ R^4 \ \beta\_t^2} \text{ or } \eta\_{sys}\ \eta\_{atm}$$

This formulation considers the following influencing factors: The receiver aperture diameter *Dr*, the range between the sensor and the target *R*, the laser beam divergence β<sup>2</sup> *<sup>t</sup>* , the backscatter cross section σ as well as losses occurring due to the atmosphere or in the laser scanner system itself, i.e., a system and atmospheric transmission factor η*sys* and η*atm* respectively. The backscatter cross section, which is given by the formula below, combines all target parameters such as the laser foot area *Ai* (i.e., the size of the area illuminated by the laser beam), the reflectivity ρ, and the directionality of the scattering of the surface Ω.

$$
\sigma = \frac{4\,\pi}{\Omega} \,\rho \, A\_i
$$

All parameters which are unknown, but can be assumed to be invariant during one ALS campaign, can be summarized in one constant, the so-called calibration constant *Ccal*. This is eventually mathematically described as:

$$\mathcal{C}\_{\rm cal} = \frac{\beta\_t^2}{P\_t \ D\_r^2 \ \eta\_{\rm sys}}$$

Combining the first and third equation the formula for the backscatter cross section σ can be written as:

$$
\sigma = \frac{\mathcal{C}\_{\text{cal}} \not\equiv \mathcal{R}^4 \not\stackrel{\mathcal{P}\_{i\mathcal{S}\_{p\bar{i}}}}{\eta\_{\text{atm}}}}{\eta\_{\text{atm}}}
$$

It is mentioned that in the above equation the actual received power *Pr* is replaced by the term *<sup>P</sup>*ˆ*iSp*,*<sup>i</sup>* . This is due to the fact that received power *Pr* is proportional to the product of the amplitude *P*ˆ*<sup>i</sup>* and the echo width *Sp*,*i*.

Source: [15,16,29].

#### **References**


© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
