Radiometric Correction of Landsat-8 and Sentinel-2A Scenes Using Drone Imagery in Synergy with Field Spectroradiometry

Abstract

The main objective of this research is to apply unmanned aerial system (UAS) data in synergy with field spectroradiometry for the accurate radiometric correction of Landsat-8 (L8) and Sentinel-2 (S2) imagery. The central hypothesis is that imagery acquired with multispectral UAS sensors that are well calibrated with highly accurate field measurements can fill in the scale gap between satellite imagery and conventional in situ measurements; this can be possible by sampling a larger area, including difficult-to-access land covers, in less time while simultaneously providing good radiometric quality. With this aim and by using near-coincident L8 and S2 imagery, we applied an upscaling workflow, whereby: (a) UAS-acquired multispectral data was empirically fitted to the reflectance of field measurements, with an extensive set of radiometric references distributed across the spectral domain; (b) drone data was resampled to satellite grids for comparison with the radiometrically corrected L8 and S2 official products (6S-LaSRC and Sen2Cor-SNAP, respectively) and the CorRad-MiraMon algorithm using pseudo-invariant areas, such as reflectance references (PIA-MiraMon), to examine their overall accuracy; (c) then, a subset of UAS data was used as reflectance references, in combination with the CorRad-MiraMon algorithm (UAS-MiraMon), to radiometrically correct the matching bands of UAS, L8, and S2; and (d) radiometrically corrected L8 and S2 scenes obtained with UAS-MiraMon were intercompared (intersensor coherence). In the first upscaling step, the results showed a good correlation between the field spectroradiometric measurements and the drone data in all evaluated bands (R² > 0.946). In the second upscaling step, drone data indicated good agreement (estimated from root mean square error, RMSE) with the satellite official products in visible (VIS) bands (RMSE_VIS < 2.484%), but yielded poor results in the near-infrared (NIR) band (RMSE_NIR > 6.688% was not very good due to spectral sensor response differences). In the third step, UAS-MiraMon indicated better agreement (RMSE_VIS < 2.018%) than the other satellite radiometric correction methods in visible bands (6S-LaSRC (RMSE < 2.680%), Sen2Cor-SNAP (RMSE < 2.192%), and PIA-MiraMon (RMSE < 3.130%), but did not achieve sufficient results in the NIR band (RMSE_NIR < 7.530%); this also occurred with all other methods. In the intercomparison step, the UAS-MiraMon method achieved an excellent intersensor (L8-S2) coherence (RMSE_VIS < 1%). The UAS-sampled area involved 51 L8 (30 m) pixels, 143 S2 (20 m) pixels, and 517 S2 (10 m) pixels. The drone time needed to cover this area was only 10 min, including areas that were difficult to access. The systematic sampling of the study area was achieved with a pixel size of 6 cm, and the raster nature of the sampling allowed for an easy but rigorous resampling of UAS data to the different satellite grids. These advances improve human capacities for conventional field spectroradiometry samplings. However, our study also shows that field spectroradiometry is the backbone that supports the full upscaling workflow. In conclusion, the synergy between field spectroradiometry, UAS sensors, and Landsat-like satellite data can be a useful tool for accurate radiometric corrections used in local environmental studies or the monitoring of protected areas around the world.

1. Introduction and Objectives

Knowledge of the Earth’s surface is essential for properly planning and managing resources and making concrete progress in effective and durable sustainability, in which people and other living beings can develop in optimal ways with the necessary means [1]. Remote sensing (RS) plays a key role in this goal by contributing with a vast amount of data towards environmental monitoring [2], not only at global scales, but also regional and local scales [3]. International associations and governments have been creating an increasingly efficient Earth observation (EO) network in order to guarantee access to high-quality data systems [4]. For these purposes, an increasing temporal availability of EO resources is arising from political will and collaboration.

To maximize efforts in terms of scientific, social, and economic benefits, spatial agencies advance together through initiatives as the Global Earth Observation System of Systems (GEOSS) [5]. In GEOSS, the National Aeronautics and Space Administration (NASA) and the European Space Agency (ESA), among others, collaborate in the acquisition of high spatial resolution (SR) multispectral imagery, with the Landsat Data Continuity Mission (LDCM) [6] and the Copernicus Sentinel-2 (S2) [7] programs, respectively. The synergistic use of different satellite missions implies that images must have similar characteristics (spatial, spectral, radiometric, etc.). In fact, the design of several missions in recent years has gone in this direction, such as the aforementioned LCDM and Copernicus, currently with two Landsat (Landsat-7 (L7) and Landsat-8 (L8)) and two Sentinel-2 (Sentinel-2A (S2A) and Sentinel-2B (S2B)) satellites in orbit. The Sentinel constellation was designed to support the Landsat Program and provide continuity to the SPOT program [8] by embedding sensors with similar spatial and spectral resolution features. The temporal orbit cycle of Sentinel-2 platforms was specifically thought to maximize combined consecutive observations with Landsat [9], thereby increasing their joint monitoring possibilities.

Nevertheless, these efforts are still considered insufficient as long as the sensed raw data cannot be ensured to receive radiometric correction that allows the combined use of the different platforms and sensors intended to act in a coordinated way. The backbone of the radiometric correction of satellite data is the proper modeling of the atmospheric parameters affecting each captured image [10]. Perhaps the most interesting physical magnitude to be obtained after an atmospheric correction of remote sensing optical imagery is spectral surface reflectance (ρ_sλ; dimensionless), which can be defined as the proportion of reflected radiant flux over the incident radiant flux in a given wavelength [11]. Indeed, surface reflectance is an essential variable for characterizing main land cover types with maximum accuracy. Besides atmospheric correction, topographic correction is also interesting for properly managing terrain incident angles, shadows, and so forth. [12]. Therefore, to this end, it is implemented in several radiometric correction processors (e.g., Sen2Cor-SNAP, CorRad-MiraMon) using a digital elevation model (DEM) as auxiliary data. The DEM can also be used for other purposes, such as modeling the transmittance as a function of the terrain height, assuming a higher atmospheric transmittance at higher altitudes (e.g., 6S-LaSRC).

Official agencies provide surface reflectance products of Landsat and Sentinel-2 data by using specifically designed radiometric correction processors. For L8, the United States Geological Survey (USGS) uses an adaption of the Second Simulation of the Satellite Signal in the Solar Spectrum (6S) [13], namely Landsat Surface Reflectance Code (6S-LaSRC) [14]. For S2, ESA uses an adaption of the Atmospheric Correction (ATCOR3) [15] processor, also implemented in Sentinel’s Application Platform (Sen2Cor-SNAP) [16,17]. Additionally, both official agencies work together towards harmonizing surface reflectance products to support the creation of time series combining L8 and S2 data (Harmonized Landsat and Sentinel (HLS)) [18]. Other works in this field have analyzed the geometric adaption of both grid origins and different spatial resolution grids [19] and the combined use of L8 and S2 data [20]. In recent years, international projects have been focusing on the radiometric coherence between sensors, and particularly between L8 and S2 products, developing and comparing atmospheric correction processors (Atmospheric Correction Intercomparison Exercise (ACIX) [21,22]). Thanks to these efforts, radiometric correction methods today provide an acceptable level of coherence, as can be validated by field spectroradiometric data captured almost simultaneously with satellite data [23]. Reaching a higher level of precision, a typical Landsat-like optical sensor measures spectral radiance (L_sλ; W m⁻² sr⁻² µm⁻¹), including both the land surface-reflected radiance, which is our main interest, and the atmospheric spectral radiance (L_atmλ), composed by the upwelling spectral radiance and the reflection of the downwelling spectral radiance of the atmosphere. Atmospheric spectral total optical depth (τ_0λ; dimensionless) weakens the downwelling solar spectral irradiance (E_0λ; W m⁻² µm⁻¹) and the upwelling surface-reflected radiance [24].

For years, satellite remote sensing has been attempting to relate the images with field spectroscopy data [25]. Field spectroradiometric measurement of surface reflectance at the satellite overpass is practically unaffected by atmospheric effects. Moreover, it is a valuable source of information for assessing and reinforcing the accuracy of radiometric correction [26,27]. Currently, in situ measurements are carried out by operators using portable field spectroradiometers, sometimes with elevator devices [28]. However, they face various challenges, such as trying to access high tree canopies or water surfaces, or the impossibility of measuring in dangerous or restricted areas. Another important restriction of conventional field spectroradiometry is the specific covered area and the sampling method themselves, since during the acceptable time of 20 min before and after the satellite overpass, a skilled operator can scan up to a reasonable area equivalent to 20 pixels of OLI 30 m, 45 pixels of MSI 20 m, or 180 pixels of MSI 10 m. Moreover, it is advisable to adapt sampling strategies to a predetermined satellite pixel size and origin, followed by a resampling of the pointer field data to the different satellite pixel sizes, which is neither a trivial nor easy task by any means. For example, although ground-sampling data is usually georeferenced in an accurate way, in some situations, such as the forest understory, it can be difficult to achieve accurate georeferencing [29].

In this context, imagery acquired with hyperspectral or multispectral sensors onboard unmanned aerial systems (UAS) is expected to complement field data, as the UAS platform can access those hardly accessible areas in a faster way than human operators and systematically sample the surface with higher spatial accuracy. The use of ultrahigh spatial resolution data acquired with UAS for quantitative remote sensing spectroscopy purposes has evolved enormously over the years, as comprehensively reviewed in [30]. In 2014, an international initiative (OPTIMISE) [31] was engaged in exploring new opportunities that UAS platforms and sensors could potentially provide to near-ground Earth observations. MacArthur and Robinson [32] presented a review of field spectroscopy and advanced the impact of new technologies, such as UAS, particularly in supporting the empirical calibration and validation of satellite observations. At the same time, they emphasize the importance of accurately processing UAS imagery to ensure its reliability.

Other works propose the use of reference panels in synergy with field spectroscopy data and UAS-acquired imagery; for example, using a band-by-band linear relationship between in situ field measurements and UAS sensor data [33], or the use of downwelling light sensors [34]. Previous studies compared UAS, airborne, and satellite data [35], since commercial sensors are now light enough to be embedded in low-weight platforms [36,37]. Zabala [38] compared VIS and NIR (VNIR) S2 bands with a multispectral consumer grade 2D-imager multicamera, such as MicaSense RedEdge [39], thereby concluding that a good correlation exists between derived vegetation indices on both sensors. MicaSense RedEdge was designed to be embedded in low-weight UAS platforms and provides acceptable geometric and radiometric data [40]. Parrot Sequoia [41] is a similar frame sensor, whose imagery radiometric performance was analyzed and compared with WorldView satellite data in [42] and was found to fit well between the corresponding bands. Apart from radiometric issues, the geometric benefit of accurate georeferencing of UAS data is expected to result in a better comparison of the ground-truth plotted area and the satellite pixel grid. In that sense, as drone data can systematically (raster, quasi-continuous) sample the ground at a distance below 10 cm and with a geometric error comparable to its pixel size [43,44], the drone data can be rigorously resampled to any satellite pixel size and grid.

The primary objective of this paper is to demonstrate that UAS data in synergy with field spectroradiometry can be useful for the accurate radiometric correction of L8 and S2 imagery. The underlying hypothesis is that with the appropriate upscaling methodology, radiometrically corrected UAS imagery through field spectroradiometric measurements [45] can be used as a reference to correct Landsat-like images in single dates with high accuracy (which we referred to as the UAS-MiraMon method), in a similar way as pseudoinvariant areas are used as references to correct Landsat-like time series in other approaches [23,46] (which we referred to as the PIA-MiraMon method). Additionally, a second objective is to show how the corrected drone data can be used to validate preexisting satellite surface reflectance methods, such as PIA-MiraMon, 6S-LaSRC, or Sen2Cor-SNAP. To realize this aim, we used near-coincident observations of L8 and S2 data, drone-acquired data, and field spectroscopy data.

These objectives were achieved by following a workflow consisting of the following steps (explained in detail in “Section 3 Material and methods”):

• Radiometrically correcting drone data by using calibration panels and field spectroradiometric measurements;
• Using drone corrected and upscaled data to validate different radiometric products (PIA-MiraMon, 6S-LaSRC and Sen2Cor-SNAP). This point fulfills the second main objective;
• Splitting the drone corrected and upscaled data into fitting pixels and test pixels. Fitting pixels are used to feed the UAS-MiraMon satellite radiometric correction approach based on field spectroradiometric data in synergy with drone data. Test pixels are used to validate the new method. This point fulfills the first main objective. Additionally, as a second validation we compare the results of UAS-MiraMon with the current 6S-LaSRC, Sen2Cor-SNAP and PIA-MiraMon methods. This complements the second main objective;
• By taking advantage of the almost simultaneous acquisitions of our experimental design, we evaluate the L8 and S2 inter-sensor fitting for the UAS-MiraMon method. This is done at different radial distances from the study area, as we hypothesized that the accuracy of the new UAS-MiraMon method could be better in the drone area. However, we also hypothesized that the accuracy could decrease with increasing distance (local overfitting). This point complements the first main objective.

2. Study Area

The data acquisition was carried out on April 21, 2018, at Can Gelabert farm (41°56′59.3″N–1°32′1.4″E) in the region of Catalonia, Spain, where the overpassing of L8 and S2A were scheduled at 10:35:37 Coordinated Universal Time (UTC) and at 10:56:29 UTC, respectively. The region of interest (ROI) covered 5.26 ha and included 51 OLI 30 m pixels, 143 MSI 20 m pixels, and 517 MSI 10 m pixels (Figure 1). The landscape is a Mediterranean mosaic dominated by a forest environment. The altitude of the terrain at the study area is between 683 m and 739 m. The topography is not absolutely flat, as there are some abandoned vineyard terraces and there are trees of less than 15 m in height. The terrain slope ranges from 0° to 56°, being especially high in the steps between terraces. In the central zone, there is a field of cereal crops that, at the time of the flight, was in the growing season. The woody vegetation of the area is mainly formed by Quercus humilis, Quercus ilex ssp. ballota, Ulmus minor, and Pinus nigra ssp. salzmanii. The man-made elements in the area are the house of Can Gelabert (650 m²), a pond, and a small well in the southeastern corner of the central cultivation field. We considered appropriately selecting a mixed natural and rural area to be representative for land cover mapping and monitoring purposes.

3. Materials and Methods

The materials we used involved three geographical scales: (a) field spectroradiometric data (Section 3.1) covering surface reflectance measurements, (b) UAS data (Section 3.2) covering drone imagery, and (c) satellite data (Section 3.3) covering L8 and S2A imagery. The methodology applied and the overall upscaling workflow are both summarized in the final subsection (Section 3.4).

3.1. Field Spectroradiometric Data

Ground reflectance measurements were carried out with an OceanOptics USB2000+ spectroradiometer [47]. The spectral range covered the VIS and NIR 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. The light was conducted on the portable instrument through a 50 µm diameter optical fiber of 25° field of view (FOV); the acquisitions were operated with the SpectraSuite software [48] (

Table 1. Main characteristics of the OceanOptics USB2000+ field spectroradiometer instrument. FOV: Field of View; CCD: Charge-Coupled Device; FWHM: full width at half maximum.

Manufacturer-Model	Size (mm)	Weight (g)	Design	Fiber Optic FOV (°)
OceanOptics USB2000+	89.1 × 63.3 × 34.4	190	Czerny-Turner	25
Sensor manufacturer and type	Sensor CCD samples	Sampling interval (nm)	Input Focal Length (mm)	Fiber optic diameter (µm)
CCD Sony ILX511B	2048	0.3	42	50
Raw radiometric resolution (bits) and expanded dynamic range (DN)	Grating #2 Spectral range (nm)	FWHM (nm)	Signal to Noise Ratio @ 50 ms	Minimum Integration time (ms)
8 (0–65,535)	340–1030	1.26	250:1	1

).

The dark current of the instrument was measured before each measurement. The optical fiber head was located in the nadir position at 20 cm, projecting an 8.87 cm circular footprint over the sampled surfaces. A Spectralon reference panel was used to approximate the incoming irradiance. Once the reference was taken, the individual target reflectance of 17 rectangular ethylene–vinyl acetate (EVA) foam panels of different colors was measured, averaging 100 readings to obtain a robust measurement. The set of EVA panels was composed of two sizes, both covering at least 6 pixels of the drone imagery: a group of 50 cm × 50 cm × 2 cm included black, violet, blue, dark green, yellow, orange, and red colors (namely EVA in Figure 2), and another group of 40 cm × 60 cm × 0.2 cm included black, blue, sky blue, grey, brown, red, green, dark green, yellow, and pink colors (namely foam in Figure 2), to demonstrate the possible work with thinner plates and a larger set of colors. Targets were separated by 10 cm to minimize adjacency effects and they were located over a black matte carpet to minimize background contributions. The EVA foams present a highly Lambertian [49,50] behavior and a high horizontal homogeneity [51]; the panels are light, easy to transport, and easy to clean with compressed air, making them a good material to be used for reflectance reference when calibrating drone imagery in synergy with field spectroradiometric measurements. The mean, minimum, maximum, and standard deviation statistics of the 100 readings were calculated to analyze and detect possible procedural errors or spectral variabilities. The spectral signature of each panel was integrated within the drone sensor’s relative spectral response function (RSRF) to be used in the empirical line correction [30]. The main motivation for using more than 15 EVA panels was to achieve a solid linear correlation between the ground-truth and the UAS sensor, with an extensive set of radiometric references distributed across the spectral domain (Figure 2).

3.2. Unmanned Aerial System Data

Two drone flights were planned at 80 m above ground level (AGL) with a side overlap of 70% and a front overlap of 80%, resulting in a pixel size of 6 cm. The first flight time was between 10:31:52 UTC and 10:41:42 UTC, which corresponded to the L8 overpass. The second flight time was between 11:23:10 UTC and 11:32:18 UTC, which corresponded to the S2A overpass. Both acquisitions were performed under excellent illumination and meteorological conditions without cloud cover and with wind speed below 2 m/s.

The UAS platform was a DJI S900 vertical take-off and landing 6-rotor multicopter coupled by the HEMAV Company to provide additional security measures. The drone used was 90 cm wide and had an 8.2 kg maximum take-off weight, a top speed of 9.7 m/s (35 km/h), and 15-min flight autonomy when carrying the payload. The UAS optical sensor was a multilens MicaSense RedEdge frame sensor [39] with a focal length of 5.5 mm and a lens field of view of 47.2°. The complementary metal-oxide semiconductor (CMOS) is 1280 × 960 pixels per frame (

Table 2. Main characteristics of the MicaSense RedEdge sensor. CMOS: Complementary Metal–Oxide–Semiconductor.

Manufacturer-MODEL	Size (mm)	Weight (g)	Raw Radiometric Resolution (Bits) and Expanded Dynamic Range (DN)
MicaSense RedEdge	121 × 66 × 46	150	12 (0–65,535)
Sensor type	Sensor size (pixels)	Pixel size (µm)	Focal Length (mm)	Field of View (°)
CMOS	1280 × 960	3.75 × 3.75	5.5	47.2
#1 Blue FWHM (nm)	#2 Green FWHM (nm)	#3 Red FWHM (nm)	#4 Red-edge FWHM (nm)	#5 NIR FWHM (nm)
468–491	548–568	666–676	712–723	814–865

). The MicaSense RedEdge is sensible in 5 solar spectrum bands, which, following the typical RSRF provided by the manufacturer [52], have a bandpass FWHM of 468 nm–491 nm (blue, #1), 548 nm–568 nm (green, #2), 666 nm–676 nm (red, #3), 712 nm–723 nm (red-edge, #5) and 814 nm–865 nm (near-infrared, #4) (Figure 3). The radiometric resolution is 12 bits, although the output is in the TIFF-16 bit format. The pixel size, which accounted for the focal length and the AGL distance, was rounded up to 6 cm.

UAS data was geometrically and radiometrically processed using the photogrammetric software Agisoft Photoscan v1.4.1 [53] and the remote sensing and geographic information system software (RS-GIS) MiraMon v8.1f [54]. Indirect georeferencing [43,44] was carried out using 7 ground control points (GCPs) measured with a post-processed Differential Global Positioning System (DGPS) Leica Viva GS14^® GNSS [55] static receiver and distributed across the ROI. Although all the available GCPs were used in the orthomosaic georeferencing, a cross validation (leave-one-out) was carried out to assess an independent accuracy test. GCPs were materialized in concrete bases and painted with a checkered template to display a clear visualization of the central point (Figure 4) when locating the markers in the photograms. The RMSExy of the first flight was 0.053 m (cross validation), and the RMSExy of the second flight was 0.044 m (cross validation). The geometric accuracy of both orthomosaics (6 cm pixel size) met the requirements of the American Society for Photogrammetry and Remote Sensing (ASPRS) Standards for Digital Geospatial Data [56], which demand a minimum RMSExy of 0.071 m for a pixel size of 5 cm.

Once the bundle block was adjusted and georeferenced, an orthomosaic of the surface reflectance values of each flight was generated. In this first radiometric correction processing step, we used the method based on the UAS camera manufacturer’s method and implemented it in the photogrammetric software [57], calibrating individual photograms thorough a grayscale reference panel (to account for the downwelling solar spectral irradiance) [33] together with exchangeable image file format (EXIF) information that considers the camera features (e.g., vignette correction factor) and the acquisition parameters (e.g., dark current, ISO, exposure time) [40]. In a second radiometric correction processing step, we used the field spectroradiometric measurements to evaluate the radiometric accuracy of the orthomosaic generated with the photogrammetric software, and furthermore, to improve the results by using an empirical line correction and the RS-GIS software [46]. The spectroradiometer pointer, located at 20 cm over the center of each target, projected a circle of 8.87 cm in diameter. In the image, a buffer with the appropriate radius was traced from the center of each of the 17 radiometric EVA targets in order to calculate the area-weighted average value of the image pixels that corresponded to the field spectroradiometric measurements. The image average reflectance of the 17 reflectance targets was correlated band-by-band with the field spectroradiometric measurements (integrated in the UAS sensor band RSRF). The correlation gave a coefficient of determination (R²) and a lineal function that, when applied to the UAS data, adjusted the UAS sensor reflectance values to the field spectroradiometric ground-truth reflectance, thus improving the surface reflectance reliability by applying an empirical line correction to the orthomosaic [45]. In empirical line calibration, several targets can be used. For example, when using two, the targets should be dark and bright, respectively. When several targets are used, they should cover the reflectance range well, but the reflectance spectra should not have intersections; otherwise, a large number of spectra is recommended to be used [15], and thus, we used a large number of targets.

3.3. Satellite Data

Two satellite multispectral push-broom sensors were used: The OLI onboard L8 platform and the MSI onboard S2A platform.

• L8 OLI: The L8 satellite was launched on February 11, 2013, and tracks a sun-synchronous polar-orbit at an approximate altitude of 705 km, providing a temporal resolution of 16 days. L8 has two main instruments on board, the OLI and the Thermal Infrared Sensor (TIRS), both with an approximate swath of 180 km and a radiometric resolution of 12 bits. The OLI sensor captures images in the visible (VIS), near-infrared (NIR), and shortwave infrared (SWIR) spectral regions, through 9 spectral bands of 30 m spatial resolution (SR) and an additional panchromatic band of 15 m SR; band RSRFs are available in [58]. On the other hand, although it is not the focus of this study, the TIRS obtains data in the thermal infrared spectral region (8 µm–12 µm) with two spectral bands that have a 100 m SR (

  Table 3. Main characteristics of the Landsat-8 and Sentinel-2A imagery used.

  Acquisition Date	21 April 2018
  Study Area	1°32′1.4″E, 41°56′59.3″N
  Platform	Landsat-8	Sentinel-2A
  Sensor	Operational Land Imager (OLI)	Multi Spectral Imager (MSI)
  Bands	10	13
  Aerosols	B1: 435.0 nm–451.0 nm (30 m)	B1: 432.7 nm–452.4 nm (60 m)
  Blue	B2: 452.0 nm–512.1 nm (30 m)	B2: 459.8 nm–524.0 nm (10 m)
  Green	B3: 532.7 nm–590.1 nm (30 m)	B3: 542.8 nm–577.6 nm (10 m)
  Red	B4: 653.9 nm–673.3 nm (30)	B4: 649.3 nm–679.9 nm (10 m)
  NIR(RedEdge)		B5: 697.3 nm–711.3 nm (20 m)
  NIR(RedEdge)		B6: 733.6 nm–747.3 nm (20 m)
  NIR(RedEdge)		B7: 773.5 nm–792.5 nm (20 m)
  NIR (wide)		B8: 782.5 nm–887.3 nm (10 m)
  NIR	B5: 850.5 nm–878.8 nm (30 m)	B8A: 854.5 nm–875.0 nm (20 m)
  Water Vapor		B9: 935.4 nm–954.9 nm (60 m)
  Cirrus	B9: 1363.0 nm–1384.0 nm (30 m)	B10: 1359.0 nm–1388.1 nm (60 m)
  SWIR1	B6: 1566.5 nm–1651.2 nm (30 m)	B11: 1568.7 nm–1658.3 nm (20 m)
  SWIR2	B7: 2107.4 nm–2294.1 nm (30 m)	B12: 2112.9 nm–2286.5 nm (20 m)
  Panchromatic	B8: 503.0 nm–676.0 nm (15 m)
  Thermal1	B10: 10.60 µm–11.19 µm (100 m)
  Thermal2	B11: 11.50 µm–12.51 µm (100 m)
  Acquisition time	10:35:49 (UTC)	10:56:29 (UTC)
  Orbit-Scene	198-031	R051-TCG
  Scene center	41°45′19.2″, 1°16′28.2″	41°56′41.6″, 1°14′58. 2″
  Solar elevation	55.67° (Scene center) 56.14° (Study area)	57.73° (Scene center) 57.57° (Study area)
  Solar Azimuth	144.75° (Scene center) 147.32° (Study area)	154.56° (Scene center) 153.80° (Study area)
  View Zenith Angle	7.26° (Study area)	1.23° (Study area)

  ). Official products were obtained from the USGS website via the EarthExplorer tool [59], which are distributed in the World Reference System-2 tiling system (WRS-2). In our study, we used images from orbit 198, row 031 (Figure 5) and two processing levels. The L1T processing level product (Level 1 T: precision terrain-corrected to UTM31N, WGS84, top of atmosphere (TOA) radiances product, Tier 1) was downloaded so that we could perform radiometric correction on it with the CorRad module of MiraMon, using, on the one hand, pseudo-invariant areas (PIA-MiraMon) [46], and on the other, the UAS imagery as the reflectance reference (UAS-MiraMon). The L1T product was also converted to TOA reflectance using the metadata parameters and the center-of-scene sun elevation angle [60] for comparison with the radiometric correction results and the drone data. The L2A processing level product (Level 2 A: geometry as L1T and surface reflectance corrected for atmospheric effects with a cloud mask and a cloud shadow mask, as well as a water and snow mask), atmospherically corrected with the Landsat Surface Reflectance Code (6S-LaSRC) [14], was downloaded to have an independent reference of the results obtained with the MiraMon algorithm and for comparison with the official product of S2. On 21 April 2018, L8 overpassed the study area at 10:35:49 UTC (Table 3).
• S2A MSI: The S2A satellite was launched on 23 June 2015 and is part of the European Space Agency (ESA) Copernicus space program. The S2A twin S2B satellite was put into orbit on March 7, 2017. Both track a sun-synchronous polar-orbit at an approximate altitude of 786 km, providing an individual temporal resolution of ten days and a combined temporal resolution of five days. S2A has a main instrument on board, the MSI, which has an approximate swath of 290 km and a radiometric resolution of 12 bits. The MSI sensor captures images in the visible (VIS), near-infrared (NIR), and shortwave infrared (SWIR) spectral regions through 4 spectral bands of 10 m SR, 7 bands of 20 m SR, and 3 bands of 60 m SR (Table 3); band RSRFs are available in [61] and are slightly different for the MSI onboard S2A and the MSI onboard S2B. Official products were obtained via the ESA—Copernicus Scientific Data Hub tool [62], which are distributed in the Military Grid Reference System tiling system (MGRS). In our study, we used images from the orbit R051 granule TCG and two processing levels: L1C (Level 1C: product results from using a digital elevation model to project the image in cartographic geometry (UTM31N, WGS84), where per-pixel radiometric measurements are provided in TOA reflectance); and L2A (Level 2A: Bottom of atmosphere (BOA reflectance in cartographic geometry (as for L1C)) (Figure 4). L1C was downloaded to be radiometrically corrected using PIA-MiraMon and UAS-MiraMon. The L2A product, radiometrically corrected with the Sen2Cor processor [17], was downloaded to have an independent reference of the results obtained with the MiraMon algorithm and be compared with the L8 official product. On April 21, 2018, the S2A overpassed the study area at 10:56:29 UTC.

The radiometric performance requirements of the Landsat-8 [63] and Sentinel-2 [64] missions require producing data with a radiometric uncertainty below 3%. This is the threshold we used to discuss if the radiometric correction methods are compliant with the requirements, considering that a RMSE is good if it is below the 3% threshold and is bad if it overpasses this threshold. Other works use a 2% threshold [22,65], below which we could consider it to be an excellent result.

The satellite bands involved were OLI 452 nm–512 nm (blue, #2), 533 nm–590 nm (green, #3), 636 nm–673 nm (red, #4), and 851 nm–879 nm (near-infrared, #5); and MSI 460 nm–524 nm (blue, #2), 543 nm–578 nm (green, #3), 649 nm–680 nm (red, #4), and 855 nm–875 nm (near-infrared, #8a) (Figure 6). Note that all OLI bands are 30 m SR, MSI visible bands are 10 m SR, and the NIR band (#8a) is 20 m SR. The official RSRF was used to select the matching bands between both sensors and determined the bands to be compared with the UAS sensor. Also note that the OLI NIR band (#5) is very similar to the MSI NIR band #8a (designed to match the aforementioned OLI NIR band) (Figure 6), although the MSI NIR band has a wider NIR band (#8) in the same spectral region (which was designed to match the NIR band of the ETM+ of Landsat-7).

The corresponding bands of the satellite sensors and UAS are located in the blue, green, red, and near-infrared regions of the electromagnetic spectrum, but have different spectral overlapping. In the blue band, the OLI and the MSI response are similarly wide, but the MSI band is displaced by about 20 nm to longer wavelengths. In the green band, the MSI band is narrower than the OLI band. In both the blue and green bands, the UAS sensor response is reasonably within the same spectral region and centered in the RSRF of the satellite sensors. In the red region, the MSI FWHM is narrower than the OLI and is displaced to lower frequencies, while the UAS sensor response matches the upper boundary of the satellites’ RSRF, being decentered from the satellite RSRFs. The OLI NIR band is slightly wider than the MSI band; nevertheless, both reasonably match the same spectral region. At any rate, the UAS NIR band peak is located at about 20-nm higher frequencies and only matches the satellite bands within a 10-nm spectral region where the response decreases such that, a priori, a weaker correlation is expected between the satellite and drone imagery along these bands (Figure 7).

Satellite imagery was masked by clouds and cloud shadows, accounting for both the clouds in the L8 image and the clouds in the S2A image, because cloud movement and growth was notable during the 20-min delay period between the overpass of both platforms, and as such, could cause a bad interpretation of the intercomparison correlations. The total matching area of the L8 198031 scene and the S2 R051 TCG tile was 1,204,945 ha, and the masked area was 111,562 ha (9.26% of the common area) (Figure 8).

The CorRad algorithm was designed in 1994 to correct the topographic and atmospheric effects of Landsat-like imagery [24]. The topographic corrections followed a cosine-correction approach. The atmospheric correction followed a dark object subtraction (DOS) [66] approach to resolve the atmospheric path radiance, while the atmospheric total optical depth was assumed to be constant (wavelength-dependent). In 2014, Pons et al. [46] presented an automatic method designed to massively correct Landsat-like imagery, thus solving the path radiance and the atmospheric total optical depth by fitting these unknowns to radiometric reflectance references obtained from Moderate-Resolution Imaging Spectroradiometer (MODIS) pixels that remained under a standard deviation threshold during the whole MODIS time series, namely pseudo-invariant areas (PIA). When a given Landsat-like image is processed with the PIA-MiraMon method, in each PIA, the CorRad algorithm resolves the optimal path radiance and the atmospheric total optical depth needed to obtain the reference value, and this process is carried out for all the PIA within the processed image extent. The PIA-MiraMon gets a band-by-band pairing of path radiance and atmospheric total optical depth values to correct the full scene.

However, the focus of this study was to use resampled UAS imagery as radiometric reflectance references to process a given Landsat-like image with the CorRad-MiraMon algorithm (UAS-MiraMon), thereby aiming to resolve the optimal path radiance and the atmospheric total optical depth in the study area in a still finer way by taking advantage of the current availability of drones and light multispectral sensors. The L1T L8 and S2 scenes were corrected with these atmospheric parameters, expecting to correct the satellite image with great accuracy near the UAS drone flight and also anticipating a possible worse-case adjustment in distant areas.

3.4. Summary of Upscaling Workflow and Validation from Field Measurements to Satellite Data

First Upscaling Step, from Ground-Truth to UAS:

The field spectroradiometric data at <1 nm spectral sampling was resampled to the drone sensor bandwidth, accounting for the RSRF provided by the manufacturer. Spectrally resampled ground-truth data was correlated band-by-band to the drone data (Figure 2).

Second Upscaling Step, from UAS to Satellite:

The weighted average value was calculated to resample drone data (6 cm) to the satellite pixel size (30 m in the case of L8, and 10 m or 20 m in the case of S2) using a 2D Gaussian procedure described in [67], similar to the one used in [68] to resample Landsat imagery to the MODIS pixel size. There were 250,000 drone pixels in a L8 30 m pixel, 111,111 drone pixels in a S2 20 m pixel, and 27,776 drone pixels in a S2 10 m pixel. We validated the overall consistency of the UAS data calculating the band-by-band RMSE between all the UAS resampled pixels and the corresponding satellite L1T product, the L2A official product (6S-LaSRC in the case of L8 and Sen2Cor-SNAP in the case of S2), and the PIA-MiraMon surface reflectance product.

Third Upscaling Step, Using UAS Data to Radiometrically Correct Satellite Data:

Of the total drone-resampled pixels, 40% were reserved for validation purposes as test areas, while 60% were used as radiometric references to perform the radiometric correction of satellite data using the CorRad algorithm (UAS-MiraMon). The RMSE between test areas and the L8 or S2 VNIR bands was calculated for the UAS-MiraMon, the PIA-MiraMon, and the two official products.

Fourth Upscaling Step, Radiometric Corrections Intersensor Comparison:

The availability of near-coincident observations of the L8 and the S2A, being less than 22 min apart, allowed us to correlate not only the UAS-MiraMon surface reflectance provided for the L8 and the S2A, but also the surface reflectance provided for the L8 official product and the S2 official product, and the PIA-MiraMon surface reflectance provided for the L8 and the S2A. To perform the correlation between these products with different pixel sizes and a misaligned grid, the MSI grid was adapted to the 30 m OLI grid by calculating the area-weighted pixel contribution of the 10 m or 20 m pixels to the 30 m resampled pixel, since the number of MSI pixels inside an OLI pixel does not exceed 16 units (Figure 9).

Once the S2 images were adapted to the L8 spatial grid (30 m pixel size), the band-by-band RMSE in reflectance values was calculated to intercompare the UAS-MiraMon radiometric correction method in near time-coincident acquisitions. The differences were calculated in concentric rings from the center of the study area in order to search for a possible degradation of L8–S2 correlation in the UAS-MiraMon method. The analyzed zones were from 0 km to 15 km (n = 681,950 pixels), 15 km to 30 km (n = 2,203,730), and additionally, the full TCG scene extent (n = 12,148,700) (Figure 10).

The overall methodology and the results can be summarized in (a) the fitting between the ground-truth and drone data, (b) the fitting between the drone data and satellite data, and (c) the fitting between the radiometric corrections of almost simultaneously acquired data with different satellite sensors. The full upscaling workflow from field measurements to satellite data can be seen in Figure 11.

4. Results

4.1. Fitting between Ground-Truth and Drone Data

The UAS imagery acquired by the overpass of L8 showed an excellent fitting to ground-truth measurements. VIS bands showed a high R² (R²_BLUE > 0.955; R²_GREEN > 0.967; R²_RED > 0.982) and relatively low negative bias (−3.346%; −4.067%; −3.762%), meaning that the UAS reflectance values were slightly but systematically greater than those measured with the field spectroradiometric measurements. In the NIR band, this effect was relevant in dark targets, where the radiometer-measured surface reflectance was below 5%, while the UAS sensor-measured reflectance was between 15% and 20%, resulting in a relatively high negative bias (−16.956%), albeit also showing a high coefficient of determination (R²_NIR > 0.942) (Figure 12a).

The UAS imagery acquired by the overpass of S2A obtained results that are very similar to the precedent L8. It also showed an excellent fitting to ground-truth measurements. VIS bands showed a high R² (R²_BLUE > 0.971; R²_GREEN > 0.969; R²_RED > 0.975) and relatively low negative bias (−2.534%; −3.652%; −2.454%), also meaning that the UAS reflectance values were slightly but systematically greater than those measured with the field spectroradiometric measurements. In the NIR band, the effect was the same as those mentioned in the dark targets, where the radiometer-measured surface reflectance values were below 5%, while the UAS sensor-measured values were between 15% and 18% in this case, resulting in a relatively high negative bias (−17.921%), albeit also showing a high coefficient of determination (R²_NIR > 0.946) (Figure 12b).

It is worth noting that the correlations in all the cases are significant, with p-values of <0.05.

When the linear function was applied, the bias was corrected and the coefficient of determination was maintained, thus accurately fitting the image values to the field spectroradiometric reflectance. The differences between both UAS images were in practice very low, since a delay of less than 22 min from one flight to the other only affected the sun position change, the small portion of atmosphere between the panels and the sensors, and the slight view difference.

4.2. Fitting between the Drone Data and the Satellite Data

Drone data (6 cm) resampled to 51 OLI grid pixels (30 m), 143 MSI grid pixels (20 m), and 517 MSI grid pixels (10 m) was compared band-by-band and pixel-by-pixel with official surface reflectance products (6S-LaSRC, Sen2Cor-SNAP) and with the PIA-MiraMon method. Results showed that the UAS data and satellite data follow the same trend, both for the OLI imagery (Figure 13a) and the MSI imagery (Figure 13b), although obviously the higher spatial resolution of the MSI increased the pixel variability over a smoother trend in the OLI sensor.

Focusing on the blue band, the PIA-MiraMon radiometric correction values were generally closer to the UAS data (RMSE_OLI ≤ 0.703%, RMSE_MSI ≤ 1.549%) than the 6S-LaSRC (RMSE_OLI ≤ 2.054%) or the Sen2Cor-SNAP (RMSE_MSI ≤ 2.192%) data (

Table 4. Band-by-band RMSE between the UAS resampled imagery, the satellite data corrected with the PIA-MiraMon method and the official surface reflectance products for the OLI (6S-LaSRC) and the MSI (Sen2Cor-SNAP). Bold values are the lowest difference in the corresponding band and sensor, being below the 3% accuracy requirement. Units: Reflectance (%).

	BLUE		GREEN		RED		NIR
RMSE UAS Data vs. Data Source	OLI #2 (n = 51)	MSI #2 (n = 517)	OLI #3 (n = 51)	MSI #3 (n = 517)	OLI #4 (n = 51)	MSI #4 (n = 517)	OLI #5 (n = 51)	MSI #8a (n = 143)
TOA reflectance	7.743	6.981	3.684	3.674	3.533	3.101	6.847	8.050
PIA-MiraMon	0.703	1.549	1.745	2.477	3.017	2.802	9.177	11.349
Official products	2.054	2.192	1.623	2.058	2.484	2.361	6.688	8.783
RMSE PIA-MiraMon vs. Official products	1.882	0.818	0.746	0.715	1.061	0.780	3.693	2.757

). The UAS surface reflectance data was always below that of the satellite, except in a bright pixel located over the house’s roof. The MSI blue band showed the same trend as the OLI band, whereby it exhibited a good agreement between UAS and satellite surface reflectance data, despite the higher interpixel variability of the MSI blue band due to its higher spatial resolution (Figure 13).

Focusing on the green band, the PIA-MiraMon and official products’ surface reflectance values were both close to the UAS data, as the pixel-by-pixel plot indicates (Figure 13), but the 6S-LaSRC (RMSE_OLI ≤ 1.623%) and the Sen2Cor-SNAP (RMSE_MSI ≤ 2.058%) obtained slightly better adjustments than the PIA-MiraMon (RMSE_OLI ≤ 1.745%, RMSE_MSI ≤ 2.477%) in this spectral region (Table 4).

Focusing on the red band, the differences in radiometric corrections over the UAS data were slightly higher for the PIA-MiraMon (RMSE_OLI ≤ 3.017%, RMSE_MSI ≤ 2.802%), LaSRC (RMSE_OLI ≤ 2.484%), and Sen2Cor (RMSE_MSI ≤ 2.361%). The pixel-by-pixel comparison showed a high variability of the UAS red band, with a range of values wider than in the green and blue bands (Figure 13).

Focusing on the NIR band, the differences were dramatically high between the UAS data and the satellite data (Table 4) in both the PIA-MiraMon (RMSE_OLI ≤ 9.177%, RMSE_MSI ≤ 11.349%) and the official products (RMSE_OLI ≤ 6.688%, RMSE_MSI ≤ 8.783%). Although the overall trend was coherent, these differences were expected due to the RSRF drone sensor and the satellite sensors mismatching in this band (Figure 13).

With these results, we were concerned about the high differences between the UAS surface reflectance data and the satellite surface reflectance data in the NIR band. Thus, we compared the satellite surface reflectance products to investigate if there were higher differences than between the UAS data and satellite data. The differences between the PIA-MiraMon and official products (last row of Table 4) in the blue band (RMSE_OLI ≤ 1.882%, RMSE_MSI ≤ 0.818%), green band (RMSE_OLI ≤ 0.746%, RMSE_MSI ≤ 0.715%), and red band (RMSE_OLI ≤ 1.061%, RMSE_MSI ≤ 0.780%) indicated a general agreement of both radiometric correction methods, thus either confirming the differences below, or that they are equal to 3.0% between the UAS data and the satellite data, which were radiometrically corrected using the other methods. In the NIR band, those differences were below 4% between the PIA-MiraMon and the 6S-LaSRC (RMSE_OLI ≤ 3.693%) and also between PIA-MiraMon and Sen2Cor-SNAP (RMSE_MSI ≤ 2.757%), which definitely points to a disagreement between the UAS data and the satellite data in this band (Figure 8).

4.3. Fitting of the UAS-Based Satellite Radiometric Corrections with the Test Areas

Comparing the test pixels (40% of 51) in the corrected L8 image (30 m) with the UAS-MiraMon method, for band-by-band reflectance, the RMSE was 0.712% in the blue band, 1.000% in the green band, 1.871% in the red band, and 7.530% in the NIR band, as expected based on the previous section’s results. Comparing the test pixels (40% of 517) in the corrected S2 (10 m) image with the UAS-MiraMon method, for band-by-band reflectance, the RMSE was 0.745% in the blue band, 1.230% in the green band, and 2.018% in the red band; when comparing the test pixels (40% of 143) in the NIR band (20 m), the RMSE was 6.760%, also in consonance with previous results (

Table 5. Band-by-band RMSE between the test UAS areas and the satellite data radiometrically corrected by using: fitting UAS references (UAS-MiraMon), MODIS references (PIA-MiraMon) and the official products (6S-LaSRC for OLI and Sen2Cor-SNAP for MSI). Bold values are the lowest difference in the corresponding band and sensor, being below the 3% accuracy requirement. Units: Reflectance (%).

	BLUE		GREEN		RED		NIR
RMSE Test Areas vs. Data Source	OLI #2 (n = 20)	MSI #2 (n = 207)	OLI #3 (n = 20)	MSI #3 (n = 207)	OLI #4 (n = 20)	MSI #4 (n = 207)	OLI #5 (n = 20)	MSI #8a (n = 57)
UAS-MiraMon	0.712	0.745	1.000	1.230	1.871	2.018	7.530	6.760
PIA-MiraMon	0.585	1.366	1.746	1.867	3.130	2.373	8.528	7.281
Official Products	2.197	2.149	1.653	1.806	2.680	2.192	5.841	7.690

).

When calculating the RMSE between test areas and the other radiometric correction methods, in the blue band, they were higher in the official products (RMSE_OLI ≤ 2.197%, RMSE_MSI ≤ 2.149%) than with the PIA-MiraMon correction method (RMSE_OLI ≤ 0.585%, RMSE_MSI ≤ 1.366%). In the green band, the official products (RMSE_OLI ≤ 1.653%, RMSE_MSI ≤ 1.806%) were slightly better than with the PIA-MiraMon method (RMSE_OLI ≤ 1.746%, RMSE_MSI ≤ 1.806%), as in the red band, where 6S-LaSRC and Sen2Cor-SNAP (RMSE_OLI ≤ 2.680%, RMSE_MSI ≤ 2.192%) obtained lower differences when compared to test areas than with the PIA-MiraMon method (RMSE_OLI ≤ 3.130%, RMSE_MSI ≤ 2.373%) (Table 5).

4.4. Fitting between the Radiometric Correction of Almost Simultaneously Acquired L8 and S2 Data

The surface reflectance of the L8 and S2A images acquired at almost at the same time showed acceptable differences in the visible bands (RMSE_VIS < 1.147%) and higher ones in the NIR band (RMSE_VIS < 1.707%) (

Table 6. Band-by-band correlation between almost simultaneous acquired images of L8 and S2 (resampled to 30 m pixels). Analysis performed in concentric radial distances from the center of the study area and in the R051 TCG tile extent. The radiometric correction was derived from data using UAS references (UAS-MiraMon).

RMSE (%) between Near Coincident L8 and S2 Images Corrected with the UAS-MiraMon Method
Distance to UAS Central Scene (Number of Inter-Compared Pixels)	BLUE	GREEN	RED	NIR
0–15 km (n = 681,950)	0.621	0.741	0.919	1.664
15–30 km (n = 2,203,730)	0.739	0.815	0.984	1.666
TCG scene (n = 12,148,700)	0.852	0.940	1.147	1.707

). The blue band had the lowest intersensor differences among all the performed comparisons, either when focusing 15 km around the study area (0.621%), in the ring between 15 km and 30 km (0.739%), or in the full TCG scene (0.852%). The green band presented a good fitting between both sensors, with differences lower than 1% in the three distances (0.741%, 0.815%, and 0.940%). In the red band, the fittings showed slightly higher differences (0.919%, 0.984%, and 1.147%) than in the blue and green bands. In the NIR band, the differences increased with respect to the visible bands (1.664%, 1.666%, 1.707%), but they followed the same distance pattern.

5. Discussion

First of all, we would like to present that obtaining the field spectroradiometric measurement with the hand-held instrument and the reflectance reference EVA panels was operational and provided very useful data for this research. Analyzing the results of the fitting between the ground-truth and the 6-cm drone data described in the previous section, we found that the photogrammetric software reflectance calibration using the calibration target provided by the MicaSense RedEdge manufacturer and EXIF information obtains good fitting with in situ EVA panel measurements (R² > 0.942), resulting in a good performance of the MicaSense RedEge sensor and quite good radiometric treatment of the software. Although we found very low differences between the fitting of the in situ data and the two UAS flights carried out by the L8 overpass and by the S2A overpass, we should highlight a common negative bias, which means that the automatic UAS reflectance values are systematically higher than the field spectroradiometric reflectance values. In summary, in this first upscaling step, we were able to determine that using multiple reflectance targets distributed in the spectral domain, covering a wide range of intensities and characterized with field spectroradiometry measurements, facilitates the correction of drifts introduced by the simpler usage of the calibration target provided by the drone sensor manufacturer, all of which are consistent with very recent studies [30]. It is worth noting that two or three high-quality targets should be enough for calibration purposes, and the remaining targets could be used for the validation of results and accuracy analysis, but according to our experience, more panels could be used to guarantee an appropriate distribution of intensities and spectral references, and also to have some alternatives if some color is so intense that it saturates the UAS data. Nevertheless, as this research did not attempt to determine the minimum number of panels or the more appropriate colors, we have used a fairly large number of panels (17) and a wide set of colors in order to be sure that no problem occurred because of that. Further research will determine the appropriateness of the number of panels and the concrete colors we have suggested after more intensive tests focused on this objective.

Analyzing the results of the fitting between the drone data and the satellite data, the RMSE between the UAS and the satellite data was acceptable in visible bands (RMSE < 3.017%), as supported in the recent work of [38]. The spatial resampling of the drone data from 6 cm to 30 m, 20 m, and 10 m using a 2D Gaussian weighting function provided synthetic values coherent with the values obtained from the satellite imagery that was radiometrically corrected using the official products and a well-proven method, such as the PIA-MiraMon. In this second upscaling step, we have to highlight three effects among the comparisons in all the bands, both in the OLI and the MSI. Firstly, the UAS surface reflectance values were systematically lower than satellite values (this effect was most visible in the blue and green bands). Secondly, the UAS data was more diverse than the satellite data; also, the 10 m and 20 m data were more diverse than the 30 m data; this is expected due to a larger pixel size (Figure 13). Third, and most importantly, the NIR band presented very high reflectance differences between the UAS data and the satellite data (RMSE > 6.688%), both in the OLI and the MSI and in all the radiometric correction methods, thus confirming that the band mismatch in the spectral configuration is too significant to be ignored. The spectral mismatching is the main source of this high difference, but in further works where the UAS band and the satellite NIR band have a good spectral matching, it will be necessary to consider if the signal-to-noise ratio (SNR) is the main significant source of error, because in OLI and MSI sensors, the SNR is poorer than in the VIS bands (due to the lower radiance signal in IR wavelengths), as stated in the Landsat-8 Data users handbook [63] and in the Sentinel-2 Technical Guide [64].

Analyzing the results of the fitting between the drone-based satellite radiometric corrections and test areas, we found that in the visible bands, the UAS-MiraMon radiometric correction achieves similar-quality results to the 6S-LaSRC, Sen2Cor-SNAP, or PIA-MiraMon alternatives (in the NIR band, the results are highly biased and cannot be compared), meeting the radiometric accuracy requirements. Indeed, in the visible bands of the OLI and MSI, the UAS-MiraMon achieves better fittings (RMSE < 2.018%) than the other tested methods (6S-LaSRC (RMSE < 2.680%), Sen2Cor-SNAP (RMSE < 2.192%), and PIA-MiraMon (RMSE < 3.130%)). This could be explained by the fact that the 6S-LaSRC, Sen2Cor-SNAP, and PIA-MiraMon use information for the whole scene, whereas in the UAS-MiraMon case, the atmospheric parameters of the CorRad-MiraMon model (L_atmλ, τ_0λ) are fitted using the local drone references. Consequently, using the UAS references has improved local accuracy in radiometric correction, arising as a new methodology to radiometrically correct satellite imagery focusing in a special area of interest, rather than correcting a full scene with atmospheric parameters obtained from other sources, often out of the area of interest, as is done by the official products and other alternatives. Moreover, the presented method allows a much more robust usage of the point spectroradiometric measurements with regard to their direct comparison to satellite imagery. On the other hand, the UAS-MiraMon method should be used for correcting satellite images acquired almost simultaneously to UAS data.

Consistent with other findings [23], when analyzing the results of the surface reflectance differences between almost simultaneously acquired L8 and S2 data radiometrically corrected with UAS-MiraMon, we found that the method obtains good fittings (RMSE < 1.7%). As we hypothesized, the accuracy of the new UAS-MiraMon method is better in the drone area and decreases with increasing distance (local overfitting), since the solution for the atmospheric parameters is fitted to the study area surface reflectance references. However, the results for the coincident S2 and L8 tiles extent are also good (RMSE_VIS < 1%), thus implying a good consistency of the UAS-MiraMon experimental approach. This high intersensor coherence, added to the improved radiometric correction accuracy in the area of interest, is an advantage of the UAS-MiraMon method for the generation of consistent temporal series using data from different sensors in specific areas of interest that need to be monitored in an accurate way.

6. Conclusions

In this research, we demonstrated the operational capabilities of unmanned aerial system (UAS) imagery in synergy with field spectroradiometric measurements to perform radiometric correction of satellite imagery. The use of UAS imagery as radiometric references in combination with the CorRad-MiraMon radiometric correction algorithm (UAS-MiraMon method) resulted in improved local accuracy of surface reflectance values of the satellite imagery. Radiometrically corrected official products (6S-LaSRC and SNAP-Sen2Cor) and other alternatives (PIA-MiraMon) met the mission radiometric accuracy requirements as well (RMSE < 3%), but UAS-MiraMon obtains better accuracy because it focuses on the retrieval of atmospheric parameters in the area of interest, rather that correcting a full scene with atmospheric parameters obtained from other sources, often out of the area of interest. Also, the UAS imagery covers a larger area than conventional field spectroscopy measurements with reliable radiometric quality, and can be useful to support the validation of surface reflectance values of satellite imagery in specific areas of interest as well as their continuous monitoring with higher radiometric accuracy.

Correlating the drone data with the field spectroradiometric measurements was evaluated positively. However, a linear fitting using reference panels of different spectral response and intensity is required to improve the simpler correction based on the manufacturer’s calibration panel. The UAS data provided an ultrahigh spatial resolution (6 cm) and sampled areas not accessible by human operators in previous field spectroradiometric campaigns; the results indicated consistent surface reflectance accuracy. One of the best methodological advantages of using the UAS instruments is the high sampling capacity over space, allowing an easy and efficient resampling to any satellite grid; in this study, we covered 517 MSI 10 m pixels in ±5 min with the satellite overpass.

Altogether, the results showed that when the spectral response of the drone sensor reasonably matches the response of the satellite sensors, the results are good. However, a mismatch between them is a problematic issue, and we encourage the design of UAS sensors with spectral responses more similar to OLI or MSI sensors, especially in the NIR region. The UAS data resampled to satellite pixel sizes proved to be highly coherent with the official surface reflectance products in the visible bands of the L8 and the S2A. Nevertheless, use of the UAS data as a radiometric reference to fit the atmospheric parameters of the radiometric correction in combination with the CorRad-MiraMon algorithm revealed very good results in the study area. Moreover, although the surface reflectance accuracy slightly degrades with increasing distance from the study area, the comparison of near-coincident L8 and S2 images corrected with the UAS-MiraMon also provided an excellent intersensor fitting (RMSE_VIS < 1%), thus suggesting a robust way for multisensor radiometric correction and the generation of more coherent time series focused on a specific study area that needs to be accurately monitored.

Field spectroradiometric measurement remains as the ground-truth reference source due to the proximity of the sensor to the surface, the intensity of scans, the control of the measurements, and the instrument features. Thus, the UAS-acquired data can benefit from conventional field spectroradiometry. Moreover, this type of synergy can provide valuable surface reflectance reference data.

Radiometric correction of Landsat-like images using drone imagery as reflectance references provides good results in the visible bands (RMSE_VIS < 2.018%), meeting the mission radiometric requirements. This method is especially useful for the surrounding areas of the drone flight, such as in the studying and monitoring of protected areas. This contribution opens the door to the operational usage of spatially disperse drone flights campaigns, in synergy with ground-truth reflectance references, for radiometrically correcting satellite images or validating radiometric corrections results, especially when using remote sensors with a similar spectral configuration.