remote sensing Article An Operational Radiometric Correction Technique for Shadow Reduction in Multispectral UAV Imagery Xavier Pons and Joan-Cristian Padró * Grumets Research Group, Departament de Geografia, Edifici B, Universitat Autònoma de Barcelona, 08193 Bellaterra, Catalonia, Spain; Xavier.Pons@uab.cat * Correspondence: JoanCristian.Padro@uab.cat; Tel.: +34-935-814-343 Abstract: This study focuses on the recovery of information from shadowed pixels in RGB or multispectral imagery sensed from unmanned aerial vehicles (UAVs). The proposed technique is based on the concept that a property characterizing a given surface is its spectral reflectance, i.e., the ratio between the flux reflected by the surface and the radiant flux received by the surface, and this ratio is usually similar under direct-plus-diffuse irradiance and under diffuse irradiance when a Lambertian behavior can be assumed. Scene-dependent elements, such as trees, shrubs, man-made constructions, or terrain relief, can block part of the direct irradiance (usually sunbeams), in which part of the surface only receives diffuse irradiance. As a consequence, shadowed surfaces comprising pixels of the image created by the UAV remote sensor appear. Regardless of whether the imagery is analyzed by means of photointerpretation or digital classification methods, when the objective is to create land cover maps, it is hard to treat these areas in a coherent way in terms of the areas receiving direct and diffuse irradiance. The hypothesis of the present work is that the relationship  between irradiance conditions in shadowed areas and non-shadowed areas can be determined by  following classical empirical line techniques for fulfilling the objective of a coherent treatment in Citation: Pons, X.; Padró, J.-C. An both kinds of areas. The novelty of the presented method relies on the simultaneous recovery of Operational Radiometric Correction information in non-shadowed and shadowed areas by the in situ spectral reflectance measurements Technique for Shadow Reduction in of characterized Lambertian targets followed by smoothing of the penumbra area. Once in the lab, Multispectral UAV Imagery. Remote firstly, we accurately detected the shadowed pixels by combining two well-known techniques for Sens. 2021, 13, 3808. https://doi.org/ 10.3390/rs13193808 the detection of the shadowed areas: (1) using a physical approach based on the sun’s position and the digital surface model of the area covered by the imagery; and (2) the image-based approach Academic Editor: using the histogram properties of the intensity image. In this paper, we present the benefits of Fernando Vicente-Guijalba the combined usage of both techniques. Secondly, we applied a fit between non-shadowed and shadowed areas by using a twin set of spectrally characterized target sets. One set was placed Received: 29 August 2021 under direct and diffuse irradiance (non-shadowed targets), whereas the second set (with the same Accepted: 21 September 2021 spectral characteristics) was placed under diffuse irradiance (shadowed targets). Assuming that the Published: 23 September 2021 reflectance of the homologous targets of each set was the same, we approximated the diffuse incoming irradiance through an empirical line correction. The model was applied to all detected shadowed Publisher’s Note: MDPI stays neutral areas in the whole scene. Finally, a smoothing filter was applied to the penumbra transitions. The with regard to jurisdictional claims in presented empirical method allowed the operational and coherent recovery of information from published maps and institutional affil- shadowed areas, which is very common in high-resolution UAV imagery. iations. Keywords: shadow detection; shadow reduction; field spectroradiometry; radiometric correction; UAS; UAV; drone Copyright: © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article 1. Introduction distributed under the terms and conditions of the Creative Commons Since the beginning of optical remote sensing (RS), much of the scientific literature Attribution (CC BY) license (https:// has been devoted to shadow hindrances: first, in aerial; later, in satellites; and recently, creativecommons.org/licenses/by/ in drone-acquired imagery. The shadow problematics are different at these three scales 4.0/). mainly due to the spatial resolution (SR), but they all share the common trait that the Remote Sens. 2021, 13, 3808. https://doi.org/10.3390/rs13193808 https://www.mdpi.com/journal/remotesensing Remote Sens. 2021, 13, 3808 2 of 16 shadowed pixels are darker areas where the surface radiance is lower due to a lack of direct illumination, which often causes misinterpretations in reflectance estimation [1]. In an ideal de-shadowing process, a surface located in non-shadowed conditions and the same surface located in shadowed conditions should result in the same surface reflectance value. This is because the reflectance factor is the ratio of the radiant flux actually reflected by a sample surface to that which would be reflected into the same reflected-beam geometry by an ideal Lambertian surface [2,3]. Topographic shadows [4–8] and cloud shadows [9,10] are the main focus of the shadow problematic correction both in satellite optical RS imagery with moderate SR (coarser than 100 m), such as MODIS [11] or AVHRR [12], as well as in RS with high SR (coarser than 10 m), such as Landsat [13] or Sentinel-2 [14]. In very-high SR imagery (coarser than 0.25 m) from either satellite platforms, such as WorldView [15], or airborne sensors, such as AHS-160 [16,17], CASI [18] or DMCIII [19], the focus has additionally been the shadows projected by superficial elements such as buildings and trees [20–23]. In drone (also known as unmanned aerial vehicle (UAV)) imagery, new challenges are added to those previous problematics due to the ultra-high SR (up to, typically, 0.01 m), such as that provided by sensors such as MicaSense RedEdege [24], Parrot Sequoia [25], Mapir [26] or MAIA [27], making it necessary to adapt inherited airborne and satellite shadow detection and reduction techniques to the unprecedented image detail [28]. The partial or total occlusion of direct light from a source of illumination [23] causes the cast shadow of an object. There are two main classes of shadows [23,29,30]: (a) cast shadows, which are projected by the object towards the direction of the light source; and (b) self-shadows, which are in the areas of the object not illuminated by direct light (pixels with incident angle ≥ 90◦ ). Within a cast shadow, the area called umbra is where direct light is completely blocked by the object, and the area called penumbra is where direct light is partially blocked by the object (the transition zone between umbra and not-umbra). The width of the projected penumbra can be delimited using trigonometry, accounting for the angular width of the illumination source and the distance between the object and the surface onto which the shadow is projected; this is especially relevant for the top shadows of tall elements (e.g., skyscrapers) [21]. However, only umbra-cast shadowed pixels are tackled by most shadow detection and reduction methods because ultra-realistic object modelling is usually required for the accurate delimitation of the penumbra-cast shadow [23]. It is worth noting that material transmissivity is an additional problem in shadow removal [1], especially in UAV ultra-high SR. The intrinsic complexity of each single object light transmission is avoided in this work, although it has been approached in specific studies at leaf scale [31]. The radiometric theoretical basis of optical multispectral remote sensing in the visible and near-infrared regions (VNIR), i.e., between 400 nm and 900 nm, is based on the measurement of the spectral surface reflectance (ρλ ), in order to be comparable when imagery is acquired on different dates and/or times, by different sensors, and in different illumination conditions. Imagery is converted to ρλ at pixel level, i.e., the spectral radiance (Lλ , W·m−2 ·sr−1 ·µm−1 ) reflected by a surface located in a given pixel with respect to the total spectral irradiance (Iλ , W·m−2 ·µm−1 ) received [32], which are both wavelength- dependent. The total spectral irradiance in non-cloud-covered conditions reaching a given surface, which is composed of the direct solar irradiance (Idirectλ ), a less important component of diffuse irradiance (skylight) scattered by the atmosphere (Idiffuseλ ), and a still much lower component of irradiance coming from the reflection by surrounding elements (Ireflectedλ ), is very difficult to estimate due to the intrinsic characteristics of the unknown materials around the given surface. The estimation of diffuse radiance is key for the correction of shadowed pixels, where it is the main source of irradiation [33]. In UAV imagery, hyperspectral and multispectral sensor systems implement specific downwelling light sensor (DLS) devices [34], which obtain direct information to calculate the at-sensor reflectance. Moreover, UAV sensor manufacturers provide calibration panels to perform empirical line corrections [35] to obtain spectral surface reflectance. Both methods are Remote Sens. 2021, 13, 3808 3 of 16 commonly used in UAV frame cameras to obtain dimensionless reflectance values. Field spectroradiometric measurements are usually the ground truth for validating the spectral surface reflectance imagery resulting from radiometric corrections of satellite, airborne or UAV optical data [36–40]. The spectral surface reflectance of a given material does not change when located in shadowed conditions (receiving only diffuse irradiance) or in non-shadowed conditions (receiving mainly direct irradiance): what changes is the reflected spectral radiance because the surface receives less total spectral irradiance, thereby appearing darker in shadowed conditions. The hypothesis of the present work is that the relationship between irradiance condi- tions in shadowed areas and non-shadowed areas can be determined by following classical empirical line techniques [41–43] for fulfilling the objective of a coherent treatment in both kinds of areas. To this end, we present a simple and operational radiometric correction technique based on the usage of twin sets of radiometrically characterized panels: one set located in shadowed conditions, and the other set in non-shadowed conditions. By first determining the response values of the targets in non-shadowed areas, the procedure seeks, on the one hand, to obtain the same reflectance values in the homologous twin targets located in the shadowed areas after image de-shadowing, and on the other hand, to expand the correction model to all the detected shadows in the scene. The method includes three sequential blocks: (1) shadow detection; (2) shadow reduction; and (3) edge correction filtering to smooth the transition between de-shadowed areas and non-shadowed areas. 2. Materials and Methods 2.1. Experiment Conceptualization and Materials The main concept is that the spectral reflectance of a sample surface is the same in direct- plus-diffuse illumination conditions as in diffuse illumination conditions (what is variable is the spectral irradiance) [2,3]. Additionally, there is a need to account for atmospheric radiance contributions to the at-sensor spectral radiance [4,5]. These principles can be applied using an empirical line approach, among others [41–43]. Then, the information retrieval from shadowed pixels can be achieved in three stages: (1) shadow detection; (2) Remote Sens. 2021, 13, x FOR PEER REVIEWshadow reduction in all shadowed scene pixels using empirical relationships between 4 of 17 the characterized set of panels located in shadowed and non-shadowed areas; and (3) edge filtering to smooth the transition between de-shadowed and non-shadowed areas (Figure 1). Figure 1. (a) Basic concepts of the nature of shadows, including the three steps of the presented methodology: (1) shadow Figure 1. (a) Basic concepts of the nature of shadows, including the three steps of the presented methodology: (1) shadow detection; (2) shadow correction; (3) edge filtering. (b) General workflow for obtaining de-shadowed images from UAV detection; (2) shadow correction; (3) edge filtering. (b) General workflow for obtaining de-shadowed images from UAV orthomosaics with orthomosaics with shadowed shadowed areas. areas. The objective was to test whether the de-shadowing workflow was successful and operational in a real scenario. The method for checking the hypothesis is presented in an applied scenario (forested area). The case study was carried out in a Mediterranean for- Remote Sens. 2021, 13, 3808 4 of 16 The objective was to test whether the de-shadowing workflow was successful and operational in a real scenario. The method for checking the hypothesis is presented in an applied scenario (forested area). The case study was carried out in a Mediterranean forestry environment; specifically, within a restored quarry in Catalonia (X = 411,713 m; Y = 4,616,052 m (UTM31N-WGS84)) with grasslands, scrublands, and Pinus hapelensis trees of several heights due to the afforestation carried out in the frame of restoration actuations (Figure 2a). The authors previously used this area for UAV vegetation mapping studies, geometric and radiometric calibration experiments [38–40,44]. The flight was performed on 27 April 2018 between 10:36 and 10:47 UTC, at 90 m above ground level, designed to cover a 64,386 m2 study area (Table 1). By applying photogrammetric and Structure from Remote Sens. 2021, 13, x FOR PEER REVIEW 5 of 17 Motion processing, an orthoimage of 6 cm SR (Figure 2b) and a digital surface model (DSM) of 9 cm of pixel size were obtained. Figure 2. (a) Conceptualization of the applied scenario for checking the applicability of the hypothesis. (b) Experimental Figure 2. (a) Conceptualization of the applied scenario for checking the applicability of the hypothesis. (b) Experimental display of display of the the EVA EVApanels panelsused usedto obtain the the to obtain empirical line correction empirical between line correction shadowed between and non-shadowed shadowed areas (bottom and non-shadowed areas image: Parrot (bottom image:Sequoia false color Parrot Sequoia composite). false color composite). Flight plan TableIn1. order main features to relate and scene highlights. the non-shadowed illumination conditions and the shadowed conditions, we used a set of 500 mm × 500 mm × 20 mm ethylene-vinyl acetate (EVA) Date 27 April 2018 panels in direct illumination conditions and a twin set of panels in a cast shadow. Panels First photogram were spectrally capture characterized time a(UTC) using hand-held OceanOptics USB2000+ 10:36:33spectroradiom- Last photogram capture time (UTC) 10:46:23 eter [45] (Table 2, Figure 3a–e) and present a Lambertian behavior [37]. The spectral range Flight height Above Ground Level (m) 90 covered the VNIRofregions Number (340 (x4 photograms nm–1030 = bands)nm), with a sampling interval of 0.3 nm and a 623 (2492) spectral resolution (full width Planned along-track at half maximum, overlapping (%) FWHM) of 1.26 nm.90A reference panel was usedPlanned to measure the incoming across-track overlappingspectral (%) irradiance. At the same85time, UAV multi- Latitude spectral imageryof the wascenter of the whole acquired with scene a Parrot(◦ ) Sequoia sensor (Table 41◦ 410 31.29 00 N 3, Figure 3c,d) on ◦ ◦ 0 00 boardLongitude of the center a DJI Phantom of the whole III quadcopter forscene ( ) the◦ applied 1 49 43.18 case in a forestry E environment. Sun azimuth at central time of flight ( ) 146.50 Sun elevation at central time of flight (◦ ) 58.37 Table 2. Field spectroradiometer main features (OceanOptics USB2000+). FWHM (nm) Size (mm) Weight (g) Design Manufacturer/Model In order to relate the non-shadowed illumination conditions and the shadowed condi- OceanOptics tions,1.26 we used a89.1 set×of63.3 34.4 × 500 190 500×mm mm × 20 mm ethylene-vinyl acetate Czerny-Turner (EVA) panels USB2000+ in direct illumination conditions and a twin set of panels in a cast shadow. Panels were Input Focal Fiber optic Sampling in- Sensor CCD Grating #2 Spectral Length (mm) FOV (°) terval (nm) samples range (nm) 42 25 0.3 2048 340–1030 Remote Sens. 2021, 13, 3808 5 of 16 spectrally characterized using a hand-held OceanOptics USB2000+ spectroradiometer [45] (Table 2, Figure 3a–e) and present a Lambertian behavior [37]. The spectral range covered the VNIR regions (340 nm–1030 nm), with a sampling interval of 0.3 nm and a spectral resolution (full width at half maximum, FWHM) of 1.26 nm. A reference panel was used to measure the incoming spectral irradiance. At the same time, UAV multispectral imagery was acquired with a Parrot Sequoia sensor (Table 3, Figure 3c,d) on board a DJI Phantom III quadcopter for the applied case in a forestry environment. Table 2. Field spectroradiometer main features (OceanOptics USB2000+). FWHM (nm) Size (mm) Weight (g) Design Manufacturer/Model 1.26 89.1 × 63.3 × 34.4 190 Czerny-Turner OceanOptics USB2000+ Input Focal Length Sampling interval Grating #2 Spectral Fiber optic FOV (◦ ) Sensor CCD samples Remote Sens.(mm) 2021, 13, x FOR PEER REVIEW (nm) range (nm) 6 of 17 42 25 0.3 2048 340–1030 Figure 3. (a) (a)Field Fieldspectroradiometry spectroradiometrymeasurements. measurements.(b)(b) Twin Twinsetset of of targets located targets in diffuse located illumination in diffuse conditions. illumination (c) conditions. UAV sensor (c) UAV Parrot sensor Sequoia Parrot spectral Sequoia bandwidth. spectral (d) UAV bandwidth. (d) vehicle configuration UAV vehicle used for configuration the applied used scenarioscenario for the applied in a forested in a area. (e) area. forested Spectral signatures (e) Spectral of the of signatures EVAthe panels usedused EVA panels to obtain the the to obtain empirical line empirical correction line correctionbetween betweenshadowed shadowed and and non-shadowed areas. non-shadowed areas. Table 3. Parrot Sequoia camera main features. Expanded Dynamic Manufactur- Size (mm) Weight (g) Sensor Type Range (DN) er-Model 0–65,535 59 × 41 × 28 72 CMOS Parrot Sequoia Raw radiometric Field of View Input Focal Pixel size (µm) Sensor size (pixels) resolution (bits) (°) Length (mm) 10 47.5 3.98 3.75 × 3.75 1280 × 960 #1 Blue FWHM #2 Green #3 Red FWHM #4 Red-edge #5 NIR FWHM (nm) Remote Sens. 2021, 13, 3808 6 of 16 Table 3. Parrot Sequoia camera main features. Expanded Dynamic Size (mm) Weight (g) Sensor Type Manufacturer-Model Range (DN) 0–65,535 59 × 41 × 28 72 CMOS Parrot Sequoia Raw radiometric Input Focal Length Field of View (◦ ) Pixel size (µm) Sensor size (pixels) resolution (bits) (mm) 10 47.5 3.98 3.75 × 3.75 1280 × 960 #4 Red-edge FWHM #1 Blue FWHM (nm) #2 Green FWHM (nm) #3 Red FWHM (nm) #5 NIR FWHM (nm) (nm) No blue band 530–570 640–680 730–740 770–810 2.2. Methods: Shadow Detection, Shadow Reduction, Edge Filtering The UAV at-sensor reflectance image (obtained with the DLS irradiance reference) was empirically fitted to the characterized panels to obtain at-surface reflectance values [38]. Afterwards, the methodology for recovering information from shadowed areas consisted of three processing blocks, previously mentioned (Figure 1) and detailed in this section: (a) shadow detection, (b) shadow reduction, and (c) edge filtering. 2.2.1. Shadow Detection Let us consider the scenario of a clear sky and direct illumination, as in the example case, wherein the sun vector hits a given surface, thereby projecting shadows. Shadow de- tection can mainly be carried out with physical and/or image-based approaches (Figure 4). Physical approaches are based on the sun’s position and the cast shadow projected by the surface objects at a given location. The sun’s position is well known for its given geographical position and the surface is modeled using a DSM, geometrically obtaining a digital cast shadow model (DCSM) and a digital illumination model (DIM). Image-based approaches are based on the image band histogram, selecting those dark pixels located under a given threshold, assuming that those pixels have low illumination and are located in shadows. Both approaches provide good results but commission errors (false positives) for non-shadow pixels classified as shadows (low signal surfaces classified as shadows [23]) and omission errors (false negatives) of shadowed pixels classified as non-shadow, often appear. False positives in the physical approach appear when the DSM is not of enough quality or is not truly 3D but 2.5D, providing false positives specially in the trunk shadows in vegetation environments. False positives in the image-based approach appear in dark surfaces that are not located in shadows, such as waters or dark soils. To overpass the commission errors of both methods, in our approach, the detailed shadow detection was obtained combining image-based and physically based methods [33]. We selected pixels achieving the joined condition of being under the histogram threshold of the first quartile in the intensity image, in the Red Edge band and in the Near-Infrared band (image-based condition), also being under the shadow of the Digital Cast Shadow Model (DCSM) and in a self-shadow of the Digital Illumination Model (DIM) (physically based condition) in order to constrain the possible histogram commissions within the geometrically candidate pixels and considering that the Digital Elevation Models (DEM) are not truly 3D but 2.5D and provide false positives, especially in the trunk shadows. der the histogram threshold of the first quartile in the intensity image, in the Red Edge band and in the Near-Infrared band (image-based condition), also being under the shadow of the Digital Cast Shadow Model (DCSM) and in a self-shadow of the Digital Illumination Model (DIM) (physically based condition) in order to constrain the possible histogram commissions within the geometrically candidate pixels and considering that Remote Sens. 2021, 13, 3808 7 of 16 the Digital Elevation Models (DEM) are not truly 3D but 2.5D and provide false positives, especially in the trunk shadows. Figure Figure 4. 4. Shadow Shadow detection detection step step methodology methodology workflow. workflow. Two Twoproblems problemsmust mustbe besolved: solved: (1) (1) when when detecting detecting the the shadows shadows with aa physical with physical approach, approach, the the typical typical DEM DEM (2.5D) (2.5D) does does not not reproduce reproduce the the trees trees silhouette silhouette and and produces produces commission commission error error (left); (2) (left); (2)when whendetecting detecting shadows shadows with with the the image image histogram histogram approach, approach, dark dark objects objects located located inin non-shadowed non-shadowed areas areas also also produce commission errors (right). The solution is the masking of only those dark pixels located in projected produce commission errors (right). The solution is the masking of only those dark pixels located in projected shadow areas. shadow areas. 2.2.2. Shadow Reduction 2.2.2.Detected Shadow shadows Reductioncan be managed by using three approaches [21]: (1) simply mask- ing, (2) multisource Detected data can shadows fusion, be (3) radiometric managed enhancement. by using In the masking three approaches approach, [21]: (1) simply shadowed pixels are simply converted to values (removed), whereas masking, (2) multisource data fusion, (3) radiometric enhancement. In the masking NODATA in the multi- ap- source data fusion approach, the shadowed pixels are filled with imagery captured in other moments, both being approaches out of the focus of this work. Shadow reduction with radiometric enhancement techniques recovers the radiometry of the shadowed pixels, in- creasing their brightness [21] to make them seamless to neighboring non-shadowed pixels. Due to the great complexity and diversity of shadows, errors can appear in pixels with too high values after shadow reduction (overcorrection errors) or pixels with excessively high values after shadow reduction (undercorrection errors). In our approach (Figure 5), once the shadowed pixels were selected, they were further processed separately. The reduction was achieved with an invariant color model method [33], which takes advantage of pairs of panels located in shadowed and non-shadowed conditions simultaneously. We correlated four panels, whereas the remaining were used for validation purposes [38]. The linear function relationship between both panels was applied to all the detected shadowed pixels in the image. Finally, the de-shadowed pixels were covered over the image, resulting almost in the final result. However, at this stage, the problem of a bad transition between the de-shadowed areas and the non-shadowed areas appeared (Figure 5). To solve this problem located in the penumbra areas, a further step is needed that involved edge filtering (Section 2.2.3). conditions simultaneously. We correlated four panels, whereas the remaining were used for validation purposes [38]. The linear function relationship between both panels was applied to all the detected shadowed pixels in the image. Finally, the de-shadowed pixels were covered over the image, resulting almost in the final result. However, at this stage, the problem of a bad transition between the de-shadowed areas and the non-shadowed Remote Sens. 2021, 13, 3808 8 of 16 areas appeared (Figure 5). To solve this problem located in the penumbra areas, a further step is needed that involved edge filtering (Section 2.2.3). Figure 5. Shadow reduction block methodology: finding the empirical line function from the characterized targets located Figure 5. Shadow reduction block methodology: finding the empirical line function from the characterized targets located in shadowed and non-shadowed areas, and further applying it in all the detected shadowed areas. in shadowed and non-shadowed areas, and further applying it in all the detected shadowed areas. 2.2.3. Edge Filtering 2.2.3. Edge Filtering Boundary pixels between de-shadowed and non-shadowed regions are problematic due Boundary pixels to the spatial between mixture de-shadowedconditions of illumination and non-shadowed regionsperfection and the limited are problematic of the due to the spatial mixture of illumination conditions and the limited perfection DSM. Thereby, these pixels tend to be more or less overcorrected if considered as shadows, of the DSM. Thereby, these or undercorrected pixels tend otherwise, to be leading more orand to artifacts lessa overcorrected visual sensation if of considered as patched de- shadows, shadowedor undercorrected areas; the extractionotherwise, leading of a boundary map to artifacts [21,33] and a visual is a possible solutionsensation of for dealing patched with these de-shadowed areas. In ourareas; the extraction approach, of a boundary the shadow-mask was map [21,33]and vectorized is a apossible solu- buffer (edge tion belt)for wasdealing with for calculated these areas. that limitInarea; our in approach, this casethe shadow-mask study, the penumbrawas area vectorized was inand the aorder bufferof (edge 10 cm,belt) was so the calculated buffer area wasfor2that limit pixels area; at (1 pixel in each this case side study, the penumbra of the shadow mask). As we have pointed out in Section 1, precise determination of that area can be performed, according to [21], through Equation (1): ! 1 1 w=H −  (1) tan e − 2ε tan e + 2ε where w is the width of the penumbra, H is the height of the object (a tree, a house, etc.), e is the elevation angle of the center of the sun, and ε is the angular width of the sun at the top of the object. Figure 6 illustrates this formula. Nevertheless, due to our non-ultra-realistic object modelling goal, a visual inspection was enough to set a convenient value. That edge belt was used as follows: firstly, the output image of the shadow reduction block was smoothed with a 3 × 3 low-pass filter kernel; secondly, this smoothed image was clipped with the edge belt, and the smoothed belt was overlayed to the patched image, thus applying the kernel filter only in the penumbra transition (Figure 7). 𝜀 𝜀 tan 𝑒 tan 𝑒 2 2 where w is the width of the penumbra, H is the height of the object (a tree, a house, etc.), e is the elevation angle of the center of the sun, and ε is the angular width of the sun at the top of the object. Figure 6 illustrates this formula. Nevertheless, due to our Remote Sens. 2021, 13, 3808 9 of 16 non-ultra-realistic object modelling goal, a visual inspection was enough to set a con- venient value. Remote Sens. 2021, 13, x FOR PEER REVIEW 10 of 17 Figure 6. Theoretical Figure 6.estimation Theoreticalofestimation the penumbra width, of the assuming penumbra width,a flat terrain without assuming other without a flat terrain objects where other the ob-shadow is jects where projected (source: authorthe 1).shadow is projected (source: author 1). That edge belt was used as follows: firstly, the output image of the shadow reduc- tion block was smoothed with a 3 × 3 low-pass filter kernel; secondly, this smoothed image was clipped with the edge belt, and the smoothed belt was overlayed to the patched image, thus applying the kernel filter only in the penumbra transition (Figure 7). Figure 7. Edge filtering step methodology workflow: the cover problem is a patching pattern of de-shadowed pixels Figure 7. Edge filtering step methodology workflow: the cover problem is a patching pattern of de-shadowed pixels over overoriginal the the original image. image. To solve To solve this,this, the the shadow shadow maskmask is vectorized, is vectorized, a two-pixel a two-pixel bufferisisapplied, buffer applied,and andused used to to clip clip a a smoothed raster of the patched image. This clip is overlayed again on the patched image to smooth the transition smoothed raster of the patched image. This clip is overlayed again on the patched image to smooth the transition between between the de-shadowed and the non-shadowed pixels. 2.3. Methods: 2.3. Methods: Statistical Statistical Data Data Analysis Analysis We will present the coefficients of determination, R2 , of the linear regression between We will present the coefficients of determination, R², of the linear regression be- the seven twin targets located both in diffuse illumination conditions and in direct illumina- tween the seven twin targets located both in diffuse illumination conditions and in direct tion conditions. If these regressions show high R2 (>0.90) and p-values lower than 0.01 in all illumination conditions. If these regressions show high R² (>0.90) and p-values lower than the spectral bands of the sensor, with a reasonably good distribution along the [0, 100] % of 0.01 in all the spectral bands of the sensor, with a reasonably good distribution along the reflectance values (also for all spectral bands), we will consider that the approach is robust [0, 100] % of reflectance values (also for all spectral bands), we will consider that the ap- and fulfils the hypothesis of our proposal. This analysis will be completed with a visual proach is robust and fulfils the hypothesis of our proposal. This analysis will be com- inspection of the corrected image to check whether most shadows have been removed. pleted with a visual inspection of the corrected image to check whether most shadows have been removed. 3. Results The shadow detection method provided good results in detecting those dark pixels located in projected shadows. The sun position was well known during the drone mis- Remote Sens. 2021, 13, 3808 10 of 16 3. Results The shadow detection method provided good results in detecting those dark pixels located in projected shadows. The sun position was well known during the drone mission (±5 min the central time of the flight), and the DSM was detailed enough (9 cm/pixel) to obtain a good modelling of the projected shadows (DCSM). However, the application Remote Sens. 2021, 13, x FOR PEER REVIEW of 11 of 17 a constraint to select only dark pixels in the projected shadows avoided false positives, particularly in the trunk shadows (Figure 8). Figure 8. (a) Original orthoimage with shadows. (b) Shadow mask obtained after the shadow detection procedure. (c) Detail Figure 8. (a) Original orthoimage with shadows. (b) Shadow mask obtained after the shadow detection procedure. (c) of some tree shadows. (d) Detail of the digital cast shadow model projection, based on a 2.5D digital surface model. (e) Detail Detail of some tree shadows. (d) Detail of the digital cast shadow model projection, based on a 2.5D digital surface model. of (e)the constraint Detail for selecting of the constraint foronly dark only selecting pixelsdark in the projected pixels in the shadows projectedthat avoided shadows thatfalse positives. avoided false positives. De-shadow reduction method, which is based on the empirical relationship between De-shadow reduction method, which is based on the empirical relationship between the shadowed panels and the non-shadowed panels, providing the ad hoc parameters for the shadowed panels and the non-shadowed panels, providing the ad hoc parameters for converting the shadowed pixels of the image to de-shadowed pixels, in each of the Parrot converting the shadowed pixels of the image to de-shadowed pixels, in each of the Parrot Sequoia bands (Green, Red, Red Edge and Near-Infrared) (Figure 9, Table 4). Sequoia bands (Green, Red, Red Edge and Near-Infrared) (Figure 9, Table 4). Remote RemoteSens. Sens.2021, 2021,13, 13,x3808 FOR PEER REVIEW 12 11of of17 16 Figure 9. Spectral-dependent linear regression plots for Parrot Sequoia bands, between twin targets (n = 7) located in diffuse Figure 9. Spectral-dependent linear regression plots for Parrot Sequoia bands, between twin targets (n = 7) located in illumination diffuse conditions illumination and direct conditions andillumination conditions.conditions. direct illumination All regressions show an R2 show All regressions higheranthan R² 0.90 a p-value andthan higher < 0.01; 0.90 and a additionally, p-value < 0.01;note the reasonably additionally, good note the distribution reasonably good along the [0, 100] distribution % of along reflectance the [0, 100] %values (also forvalues of reflectance all spectral (also bands). for all These empirical spectral line functions bands). These empiricalwere lineapplied functionsto were all theapplied detectedtoshadowed pixelsshadowed all the detected to enhancepixels their radiometry to enhance and theirobtain radi- de-shadowed ometry pixels. and obtain de-shadowed pixels. Table Table 4.4. Linear Linearregressions regressions to fittoshadowed fit shadowed reflectance reflectance panels topanels to the ofreflectance the reflectance of the the non-shadowed non-shadowed panels for the Parrot Sequoia panels for the Parrot Sequoia camera (case study). camera (case study). Parrot ParrotSequoia Bands Sequoia Bands R2 R2 Slope Slope Bias Bias Green: 530–570 (nm) 0.969 2.258 12.152 Green: 530–570 (nm) 0.969 2.258 12.152 Red: Red: 640–680 (nm) 640–680 (nm) 0.946 0.946 1.945 1.945 11.050 11.050 Red RedEdge: Edge: 730–740 (nm) 730–740 (nm) 0.983 0.983 1.533 1.533 10.771 10.771 NIR:770–810 NIR: 770–810 (nm) (nm) 0.990 0.990 1.479 1.479 3.642 3.642 The Theapplication applicationof ofthis this linear linear function function toto the the previously previously detected detected shadowed shadowed pixels pixels led led to the recovery of information that was previously hidden, thus obtaining at-surface to the recovery of information that was previously hidden, thus obtaining at-surface reflectance reflectance values both in values both indirect directillumination illuminationand anddiffuse diffuse illumination illumination conditions. conditions. How- However, ever, and asand as mentioned mentioned above, above, edgeedge filtering filtering waswas required required to smooth to smooth thethe patched patched transi- transition tion between between the de-shadowed the de-shadowed andand the shadowed the shadowed pixels pixels (Figure (Figure 10). 10). Remote Sens. 2021, 13,Remote 3808 Sens. 2021, 13, x FOR PEER REVIEW 12 of 16 13 Figure 10.Figure of Results Results10. of the the shadow shadowand detection detection shadowand shadow reduction reduction stages. This isstages. Thisresult the partial is thebefore partial result edge filtering. before edge filtering. The edge map, which is based on the mask shadow, delimited the zone with o The edge map, whichand corrections is based on the maskdue undercorrections shadow, to the delimited penumbra,theandzone with of a buffer over- two pixels corrections andapplied undercorrections duetransition to delimit the to the penumbra, zone. This and a buffer buffer of two was used pixels to clip the was smoothed im applied to delimit andthe transition covered the zone. This buffer de-shadowed wasthereby image, used to obtaining clip the smoothed reasonablyimage goodand results (Fig covered the de-shadowed 11). image, thereby obtaining reasonably good results (Figure 11). Remote Sens. 2021, 13,Remote 3808 Sens. 2021, 13, x FOR PEER REVIEW 13 of 16 14 Figure 11. Results of the radiometric correction of UAV imagery using the empirical approach Figure 11. Results of the radiometric correction of UAV imagery using the empirical approach on shadow reduction. (a) on shadow reduction. (a) At-sensor Parrot Sequoia UAV image (obtained from DLS radiometric At-sensor Parrot Sequoia UAV image (obtained from DLS radiometric correction) before de-shadowing. (b) At-surface UAV image correction) before empirically fittedde-shadowing. (b) At-surface to the characterized UAV panels and image empirically empirically de-shadowed. fitted (c)to the characterized Detail before de-shadowing. panels and empirically de-shadowed. (c) Detail before de-shadowing. (d) Detail of shadow detection (d) Detail of shadow detection results. (e) Detail of edge filtering results. (f) Detail of the image after de-shadowing. results. (e) Detail of edge filtering results. (f) Detail of the image after de-shadowing. 4. Discussion 4. Discussion The ultra-high spatial resolution of UAV imagery allows the use of twin sets of The ultra-high fordablespatial resolution and useful of in targets UAV fieldimagery campaigns. allows Images thewith use less of twin than sets of allow 10 cm/pix affordable and enough useful targets in field campaigns. Images with less than 10 cm/pix pixels inside panels of 500 mm × 500 mm for obtaining a central pixel with allow for enough pixels inside panels adjacency of 500 mm contributions. This× new 500 mm for obtaining opportunity a central for remote pixel opens sensing withoutdoors for e adjacency contributions. This new opportunity for remote sensing opens use of the invariant color model method [33] using twin panels simultaneouslydoors for easy locate use of the invariant color model shadowed method [33] using and non-shadowed twin panels conditions, simultaneously and ultra-high located insensors. spatial resolution shadowed and non-shadowed conditions, and ultra-high spatial resolution sensors. The detailed imagery and derived DSM allowed a reasonable modelling of the v The detailed imagery etation and derived structures, DSM but the allowed nature a reasonable of the 2.5D DSMmodelling of the vegeta- makes it necessary to combine tion structures, but the nature physical of the approach 2.5D with anDSM makes it approach image-based necessaryin to order combine the physical to obtain better shadow approach with an image-based tection. However, approach in order to not all shadows obtain were bettermainly detected, shadowdue detection. How- comple to the specific ever, not all shadows of canopy were detected, structures andmainly due to the specific the photogrammetric surface complexity of canopy modelling inaccuracies. structures and the photogrammetric Shadow correction surface basedmodelling on the linearinaccuracies. relationship between the targets located Shadow correction der directbased on the linear illumination relationship and diffuse between illumination ledthe targets to good located results under in the extrapolatio direct illumination and all the diffuseshadows detected illumination on theled to good scene, where results in the general extrapolation overcorrections orto all undercorrecti the detected shadows are noton the scene, observed. Onwhere general the other hand,overcorrections there were artifactsor undercorrections in the limit areasareof the detec not observed. On shadows. Thishand, the other problemtherewas tackled were withinthe artifacts theedge limitfiltering areas ofcorrection. the detected The edge m shadows. This problem was tackled with the edge filtering correction. The edge map based on the shadow mask was precise for detecting the artifacts, whereas the buffer accounted Remote Sens. 2021, 13, 3808 14 of 16 for all of them. The kernel filter, which smoothed the penumbra area, was effective for reducing the patching effects, but obviously blurred these areas. Limitations: This method is designed for ad hoc radiometric correction, as the approach is completely empirical. Therefore, the empirical line parameters are scene- dependent. The scale factor depends on the opacity of objects that causes the shadow, and the bias factor depends on the atmospheric composition (atmospheric optical depth). This is a limitation to establish generalized, fixed parameters, to relate shadowed pixels and non-shadowed pixels. However, the method is also useful for common UAV applications where the area covered is not large in scope for containing great atmospheric composition differences; this is usually reasonable when the flight height does not exceed 120 m (400 ft), as in current regulations in most European countries [46]. Although the overall enhancement of the image is evident, it is worth noting that not absolutely all of the shadows were recovered. Due to the ultra-high spatial detail of UAV imagery, some shadows remained inside the tree canopy. Additionally, not all the shadows were enhanced with the same quality because the empirical approach was fitted to the shadows of the tree canopy, but not in a closed forest situation; therefore, we followed a conservative trend to avoid overcorrection errors. The de-shadowed areas were not patched over the non-shadowed areas, and the transition between both illumination conditions is, in our humble opinion, a strength of the presented methodology. The shadow detection, shadow reduction and edge filtering workflow use GIS tools implemented in most GIS software, which can thus be automatized by using model builders. In our case, we used the free MiraMon RS and GIS software [47]. 5. Conclusions In this study, a method for shadow reduction in UAV imagery was described and tested; the method used well-known RS and image enhancement techniques from the airborne and satellite image processing world combined with the new opportunities of drone imagery in reduced scenarios. The empirical line correction using radiometric references has proven to be a good method for obtaining at-surface reflectance values from UAV imagery. In this work, it has also proven to be an interesting method for obtaining surface reflectance values in shadowed forested areas. Apart from the efficiency of the method for recovering radiometric information, it is a simple and operational method for smoothing the transition in the edge zone by applying a correction in the buffer of the limit area. The high spatial resolution of UAV imagery allows the use twin sets of affordable reference targets in field campaigns, which is a common requirement and limitation of ad hoc radiometric corrections. This new opportunity for RS opens the door to the easy use of invariant color models, as well as other well-known remote sensing techniques such as shadow detection or edge filtering smoothing. Author Contributions: Conceptualization, X.P. and J.-C.P.; methodology, X.P. and J.-C.P.; software, X.P.; validation, X.P. and J.-C.P.; formal analysis, J.-C.P.; investigation, J.-C.P.; resources, X.P.; data curation, J.-C.P.; writing—original draft preparation, J.-C.P.; writing—review and editing, X.P.; visualization, J.-C.P.; supervision, X.P.; project administration, X.P.; funding acquisition, X.P. All authors have read and agreed to the published version of the manuscript. Funding: This work was partially funded by the Catalan Government under Grant (SGR2017-1690), by the European Union through the NewLife4Drylands Project (LIFE20 PRE/IT/000007) and by the Spanish MCIU Ministry through the NEWFORLAND research project (RTI2018-099397-B-C21/C22 MCIU/AEI/ERDF, EU). Institutional Review Board Statement: Not applicable. Informed Consent Statement: Not applicable. Conflicts of Interest: The authors declare no conflict of interest. Remote Sens. 2021, 13, 3808 15 of 16 References 1. Adler-Golden, S.; Matthew, M.W.; Anderson, G.P.; Felde, G.W.; Gardner, J.A. An Algorithm for De-Shadowing Spectral Imagery. Proc. SPIE Imaging Spectrom. 2002, VIII 4816, 203–210. [CrossRef] 2. Nicodemus, F.E.; Richmond, J.C.; Hsia, J.J. Geometrical Considerations and Nomenclature for Reflectance; National Bureau of Standards, US Department of Commerce: Washington, DC, USA, 1977. Available online: http://graphics.stanford.edu/courses/cs448-05- winter/papers/nicodemus-brdf-nist.pdf (accessed on 1 August 2021). 3. Milton, E.J.; Schaepmann, M.E.; Anderson, K.; Kneubühler, M.; Fox, N. Progress in field spectroscopy. Remote Sens. Environ. 2009, 113, S92–S109. [CrossRef] 4. Teillet, P.M.; Guindon, B.; Goodenough, D.G. On the slope-aspect correction of multispectral scanner data. Can. J. Remote Sens. 1982, 8, 84–106. [CrossRef] 5. Pons, X.; Solé-Sugrañes, L. A simple radiometric correction model to improve automatic mapping of vegetation from multispectral satellite data. Remote Sens. Environ. 1994, 45, 317–332. [CrossRef] 6. Pons, X.; Pesquer, L.; Cristóbal, J.; González-Guerrero, O. Automatic and improved radiometric correction of Landsat imagery using reference values from MODIS surface reflectance images. Int. J. Appl. Earth Obs. Geoinf. 2014, 33, 243–254. [CrossRef] 7. Riaño, D.; Chuvieco, E.; Salas, J.; Aguado, I. Assessment of different topographic corrections in Landsat-TM data for mapping vegetation types. IEEE Trans. Geosci. Remote Sens. 2003, 41, 1056–1061. [CrossRef] 8. Richter, R.; Kellenberger, T.; Kaufmann, H. Comparison of topographic correction methods. Remote Sens. 2009, 1, 184–196. [CrossRef] 9. Zhu, Z.; Woodcock, C.E. Object-based cloud and cloud shadow detection in Landsat imagery. Remote Sens. Environ. 2012, 118, 83–94. [CrossRef] 10. Zhu, Z.; Woodcock, C.E. Automated cloud, cloud shadow, and snow detection in multitemporal Landsat data: An algorithm designed specifically for monitoring land cover change. Remote Sens. Environ. 2014, 152, 217–234. [CrossRef] 11. Aeronautics and Space Administration (NASA). MODIS Web. Available online: https://modis.gsfc.nasa.gov/ (accessed on 1 August 2021). 12. National Oceanic and Atmospheric Administration (NOAA). Advanced Very High Resolution Radiometer—AVHRR. Available online: https://www.avl.class.noaa.gov/release/data_available/avhrr/index.htm (accessed on 1 August 2021). 13. Aeronautics and Space Administration (NASA). Landsat Data Continuity Mission (LDCM). Available online: https://www.nasa. gov/mission_pages/landsat/main/index.html (accessed on 1 August 2021). 14. European Space Agency (ESAa). ESA Sentinel Online. Sentinel-2 Mission. Available online: http://www.esa.int/Our_Activities/ Observing_the_Earth/Copernicus/Sentinel-2 (accessed on 1 August 2021). 15. DigitalGlobe. Satellite Information. Available online: https://www.digitalglobe.com/resources/satellite-information (accessed on 1 August 2021). 16. Fernández-Renau, A.; Gómez, J.A.; de Miguel, E. The INTA AHS System. Proc. SPIE 2005. [CrossRef] 17. Jiménez, M.; Díaz-Delgado, R.; Vaughan, P.; De Santis, A.; Fernández-Renau, A.; Prado, E.; Gutiérrez de la Cámara, O. Airborne hyperspectral scanner (AHS) mapping capacity simulation for the Doñana biological reserve scrublands. In Proceedings of the 10th International Symposium on Physical Measurements and Signatures in Remote Sensing, Davos, Switzerland, 12–14 March 2007; Schaepman, M., Liang, S., Groot, N., Kneubühler, M., Eds. Available online: http://www.isprs.org/proceedings/XXXVI/7- C50/papers/P81.pdf (accessed on 28 October 2017). 18. Itres Research Ltd. Compact Airborne Spectrographic Imager (CASI-1500H). Available online: https://itres.com/wp-content/ uploads/2019/09/CASI1500.pdf (accessed on 1 August 2021). 19. Leica Geosystems. Leica DMCII Airborne Digital Camera. Available online: https://leica-geosystems.com/products/airborne- systems/imaging-sensors/leica-dmciii (accessed on 1 August 2021). 20. Chen, Y.; Wen, D.; Jing, L.; Shi, P. Shadow information recovery in urban areas from very high resolution satellite imagery. Int. J. Remote Sens. 2007, 28, 3249–3254. [CrossRef] 21. Dare, P.M. Shadow analysis in high-resolution satellite imagery of urban areas. PE&RS 2005, 71, 169–177. [CrossRef] 22. Richter, R.; Muller, A. De-shadowing of satellite/airborne imagery. Int. J. Remote Sens. 2005, 26, 3137–3148. [CrossRef] 23. Arévalo, V.; González, J.; Ambrosio, G. Shadow detection in colour high-resolution satellite images. Int. J. Remote Sens. 2008, 29, 1945–1963. [CrossRef] 24. MicaSense. MicaSense RedEdge™ 3 Multispectral Camera User Manual; MicaSense, Inc.: Seattle, WA, USA, 2015; p. 33. Available online: https://support.micasense.com/hc/en-us/article_attachments/204648307/RedEdge_User_Manual_06.pdf (accessed on 1 August 2021). 25. Parrot Drones. Parrot Sequoia Technical Specifications. 2018. Available online: https://www.parrot.com/global/parrot- professional/parrot-sequoia#technicals (accessed on 1 August 2021). 26. Mapir. Specifications. 2021. Available online: https://www.mapir.camera/pages/survey3-cameras#specs (accessed on 1 August 2021). 27. SAL Engineering. MAIA—The Multispectral Camera. 2021. Available online: https://www.salengineering.it/public/en/p/ maia.asp (accessed on 1 August 2021). 28. Markelin, L.; Simis, S.G.H.; Hunter, P.D.; Spyrakos, E.; Tyler, A.N.; Clewley, D.; Groom, S. Atmospheric Correction Performance of Hyperspectral Airborne Imagery over a Small Eutrophic Lake under Changing Cloud Cover. Remote Sens. 2017, 9, 2. [CrossRef] Remote Sens. 2021, 13, 3808 16 of 16 29. Salvador, E.; Cavallaro, A.; Ebrahimi, T. Shadow identification and classification using invariant color models. In Proceedings of the 2001 IEEE International Conference on Acoustics, Speech, and Signal Processing, Salt Lake City, UT, USA, 7–11 May 2001; Volume 3, pp. 1545–1548. [CrossRef] 30. Shahtahmassebi, A.; Yang, N.; Wang, K.; Moore, N.; Shen, Z. Review of shadow detection and de-shadowing methods in remote sensing. Chinese Geogr. Sci. 2013, 23, 403–420. [CrossRef] 31. Wu, T.; Zhang, L.; Huang, C. An analysis of shadow effects on spectral vegetation indices using a ground-based imaging spectrometer. In Proceedings of the 7th Workshop on Hyperspectral Image and Signal Processing: Evolution in Remote Sensing (WHISPERS), Tokyo, Japan, 2–5 June 2015; pp. 1–4. [CrossRef] 32. Schaepman-Strub, G.; Schaepman, M.E.; Painter, T.H.; Dangel, S.; Martonchik, J.V. Reflectance quantities in optical remote sensing—definitions and case studies. Remote Sens. Environ. 2006, 103, 27–42. [CrossRef] 33. Adeline, K.R.M.; Chen, M.; Briottet, X.; Pang, S.K.; Paparoditis, N. Shadow detection in very high spatial resolution aerial images: A comparative study. ISPRS J. Photogramm. 2013, 80, 21–38. [CrossRef] 34. Honkavaara, E.; Saari, H.; Kaivosoja, J.; Pölönen, I.; Hakala, T.; Litkey, P.; Mäkynen, J.; Pesonen, L. Processing and Assessment of Spectrometric, Stereoscopic Imagery Collected Using a Lightweight UAV Spectral Camera for Precision Agriculture. Remote Sens. 2013, 5, 5006–5039. [CrossRef] 35. Wang, C.; Myint, S.W. A Simplified Empirical Line Method of Radiometric Calibration for Small Unmanned Aircraft Systems- Based Remote Sensing. IEEE JSTARS 2015, 8, 1876–1885. [CrossRef] 36. Padró, J.C.; Pons, X.; Aragonés, D.; Díaz-Delgado, R.; García, D.; Bustamante, J.; Pesquer, L.; Domingo-Marimon, C.; González- Guerrero, O.; Cristóbal, J.; et al. Radiometric Correction of Simultaneously Acquired Landsat-7/Landsat-8 and Sentinel-2A Imagery Using Pseudoinvariant Areas (PIA): Contributing to the Landsat Time Series Legacy. Remote Sens. 2017, 9, 1319. [CrossRef] 37. Padró, J.C.; Muñoz, F.J.; Avila, L.A.; Pesquer, L.; Pons, X. Radiometric Correction of Landsat-8 and Sentinel-2A Scenes Using Drone Imagery in Synergy with Field Spectroradiometry. Remote Sens. 2018, 10, 1687. [CrossRef] 38. Padró, J.C.; Carabassa, V.; Balagué, J.; Brotons, L.; Alcañiz, J.M.; Pons, X. Monitoring opencast mine restorations using Unmanned Aerial System (UAS) imagery. Sci. Total Environ. 2019, 657, 1602–1614. [CrossRef] [PubMed] 39. Pons, X.; Padró, J.C. An Empirical Approach on Shadow Reduction of UAV Imagery in Forests. In Proceedings of the IGARSS 2019–2019 IEEE International Geoscience and Remote Sensing Symposium, Yokohama, Japan, 28 July–2 August 2019; pp. 2463–2466. [CrossRef] 40. Bareth, G.; Aasen, H.; Bendig, J.; Gnyp, M.; Bolten, A.; Jung, A.; Michels, R.; Soukkamäki, J. Low-weight and UAV-based Hyperspectral Full-frame Cameras for Monitoring Crops: Spectral Comparison with Portable Spectroradiometer. PFG 2015, 1, 69–79. [CrossRef] 41. Smith, G.M.; Milton, E.J. The use of the empirical line method to calibrate remotely sensed data to reflectance. Int. J. Remote Sens. 1999, 20, 2653–2662. [CrossRef] 42. Perry, E.M.; Warner, T.; Foote, P. Comparison of atmospheric modelling versus empirical line fitting for mosaicking HYDICE imagery. Int. J. Remote Sens. 2000, 21, 799–803. [CrossRef] 43. Karpouzli, E.; Malthus, T. The empirical line method for the atmospheric correction of IKONOS imagery. Int. J. Remote Sens. 2003, 24, 1143–1150. [CrossRef] 44. Padró, J.C.; Muñoz, F.J.; Planas, J.; Pons, X. Comparison of four UAV georeferencing methods for environmental monitoring purposes focusing on the combined use with airborne and satellite remote sensing platforms. IJAEO 2019, 79, 130–140. [CrossRef] 45. Photonic Solutions. USB200+ Data Sheet. OceanOptics; Ocean Optics, Inc.: Dunedin, FL, USA, 2006. Available online: https: //www.photonicsolutions.co.uk/upfiles/USB2000PlusSpectrometerDatasheetLG14Jun18.pdf (accessed on 1 August 2021). 46. European Union Aviation Safety Agency (EASA). Drones—Regulatory Framework Background. 2019. Available online: https: //www.easa.europa.eu/easa-and-you/civil-drones-rpas/drones-regulatory-framework-background (accessed on 1 August 2021). 47. Pons, X. MiraMon. Geographic Information System and Remote Sensing Software; Version 8.1j; Centre de Recerca Ecològica i Aplicacions Forestals, CREAF: Bellaterra, Spain, 2019. Available online: https://www.miramon.cat/Index_usa.htm (accessed on 1 August 2021)ISBN 84-931323-4-9.