Nothing Special   »   [go: up one dir, main page]

Next Article in Journal
Land Cover Classification in SubArctic Regions Using Fully Polarimetric RADARSAT-2 Data
Previous Article in Journal
A Comparison between Drifter and X-Band Wave Radar for Sea Surface Current Estimation
You seem to have javascript disabled. Please note that many of the page functionalities won't work as expected without javascript enabled.
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Technical Note

A Software Tool for Atmospheric Correction and Surface Temperature Estimation of Landsat Infrared Thermal Data

1
CESBIO, Université de Toulouse, CNES/CNRS/IRD/UPS, 18 av. Edouard Belin, bpi 2801, F-31401 Toulouse Cedex 9, France
2
CNES, 18 av. Edouard Belin, F-31400 Toulouse, France
*
Author to whom correspondence should be addressed.
Remote Sens. 2016, 8(9), 696; https://doi.org/10.3390/rs8090696
Submission received: 25 March 2016 / Revised: 1 August 2016 / Accepted: 2 August 2016 / Published: 24 August 2016
Graphical abstract
">
Figure 1
<p>Localization map of the both network (green triangles), for France (<b>left</b>) and for Tunisia (<b>right</b>). Background are true-color images of L7 (<b>left</b>) and L8 (<b>right</b>) acquired respectively on 29 May 2003 and 29 March 2013. Images are scaled and oriented. The black boxes on images represent the ASTER image footprint used on each location to validate tool (see <a href="#sec4dot2-remotesensing-08-00696" class="html-sec">Section 4.2</a>).</p> ">
Figure 2
<p>Flowchart of the LANDARTs algorithm to generate LST images at a resolution of 30 m.</p> ">
Figure 3
<p>Comparison between Landsat LST and in situ temperatures (in degree Celsius) for the Tunisian site.</p> ">
Figure 4
<p>Comparison between Landsat LST and in situ temperatures (in degree Celsius) for the French site.</p> ">
Figure 5
<p>Pixel-by-pixel comparison between the ASTER surface kinetic product and the Landsat LST map at a resolution of 270 m for the Tunisian site.</p> ">
Figure 6
<p>Pixel-by-pixel comparison between the ASTER surface kinetic product and the Landsat LST map at a resolution of 270 m for the French site.</p> ">
Figure 7
<p>Difference in degree Celsius between ASTER surface kinetic product and Landsat LST map at 270 m resolution on Tunisia (29 March 2013).</p> ">
Figure 8
<p>Disparity histogram in degree Celsius between ASTER surface kinetic product and Landsat LST map at 270 m resolution on Tunisia (29 March 2013).</p> ">
Figure 9
<p>Difference comparison in degree Celsius between ASTER surface kinetic product and Landsat LST map at 270 m resolution on France (29 May 2003).</p> ">
Figure 10
<p>Disparity histogram in degree Celsius between ASTER surface kinetic product and Landsat LST map at 270 m resolution on France (29 May 2003).</p> ">
Figure 11
<p>Parameter variability and influence on LST correction for an L7 image taken on 29 May 2003 (199,030). The square black frame represents the ASTER valid data footprint used for the spatial comparison. (<b>a</b>) Transmittance (unitless); (<b>b</b>) upwelling radiance (<math display="inline"> <semantics> <mrow> <mi mathvariant="normal">W</mi> <mo>·</mo> <msup> <mi mathvariant="normal">m</mi> <mrow> <mo>−</mo> <mn>2</mn> </mrow> </msup> <mo>·</mo> <msup> <mi>sr</mi> <mrow> <mo>−</mo> <mn>1</mn> </mrow> </msup> <mo>·</mo> <mi mathvariant="sans-serif">μ</mi> <msup> <mi mathvariant="normal">m</mi> <mrow> <mo>−</mo> <mn>1</mn> </mrow> </msup> </mrow> </semantics> </math>); (<b>c</b>) map of emissivity estimates derived from NDVI values used in both corrections (unitless); (<b>d</b>) disparity map, in degree Celsius, comparing LST estimated using spatial parameters with LST estimated using the center-scene correction parameters. White pixels are cloudy or snow-covered pixels detected with the level 2A processor MACCS [<a href="#B45-remotesensing-08-00696" class="html-bibr">45</a>].</p> ">
Figure 12
<p>Multi-annual temporal evolution of the atmospheric transmittance parameter (unitless) observed for the set of French images (path/row 199/030). The blue line represents the center-scene transmittance, and the red line is the transmittance for the Aurade site extracted from the spatial parameter image at a resolution of 30 m. The green fill represents the amplitude (min/max) of the transmittance for each Landsat scene.</p> ">
Figure 13
<p>Multi-annual temporal evolution of the atmospheric upwelling radiance parameter (<math display="inline"> <semantics> <mrow> <mi mathvariant="normal">W</mi> <mo>·</mo> <msup> <mi mathvariant="normal">m</mi> <mrow> <mo>−</mo> <mn>2</mn> </mrow> </msup> <mo>·</mo> <msup> <mi>sr</mi> <mrow> <mo>−</mo> <mn>1</mn> </mrow> </msup> <mo>·</mo> <mi mathvariant="sans-serif">μ</mi> <msup> <mi mathvariant="normal">m</mi> <mrow> <mo>−</mo> <mn>1</mn> </mrow> </msup> </mrow> </semantics> </math>) observed for the set of French images (path/row 199/030). The blue line represents the center-scene upwelling radiance, and the red line is the upwelling radiance for the Aurade site extracted from the spatial parameter image at a resolution of 30 m. The green fill represents the amplitude (min/max) of the upwelling radiance for each Landsat scene.</p> ">
Figure 14
<p>Comparison in degree Celsius between Landsat LST and Tunisia in situ temperature for the differing correction methods.</p> ">
Figure 15
<p>Comparison in degree Celsius between Landsat LST and French in situ temperature for the differing correction methods.</p> ">
Versions Notes

Abstract

:
Land surface temperature (LST) is an important variable involved in the Earth’s surface energy and water budgets and a key component in many aspects of environmental research. The Landsat program, jointly carried out by NASA and the USGS, has been recording thermal infrared data for the past 40 years. Nevertheless, LST data products for Landsat remain unavailable. The atmospheric correction (AC) method commonly used for mono-window Landsat thermal data requires detailed information concerning the vertical structure (temperature, pressure) and the composition (water vapor, ozone) of the atmosphere. For a given coordinate, this information is generally obtained through either radio-sounding or atmospheric model simulations and is passed to the radiative transfer model (RTM) to estimate the local atmospheric correction parameters. Although this approach yields accurate LST data, results are relevant only near this given coordinate. To meet the scientific community’s demand for high-resolution LST maps, we developed a new software tool dedicated to processing Landsat thermal data. The proposed tool improves on the commonly-used AC algorithm by incorporating spatial variations occurring in the Earth’s atmosphere composition. The ERA-Interim dataset (ECMWFmeteorological organization) was used to retrieve vertical atmospheric conditions, which are available at a global scale with a resolution of 0.125 degrees and a temporal resolution of 6 h. A temporal and spatial linear interpolation of meteorological variables was performed to match the acquisition dates and coordinates of the Landsat images. The atmospheric correction parameters were then estimated on the basis of this reconstructed atmospheric grid using the commercial RTMsoftware MODTRAN. The needed surface emissivity was derived from the common vegetation index NDVI, obtained from the red and near-infrared (NIR) bands of the same Landsat image. This permitted an estimation of LST for the entire image without degradation in resolution. The software tool, named LANDARTs, which stands for Landsat automatic retrieval of surface temperatures, is fully automatic and coded in the programming language Python. In the present paper, LANDARTs was used for the local and spatial validation of surface temperature obtained from a Landsat dataset covering two climatically contrasting zones: southwestern France and central Tunisia. Long-term datasets of in situ surface temperature measurements for both locations were compared to corresponding Landsat LST data. This temporal comparison yielded RMSE values ranging from 1.84 ° C–2.55 ° C. Landsat surface temperature data obtained with LANDARTs were then spatially compared using the ASTER data products of kinetic surface temperatures (AST08) for both geographical zones. This comparison yielded a satisfactory RMSE of about 2.55 ° C. Finally, a sensitivity analysis for the effect of spatial validation on the LST correction process showed a variability of up to 2 ° C for an entire Landsat image, confirming that the proposed spatial approach improved the accuracy of Landsat LST estimations.

Graphical Abstract">

Graphical Abstract

1. Introduction

Land surface temperature (LST) and land surface emissivity (LSE) are two key parameters used in many environmental studies because they are closely connected to the Earth’s surface energy balance. At a regional scale, land and sea surface temperatures (LST and SST, respectively) can be estimated with the help of thermal infrared (TIR) measurements recorded by remote sensors onboard spaceborne platforms. These estimates are used for hydrological and meteorological purposes. In particular, they provide valuable information to monitor evapotranspiration and determine the soil water status on agricultural lands [1,2,3]. To study the oceans, SST data provide information on sea warming, potential evaporation and water circulation [4]. However, atmospheric effects of absorption and emission affecting the thermal spectrum must first be removed in order to use thermal imagery in absolute temperature research application [5]. At present, very few surface temperature data products obtained through remote sensing observations are available. The MODIS products distributed by the U.S. National Aeronautics and Space Administration (NASA) are the best known. The most commonly-used global and daily land surface temperature and emissivity product is MOD11A1 [6]. MOD11A1 (Collection 5), which uses a grid spatial resolution of 1000 m, is based on the generalized split-window LST algorithm applied to the MODIS thermal multispectral bands [7]. Note that recent works proposed an operational Kalman filter strategy applied to the 3 IR SEVIRI (Spring Enhanced Visible and Infrared Imager) channels from geostationary platform plus ECMWF meteorological analysis (used also here) to directly retrieve surface temperature and the three channels emissivities [8,9]. Spatial resolution of products is around 3 km but earth is fully scanned every 15 min, and the retrieved surface emissivity and temperature accuracy is of ±0.005 and ±0.2 ° C, respectively. nary platform plus ECMWF meteorological analysis (used also here) to directly retrieve surface temperature and the three channels emissivities [8,9]. Spatial resolution of products is around 3 km but earth is fully scanned every 15 min, and the retrieved surface emissivity and temperature accuracy is of ±0.005 and ±0.2 ° C, respectively. Proposed by the Ministry of Economy, Trade and Industry (METI) of Japan in collaboration with NASA, the Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) kinetic temperature and emissivity products [10] offer a higher spatial resolution of 90 m, but with a low frequency of revisits; indeed data are potentially available only one to eight times per year.
Landsat satellites have been continuously recording multispectral images of Earth since the early 1970s. The U.S. Geological Survey (USGS) Earth Resources Observation and Science (EROS) Center provides these data free of charge. In fact, from the first Landsat 1 (launched in 1972) to the current Landsat 8 (launched in 2013), the project has made available more than 40 years of Earth land surface data to study land processes. In total, this represents more than one million scenes, each one offering at least one thermal band. Landsat 5 and 7 provide one thermal channel; Landsat 8 provides two, although it is recommended to use only one band [11,12].
In this sense, the Landsat program appears to be an interesting source of LST records. However, land and sea surface temperatures have not been made available in a free-to-use format. Nevertheless, several methods have been proposed to estimate LST from Landsat’s single TIR channel [13,14,15] by using atmospheric profiles obtained from databases or through radio sounding. Several studies [16,17] have validated the application of such methods. For example, Barsi [18] proposed a web-based correction tool that used modeled atmospheric global profiles obtained from the National Center for Environmental Prediction (NCEP) in conjunction with the commercially available radiative transfer model MODTRAN [19]. This process allowed the authors to obtain the necessary atmospheric transmission, upwelling, and downwelling radiance values to convert sensor radiance data into surface brightness temperatures (Tb). A subsequent study [20] used the same process to successfully reproduce the results with the North American Regional Reanalysis (NARR) database by using atmospheric profiles as inputs. Note that in both approaches, the computed atmospheric parameters are representative of a region of interest defined by the spatial resolution of the atmospheric profiles. Moreover, in both cases, raw user manipulations are necessary to estimate brightness temperatures, emissivities, and LST from Landsat image radiances. In most studies that used the Landsat thermal channel, LST were corrected for local atmospheric parameter values obtained from [21], or through radio sounding in conjunction with a radiative transfer model such as MODTRAN, or with a single channel method [22]. As a result, the correction process is more adapted to local than to regional applications. For the latter, radio sounding could prove to be inappropriate because it may be either non-coincidental or too far away from the image footprint, leading to errors such as discussed by Mira et al. [23]. The present study expands on Barsi’s approach and proposes an independent and open source software tool that uses the European ERA-Interim global atmospheric reanalysis database in conjunction with the MODTRAN software package. The proposed software tool, named LANDARTs, which stands for Landsat automatic retrieval of surface temperature, is able to automatically estimate surface brightness temperatures and LST for an entire Landsat thermal image with a coherent spatial correction.
Estimating LST through TIR radiance measurements obtained with remote sensors usually requires four steps of correction [5]: (1) converting spectral radiances into at-sensor brightness temperatures, which is typically done by inverting Planck’s law adapted to the spectral sensor window as described by Equation (4); (2) correcting for atmospheric absorption and re-emission, which requires knowledge of the atmospheric profile composition, especially concerning the water content and temperature; (3) correcting for surface emissivity (estimated if unknown); and (4) correcting for surface roughness or topography. This last correction is rarely taken into account at resolutions of 100 m and is neglected at the nadir view angle. Neglecting atmospheric correction results in systematic errors that are proportional to the target temperature for a given atmosphere.
Various approaches have been developed to solve or simplify the atmospheric correction process for TIR data dedicated to single-channel sensors. The available methods can be classified into four categories: the mono-window algorithm, which requires air temperature data only [17]; the single-channel algorithm, which requires total atmospheric water content values only [13]; the look-up tables approach, which supposes an exhaustive learning database [14]; and the physical approach, which requires the calculation of atmospheric parameters with a radiative transfer model and demands detailed knowledge of the atmospheric profile [15,18]. The main goal of our paper is to present and evaluate an automatic Python tool that relies on the MODTRAN model and on ERA-Interim data. The tool was designed to perform atmospheric and surface correction of Landsat (TM and TIRS) thermal data sensors and to retrieve LST and LSE values at a resolution of 30 m (the resolution of radiances provided by the USGS). First, we present the theory of atmospheric correction and the applicable processing algorithm; We then describe the process of comparing the LST reference product from ASTER with in situ temperature measurements to validate the tool’s performance; Finally, we discuss the sensitivity of the results to spatial variations in the atmospheric correction parameters.

2. Materials and Methods

2.1. Atmospheric Correction

The spectral distribution of the radiation from a blackbody is described by Planck’s law:
B λ ( T ) = C 1 λ 5 [ e x p ( C 2 λ T ) 1 ]
where B λ ( T ) is the spectral radiance ( W · m 2 · sr 1 · μ m 1 ) of a blackbody at temperature T (K) and wavelength λ ( μ m 1 ); C 1 and C 2 are two physical constants ( C 1 = 1.191 × 10 8   W · μ m 4 · sr 1 · m 2 , C 2 = 1 . 439 × 10 4 μ m·K). However, land surfaces measured by remote sensing are gray bodies with a surface emissivity ( ϵ s < 1 ). Therefore, radiance measured by a sensor of a gray body without atmospheric perturbation can be expressed as:
L λ , s e n , n o a t m o ( T ) = ϵ s , λ B λ ( T )
where L λ , s e n , n o a t m o ( T ) ( W · m 2 · sr 1 · μ m 1 ) is the spectral radiance of a gray body at temperature T (K) measured by a sensor without atmosphere, ϵ s , λ is the surface spectral emissivity at wavelength λ and B λ ( T ) ( W · m 2 · sr 1 · μ m 1 ) is the spectral radiance of the equivalent black-body at temperature T (K) as detailed in Equation (2).
When processing the thermal data of sensor such Landsat (Band 6 on L5 and L7, band 10 on L8; see Table 1), we need to take into account the atmospheric interaction and use a radiative transfer model to remove the atmospheric effects. In this goal we need to solve the radiative transfer Equation (3) in the band wavelength window of the sensor. A radiative transfer model can be used to estimate the unknowns terms provided the corresponding vertical atmospheric composition is known.
L λ , s e n ( T ) = τ λ ϵ s , λ B λ ( T ) + L u , λ , a t m o + τ λ ( 1 ϵ s , λ ) L d , λ , a t m o
where:
  • L λ , s e n = the spectral at-sensor radiance (top of the atmosphere) ( W · m 2 · sr 1 · μ m 1 )
  • B λ = the radiance of a supposed blackbody surface target at a kinetic temperature T in (K)
  • τ λ = the atmospheric transmittance (unitless) at the wavelength λ ( μ m)
  • L u , λ , a t m o = the upwelling atmospheric radiance in the wavelength window ( W · m 2 · sr 1 · μ m 1 )
  • L d , λ , a t m o = the downwelling atmospheric radiance in the wavelength window ( W · m 2 · sr 1 · μ m 1 )
  • ϵ s , λ = the surface spectral emissivity (unitless).
As a result, the derived spectral atmospheric parameters (τ, L u and L d ) and the surface emissivity ( ϵ s ) permit the conversion of the at-sensor radiance into surface radiance exempt from atmospheric effects. The obtained surface radiance can then be converted into surface kinetic temperatures (LST) by inverting a simplified Planck’s law from Equation (1) adapted to specific sensor spectral window by spectral integration according to [24]:
L S T = K 2 ln ( 1 + K 1 B ( T ) )
where:
  • L S T = equivalent to T, the land surface temperature (K)
  • B ( T ) = the surface radiance ( W · m 2 · sr 1 · μ m 1 )
  • K 1 = the pre-launch calibration constant 1 ( W · m 2 · sr 1 · μ m 1 )
  • K 2 = the pre-launch calibration constant 2 (K)
  • K 1 and K 2 are provided by the USGS for L5 and L7 [24], and for L8 [25].
Among these key parameters, the transmittance is known to be inversely proportional to the downwelling and upwelling radiance (see the second and third figures of Section 4.3). The evolution of these atmospheric parameters is primarily controlled by the total water vapor content in the vertical atmospheric column of air [17] and secondary by air temperature profile. The total water vapor content changes from hour to hour, but follows a more or less annual trend depending on latitude. This trend generally attains its maximum in summer (warm and humid atmosphere) and its minimum in winter (cold and dry atmosphere). This atmospheric impact on the surface temperature estimation, that we can associate to the difference between top of atmosphere (TOA) and top of canopy (TOC) brightness temperature (respectively T b , T O A and T b , T O C ), is proportional to the target temperature. The warmer the surface target temperature is, the higher the atmospheric impact, therefore the correction. Below a temperature limit, around the atmospheric temperature, the mean lower layer atmospheric effect is inversed and T b , T O A became higher than T b , T O C . For the two extremes: hot and cold targets, suppose brightness temperature equivalent to LST correspond to under and over estimation of LST, respectively.
The atmospheric transmittance is a very sensitive parameter for LST estimation. Qin et al. [17] shows that an error of 0.05 on transmittance value generated an error comprised between 1–4 ° C on LST estimation. In addition to daily and annual variations, the existence of spatial variability in the atmospheric composition, which results in spatial variation of the atmospheric correction parameters, can cause large errors if a constant value per image is used.

2.2. Radiative Transfer Modeling with MODTRAN

MODTRAN (moderate resolution atmospheric transmission) is a commercial atmospheric radiative transfer model developed by the U.S. Air Force [19,26]. The model covers a spectrum of 0.2–100 μ m that is from the mid-ultraviolet to the far infrared. MODTRAN can be used to estimate the three atmospheric correction parameters in the Landsat thermal spectrum. Because MODTRAN cannot separate downward contributions from upward contributions in a single run, the model was run twice in nadir observation configuration (according to Landsat acquisition mode) using the same atmospheric profile with the following two configurations:
  • The upwelling radiance (i.e., the atmospheric effect between the surface target and the sensor) and the atmospheric transmittance between the surface and the sensor were estimated through simulation by setting the sensor location to 100 km above the surface (considered the sensor altitude) and setting both the surface albedo and emissivity to zero, corresponding to a complete lack of surface reflection for the entire spectrum.
  • The downwelling radiance was estimated through a second run using a configuration in which the sensor was located at 1 m above the surface and the surface albedo was set to 1. For this configuration, we assumed that the downwelling radiance was completely reflected upward toward the sensor.
As the MODTRAN outputs are provided following model spectral resolution (0.1 cm 1 ), an integration was performed between the bounds defined by the sensor thermal spectrum response .
v a r ( λ i ) = λ i , m i n λ i , m a x v a r ( λ ) R s ( λ ) d λ λ i , m i n λ i , m a x R s ( λ ) d λ
where R s is the spectral response (unitless) of the sensor at center wavelength λ i ( μ m) of a band with spectral window from λ i , m i n to λ i , m a x . v a r ( i ) alternatively denotes the transmittance τ, the upwelling radiance L u , or the downwelling radiance L d value extracted from the corresponding column of the MODTRAN output files, between the lower λ i , m i n and upper λ i , m a x wavelengths.

2.3. Atmospheric Profiles Construction

The ECMWF (European Centre for Medium-Range Weather Forecasts) provides full access to several atmospheric reanalysis datasets. In accordance with our requirements, we selected the ERA-Interim global atmospheric reanalysis [27], which provides a variety of vertical and surface atmospheric fields from 1979 to the present (accessible with a 2-month delay) in uniform latitude/longitude grids ( 0 . 125 ° , 0 . 25 ° , 0 . 75 ° , 1 ° , 1 . 125 ° , 1 . 5 ° , 2 ° , 2 . 5 ° and 3 ° ). The parameters are interpolated from the original N128 reduced Gaussian grid using bilinear methods. LANDARTs integrates the Python API provided by ECMWF to process automatic requests sent to the database [28]. For the Landsat LST data atmospheric correction process, each request sent to the ERA-Interim dataset is defined to cover a footprint that is larger than the image. The returned query results are 2D or 3D matrices of surface and vertical atmospheric fields with a default horizontal grid resolution of 0.125 ° , which represents a grid of about 18 × 28 = 504 nodes for each Landsat image. Vertical fields are provided for 37 atmospheric pressure levels from 1–1000 hPa, which corresponds to approximate altitudes of between 100 m and 35 km, depending on the air temperature. In addition to pressure, we were interested in the geopotential, the temperature, and the relative humidity. The altitude of each level can be determined correctly by dividing the geopotential by the gravity constant ( g = 9.81 s 2 ), the air temperature is given in Kelvin, and the relative humidity is given in percent. Analogous to Barsi’s atmospheric correction parameter calculator (Atmocorr, http://atmcorr.gsfc.nasa.gov/), fields from the highest atmospheric levels (for altitudes between 35 km and 100 km) were obtained from the U.S. standard atmosphere of 1976 [29]. For the northern hemisphere, the atmosphere’s summer mid-latitude profile was used for images acquired between March and September and the winter profile was used for all other images; the inverse date ranges were used for the southern hemisphere. At this step, all the vertical profiles for each of the ERA-Interim fields of interest were described using 1000 hPa as the starting value, regardless of the landform or elevation. The surface elevation and air temperature were provided by the surface field data product. However, relative humidity values were not available in the surface field, we estimated it from 2 m dew-point temperature and the surface pressure according to Equation (6) from [30]:
R H = 100 e ( t , P ) e s ( t , P ) , with e s ( t d ) = e ( t )
where R H (%) is the relative humidity, t the air temperature, t d the dew-point air temperature, P the atmospheric pressure, e s the air saturation vapor pressure and e the air vapor pressure estimated with Magnus formulae [31]. The two datasets were merged to create a consistent meteorological vertical altitude profile extending from the Earth’s surface to the top of the atmosphere. In practice, for surface field pressures higher than or equal to 1000 hPa (first vertical level), surface fields were added at the first level of the vertical profile; for surface field pressures less than 1000 hPa (occurring mainly because of topographic features), surface fields were added to vertical profiles according to the pressure value, and all field profiles below these levels were cut off. For each case, pressure and altitude values were verified to ensure a strictly decreasing order of the vertical profile. However, the values of the vertical profile were not smoothed through interpolation. In fact, failure to respect this decreasing order created a conflict in the MODTRAN execution. Finally, the number of levels in the vertical profile may change depending on the surface elevation. The resulting vertical profile was then assigned to the MODTRAN input file. MODTRAN auto-completes the standard atmospheric profile for standard gases. The ERA-Interim profiles were adapted for each point of the requested grid.

2.4. Emissivity Estimation

Land surface emissivity (LSE) is a key parameter that describes the radiative absorption power of a surface in the longwave radiation spectrum. LSE depends on the target surface top layer composition, such as presence of soil, soil type, vegetation and density, or roughness of the surface [32,33]. Although LSE data are indispensable, they have thus far proven difficult to estimate with remote sensing applications. The sensitivity analysis performed by Olioso [34] showed that an error of ±0.01 for emissivity results in LST values with an error between 0.6 ° C and 0.9 ° C (depending on the LST). Qin et al. [17] reported the same range of values. The availability of emissivity data suitable for LST corrections also depends on the sensor spectral windows and on the viewing angles [35]. Landsat spectral windows are given in Table 1. Data were acquired at the nadir, thus avoiding the need for angular corrections. Several solutions exist to retrieve emissivity from satellite data [33]: they may be based (1) on the land cover classification and emissivity coverage tables; (2) on the exploitation of multi-spectral channels, as is the case for temperature and emissivity separation (TES) methods; or (3) on empirical relationships involving vegetation indices, which is the most commonly used method. Li et al. [36] compared six emissivity extraction methods and concluded that all six methods yielded similar results. The first method required a classification of the zone of interest, which is not always easy to obtain; The second method required at least three thermal channels, which are unavailable for Landsat 5 and 7, and the USGS recommends using only thermal band 10 for L8. Therefore, we ultimately decided to rely on the vegetation index method. Most vegetation index methods involve the normalized difference vegetation index (NDVI). The NDVI is defined as:
N D V I = N I R R E D N I R + R E D
where NIR represents the reflectance in the near infrared region, and RED is the spectral reflectance in the red region. A commonly used method for calculating mixed pixel emissivity is the estimation proposed by Sobrino et al. [16], adapted from [37]. This empirical relationship is based on an exponential function that is NDVI dependent and defined as follows:
ϵ λ = ϵ v λ ( ϵ v λ ϵ s λ ) N D V I N D V I v N D V I s N D V I v k
where λ represents the spectral band, ϵ v is the vegetation emissivity (set to 0.99 in our case), and ϵ s is the soil emissivity, which may be adapted by user for more accuracy. N D V I v is the maximum N D V I for full vegetation (set to 0.99) and N D V I s is the minimum N D V I for bare soil (set to 0.17). The soil emissivity was not set because its value depends on the soil’s composition and roughness. In our case, the soil emissivity was set to 0.96 in agreement with the mean emissivity value extracted from the JPL library database [38] corresponding to the soil composition of our selected sites. The k parameter was fixed arbitrary to 2, but note that [39] proposed explore values between 2 and 3. This parameter can be adjusted in the future version. Another difficulty in determining soil emissivity is the non-homogeneity of the soil at the pixel scale. However, for acquisitions at a resolution of 100 m, Zhang et al. [40] found that the pixels can be considered homogeneous.

2.5. Ground Measurements

The ground measurements used for the local validation of the LST estimates were retrieved from two meteorological network stations one local near Kairouan in central Tunisia ( 35 ° 18 17.14 N, 9 ° 54 56.62 E) the other located near Toulouse in southwestern France and composed of two plots ( 43 ° 54 97 N, 01 ° 10 61 E and 43 ° 49 65 N, 01 ° 23 79 E) identify hereafter as respectively Aurade and Lamasquere. A localization map for both networks is given in Figure 1. The Tunisian site is representative of a semi-arid climate, with high surface temperatures during summer. The French site is representative of a continental climate with moderate summer surface temperatures.
The Tunisian site is part of the SICMED program (http://www.sicmed.net/). The measurements considered in the present study were acquired on non-irrigated olive tree fields composed of trees spaced at regular intervals of 20 m and with crown diameters of about 5 m. The ground was a sandy bare soil. The French site was composed of two cultivated plots located at a distance of 12 km from each other and with different soils and crop rotations. Aurade is a crop rotation of wheat/rapeseed and sunflower each 5 years; Lamasquere is a rotation of wheat and silage maize. The measurement methodology for both sites follows ICOS (integrated carbon observation system) recommendations. A complete description of the agricultural practices and field instrumentation for the two French plots are available in [41,42]. For both zones of interest, two kinds of micro-meteorological data were measured and processed to estimate local LST. First, directional remote temperatures were obtained with an IR-120 thermo-radiometer (Campbell Scientific) located at 2 m height and nadir oriented with a spectral window of 8–14 μ m; Second, incoming and outgoing long-wave radiation was measured with a hemispherical CNR1-Pyranometer (Kipp Zonen) with a spectral window of 5–50 μ m. Surface brightness temperatures were measured directly with the thermo-radiometer and were then compared with Landsat LST estimations after emissivity correction. For the hemispherical pyranometer measurements, the surface brightness temperature was first estimated with the Stefan-Boltzmann law:
T b λ 1 λ 2 = L W λ 1 λ 2 σ 1 4
where T b λ 1 λ 2 is the surface brightness temperature, and L W λ 1 λ 2 is the outgoing long-wave radiation for the spectral window λ 1 λ 2 equal to 5–50 μ m.
Then, the surface brightness temperature was converted into surface temperature using the approach proposed in [34]:
L S T T b λ + ( 1 ϵ ) 4 ϵ T b λ ( 1 ϵ ) L W λ 4 ϵ f λ ( T b λ ) σ T b λ 3
where λ represents the ( λ 1 λ 2 ) band spectral window, σ is the Stefan-Boltzmann constant, L W λ the incoming long-wave radiation, f λ ( T b λ ) is a function factor corresponding to the fraction of energy emitted in the λ spectral window by a black body at temperature T relative to the energy emitted over the full spectrum as defined in [43] and approximate to unity in our configuration, and ϵ is the soil emissivity.
The data used in this study are acquired between 2011 and 2013 for France and between 2012 and 2014 for Tunisia. The available data allowed us to compare more than three years of local thermal measurements with our estimates obtained from Landsat cloud-free images.

2.6. Landsat Data

Landsat data are the main inputs for LANDARTs. The present study used data retrieved by three missions: Landsat 5 TM (L5), Landsat 7 ETM+ (L7), and Landsat 8 OLI-TIRS (L8). Both L5 and L7 provide a single-band in the thermal domain (band 6). Landsat 8 provides two thermal bands (10 and 11). However, the USGS recommends against the use of band 11 because considerable uncertainty persists regarding the quality of the acquired data [11,12]. The thermal band characteristics for the three sensors are summarized in Table 1. Image data for the studied area and for the period of interest were downloaded from the USGS website (http://glovis.usgs.gov/). Note that all Landsat thermal data provided by the USGS were re-sampled to a resolution of 30 m (Table 1). In LANDARTs data are read as digital numbers and converted into sensor radiances ( L s e n ) in accordance with the calibration coefficients given in [24] for Landsat 5 and 7 and in [25] for Landsat 8.
The validation of our atmospheric correction tool considered two main areas corresponding to the ground dataset measurements (Section 2.5). The first area was located in south-western France (near Toulouse) at the intersection of two tiles: 199/030 for the western portion, and 198/030 for the eastern portion, in accordance with USGS WRS-2 nomenclature. The second area was located in Tunisia (near Kairouan) on tile 191/035. The downloaded Landsat images for France were acquired between 1 January 2011 and 31 December 2013. For the Tunisia site the images were acquired between 5 April 2012 and 31 December 2014.

2.7. ASTER Data

ASTER thermal channels constitute a valuable dataset for our investigation because both channel 13 and 14 (Table 2) are similar to the single TIR channel of L5 and L7 (Table 1). The L8 channel 10 also covers a large part of ASTER channel 13. It is therefore possible to compare Landsat and ASTER products for certain ranges of wavelengths. The ASTER database provides several data products. However, our study only had access to one higher-level product: the surface kinetic temperature called AST_08, which is the result of the TES algorithm applied to the five thermal channels at a resolution of 90 m [44].
All cloudy images were discarded to prevent problems with automatic cloud detection. A comparison between the ASTER dataset covering southwestern France for the period from 2000–2006 and the USGS database revealed an acquisition time overlap of only a few dates. All dates after 2003 were discarded, to avoid large bands of missing data on the L7 images due to SLC failure and unfortunately corresponding to the ASTER position in the Landsat images. The ASTER dataset for Tunisia consisted of only two images, one of which was cloudy. Incidentally, an L8 image was available for the only cloud-free ASTER image date. The data selected and used in the present study are presented in Section 4.2. It was our aim to implement identical procedures to evaluate the performance of the proposed tool for both sites; we consequently considered only the surface kinetic temperature product (the only product available for the Tunisian site).

3. LANDARTs

3.1. Main Algorithm

LANDARTs is a Python tool that works (currently) only with Landsat archives in the new metadata file format (data acquired after 2012) that can be downloaded from the USGS website. Although the tool accepts a custom installation of MODTRAN, we tested the following three versions only: 4.1 on Unix, 4.3 on Windows, and 5.2 on Unix and Windows operating systems. A flowchart for the proposed tool is presented in Figure 2. The implementation process consists of three steps to produce LST maps at a resolution of 30 m using Landsat 5, 7 or 8 images. Step 1 consists of calculating the three atmospheric correction parameters. The details of the parameter calculations are described below:
  • Acquisition dates and times are read automatically according to the metadata file provided by the USGS data. The four corners of the Landsat footprint are also extracted from GeoTiff metadata.
  • The information obtained in (1) are used to create four spatial queries to the ECMWF server: for each two time samples bounding the time input, two queries are done, one for the surface and one for the pressure level dataset. The query zone is defined to be larger than the Landsat footprint to allow interpolations at the edges. Data are in the uniform latitude/longitude grids in one of the resolution proposed by ERA-Interim and selected by user.
  • As explained in the ECMWF Section 2.3, the two time profiles constructed with surface and pressure levels for the height from the Earth’s surface to an elevation of 100 km are linearly interpolated between the two times samples to give acquisition time [18]. This interpolation calculus is performed for each point of the ECMWF grid to obtain a 3D matrix (2D spatial and 1D vertical) of an equivalent atmospheric variables profiles at the acquisition time.
  • In accordance with the MODTRAN documentation [19], the 10 required input cards are written using the atmospheric profiles obtained in (3). The cards are used to create two tape5 files as MODTRAN inputs. This process is repeated for each point of the ECMWF grid.
  • MODTRAN run on each point of the grid defined by atmospheric profiles. When done, the MODTRAN predicted per-wavelength transmission, upwelling radiance and downwelling radiance are integrated over the instrument’s spectral response. According to the metadata, the 2D matrices of the three parameters are finally interpolated using the gdalwarp utility from the open source Geospatial Data Abstraction Library (GDAL, http://www. gdal.org/gdalwarp.html) to create three GeoTiff images oversampled at the resolution of 30 m. These three final parameter’s images can then be used for pixel-by-pixel processing.
Step 2 consists in converting the digital numbers of multispectral data into reflectance using the metadata and the USGS formulas. Next, Red and NIR bands are used to calculate the NDVI. The NDVI map were then converted into an emissivity map at the data resolution of 30 m according to Equation (8).
Finally, Step 3 handles the thermal data: in accordance with the USGS user guide, the digital numbers of the thermal band are converted into top of atmosphere radiance using the metadata file. The emissivity estimates derived from the NDVI in Step 2 are then used to calculate the ground radiance, and LST are converted by inverting radiative transfer Equation (3) with spatialized atmospheric parameters and simplified Planck’s law (4). The LANDARTs tool is available free of charge by contacting the corresponding author. However, note that the user needs to have his or her own official instance of MODTRAN.

3.2. Computing Requirements and Efficiency

For the present study, LANDARTs was installed on a Unix server with the following configuration:
  • Ubuntu 12.04.2 (x86_64)
  • 32 Gb of memory
  • 16 Processors (Intel(R) Xeon(R) @2.67GHz)
  • Python, Version 2.3
  • MODTRAN, Versions 4.3 and 5.1
The tool was benchmark tested using 140 tiles with only one processor, and the mean processing time was calculated. With the configuration given above, processing a single image took about 5 min 30 s. This total time is divided into three segments: the ECMWF request required 20 s to process, the MODTRAN computation took around 3 min 10 s, and final image processing took around 2 min. Obtained by benchmarking against 80 processed images, these times represent mean values, which can be adjusted by changing the request time and the hardware configuration. The processing of a single image required 3.2 GB of RAM.

4. Results and Discussion

We performed a local as well as a spatial validation of the LST estimation results obtained with the LANDARTs. For the local validation, we used a large set of 140 images acquired by L7 and L8 and obtained from the USGS database archive (Section 2.6) for the two considered regions-southwestern France (southwest of Toulouse) and central Tunisia (west of Kairouan). Atmospheric and surface corrections were performed with the LANDARTs. The obtained Landsat LST were filtered for clouds and compared with in situ surface temperature measurements for the local validation (Section 2.5). We took advantage of two ASTER data from archive (Section 2.7) acquired in quite same time than an L7 (SLC-on) and an L8 image to evaluate the performance of LANDARTs correction capabilities at the spatial scale. Finally, we performed a sensitivity analysis of the impact of atmospheric parameters on the correction process. All results discussed in this section were processed using the LANDARTs presented in Section 3.1 with version 4.3 of MODTRAN for UNIX.

4.1. In Situ Validation

The LST values derived from Landsat data and corrected with the LANDARTs were compared with in situ surface temperatures obtained from the CESBIO Laboratory ground networks (see description in Section 2.5) for both the French (near Toulouse) and the Tunisian (near Kairouan) site. To this end, all Landsat data corresponding to the time window of the ground measurements were downloaded from the USGS. A total of 80 Landsat images for southwestern France and 60 images for Tunisia were available in the archive and were subsequently processed with the MACCS Level 2A cloud detection processor from Hagolle et al. [45]. Finally, after cloud mask filtering and the exclusion of SLC-failure no-data pixels for Landsat 7, we obtained two sets of clear-sky images, 27 for France and 38 for Tunisia, with valid LST representative of the measurement sites. LST retrieved with LANDARTs were extracted to a 5 by 5 pixel window, and the mean value was computed. The 5 by 5 window dimension was chosen in order to cover more than one pixel at the acquisition resolution (100 m). On the basis of the values generated by the tool, emissivity data were also extracted and averaged to compute the in situ measurements (see Section 2.5). The comparison scatter plots for the Tunisian olive-grove site and the two French crop sites (Lamasquere and Aurade combined on the same graph) are shown in Figure 3 and Figure 4, respectively.
For the French site, the comparison with valid surface temperatures (Figure 4) yielded an RMSE of 2.55 ° C and a mean bias of −1.77 ° C with a correlation coefficient of 0.98. These values are consistent with the ground validation of spatially derived temperature products proposed by Zhou et al. [46], who compared the ASTER product with ground temperatures and obtained an RMSE close to 2.54 ° C after atmospheric correction. For the Tunisian site, using the same procedure, we obtained a global RMSE of 1.84 ° C with a positive bias of 1.24 ° C and a correlation coefficient of 0.99. Overall, this local comparison gave good results that are consistent with other work in the literature [15,23,46]. The superior quality of the results obtained for the Tunisian site is likely attributable to the site’s relatively static ground management. In fact, the French experimental sites were ploughed regularly; crop rotations and growing crops produced considerable surface modifications and resulted in a change of vegetation, which affected measurement precision. In contrast, the ground of the Tunisian site was mainly composed of sand that remained untouched over time and was planted with slow-growing olive trees. Srivastava [47] has shown that a good match of ground LST and satellite data can be obtained for homogeneous soil, whereas an error of ±2 ° C was found for less homogeneous surfaces such as are often encountered in agricultural areas.

4.2. Spatial Validation

First, estimated surface temperatures were corrected for spatial atmospheric and surface emissivity with the LANDARTs. Two ASTER surface kinetic temperature products (AST08 from Level 1B), whose acquisitions matched the two Landsat acquisitions, were selected to evaluate the spatial coherence of the resulting estimates (Table 3). Images for the comparison were selected according to the following criteria: availability of the ASTER dataset in the archive, no clouds on the scene, identical dates, and availability of images from both currently active Landsat sensors, L7 and L8. For the sake of simplicity, we excluded SLC-failure data for L7. In fact, the thermal and emissivity data products proposed through ASTER’s TES algorithm [48] can be considered a reference. We compared the two surface temperature products after resampling the data at a resolution of 270 m, that is a multiple of 30 m (Landsat) and 90 m (ASTER), to avoid alignment problems and to cover more than four pixels. In contrast to the kinetic temperature product from ASTER, which was resampled directly, the Landsat radiances and reflectances were first warped into an ASTER footprint o the resolution of 270 m before being processed by LANDARTs.
Figure 5 and Figure 6 show density scatter plots for the pixel-by-pixel comparison between the Landsat LST product and ASTER kinetic temperatures at a resolution of 270 m for the Tunisian and French site, respectively. For both sites, a good correlation was obtained with a high density near the regression line. In both figures, a positive offset indicates that ASTER temperatures were warmer on average, which is consistent with ASTER’s comparatively delayed acquisition. We determined an R M S E of 2.5 ° C for the Tunisian site and an R M S E of 3.52 ° C for the French site. Indeed, ASTER acquisitions are 14 min and 30 min later than Landsat for the Tunisian and the French location respectively. Unfortunately, no in-field measurements were available on 2003 for the French site to discuss this offset. Nevertheless, on the Tunisian site, according to the thermo-radiometer temperature measurement, we observed an increase of between 1 ° C (on tree canopy) to 1.5 ° C (on bare soil) between 10:00 and 10:30. As a result, we can estimate that 0.5 ° C–0.7 ° C of the offset observed in Figure 5 can be associated to the difference due to acquisition time for the Tunisian dataset. An offset of the same order of magnitude must be associated with the French datasets, whether between one quarter and one third of the mean bias. For their comparison between ASTER, MODIS, Landsat 7, and ground measurements, Srivastava et al. [47] established a similar mean error with a value of 2.3 ° C. They also determined an R M S E of around 3 ° C for the differences between Landsat and the two other datasets. Zhou et al. [46] compared ground measurements with data from ASTER bands 13 and 14, which were corrected separately using several algorithms. For the present study, we used the atmospheric correction algorithm ( A C ) with estimated emissivity and modeled atmospheric parameters. Zhou et al. [46] determined an R M S E of around 2.2 ° C for their comparison of atmospherically corrected ASTER band 13 and in situ temperature measurements.
Figure 7, Figure 8, Figure 9 and Figure 10 show disparity maps for the ASTER surface kinetic temperature product and the surface temperatures obtained with LANDARTs. The Tunisian maps show considerable homogeneity; most variations occurred within a ±2.5 ° C temperature window. Surface variations essentially occur because of differences in emissivity estimations. The comparison reveals a large amount of noise generated by the pixel-by-pixel methods, which is observable for the ASTER products, in particular on the French comparison map (Figure 9) and also on the emissivity map in Figure 11C. These figures also provide an explanation for the difference between the R M S E values of both sites: the noise level was much higher on the French image than on the Tunisian one.

4.3. The Importance of Spatial Corrections

The following section examines the importance of using a set of spatialized atmospheric parameters for the atmospheric correction process. As illustrated by the example in Figure 11A,B, the three atmospheric correction parameters (τ, L u , L d ) may vary for any particular Landsat image acquired on a single date with one atmospheric profile. The LST variation depicted in Figure 11D, which is representative of most of the observed cases, shows the possible amplitude of atmospheric variability for a complete image. The same dataset can therefore be used for the in situ validation to retrieve the variability of the atmospheric parameters over time. To demonstrate the spatial impact of the atmosphere at the image level, we compared four correction approach. For each free cloud acquisition date, we potentially disposed of a set of spatialized atmospheric parameters for the whole image. To evaluate the variability of parameters in image thus correction, we consider the set of parameters either
  • obtained at the center of the image
  • corresponding to the set with the minimum value of transmittance, and associated atmospheric radiances in image
  • corresponding to the set with the maximum value of transmittance, and associated atmospheric radiances in image
  • or spatiality distributed as applied in the LANDARTs tool.
Atmospheric correction was then performed on the entire image according to the four approaches on parameters application and the emissivity map shown in Figure 11C, thereby allowing us to evaluate the effect of the parameter choice on the correction process.
Figure 12 and Figure 13 respectively shows the atmospheric transmittance and the upwelling atmospheric radiance for the 27 Landsat images used in the local validation of the Aurade site. Presumably, it was a consequential decision whether or not to apply spatialized parameters because the site was located in a corner of the Landsat tile. The transmittance was extracted for an area of 5 by 5 pixels; the mean value is plotted in red. Note that the upwelling and downwelling radiances correlated negatively with the atmospheric transmittance, as shown in Figure 12 and Figure 13.
The existence of such variations within a Landsat footprint indicates that the center-scene parameters may not be representative of the entire scene. As shown in Figure 12, this variation can attain amplitudes of transmittance up to 0.2 between the minimum and maximum, which has considerable impact on the LST. A mean difference of around 0.04 for the transmittance and 0.31 for the upwelling radiance was found between the center scene and the spatialized parameters extracted from the Aurade site. These mean values were calculated for the dataset shown in Figure 12 and Figure 13. A mean difference of 0.11 for the transmittance and 0.83 for the upwelling radiance was found between the minimum and the maximum of the Landsat footprint for this dataset. Qin et al. [17] have shown that an error of 0.05 on the transmittance value leads to an error of 1–4 ° C on the LST value. Consequently, performing atmospheric correction with local parameters induces certain errors in the final outcome and decreases the accuracy of the LST map for a large satellite image. In some cases, the site and the center scene may incidentally either exhibit very similar parameters, or the center-scene parameters may be the local minimum or maximum of the entire footprint. The latter case would result in either an under- or overestimated LST. To demonstrate this, we modified the procedure described previously by adding a step to correct the radiances with the different sets of parameters. The same procedure as described in Section 4.1 was used for the validation, but it was repeated four times. Figure 14 and Figure 15 shows the same process as for the validation, but with the additional step of correcting the Landsat images for the center-scene parameters and the minimum and maximum parameters. The results obtained for the Tunisian olive grove site are shown in Figure 14. The green area in this figure shows the amplitude of the regression line obtained for the comparison of LST values with the minimum and maximum parameters.
The findings were unexpected for the following reasons:
  • For the Tunisian site, the choice of correction method had no impact on the results.
  • For the French site, the center-scene parameters provided better atmospheric correction.
Our main finding is that the temperature correction process shows significant sensitivity to the spatial variability of the parameters. The envelope of possible between the minimum and maximum values on the image as shown in Figure 14 and Figure 15 is the proof.
Table 4 presents the difference between the minimum and maximum set of parameters used. Regression line equations were used to obtain the values for each case. On the basis of these results, the error can be quantified using localized parameters for a complete Landsat image. Despite the difference in atmospheric composition between the two sites, mis-estimating any of the three atmospheric correction parameters results in the same impact on the retrieved temperature. According to the radiative transfer Equation (3) and Table 4 the point of inflection is located at around 12 ° C for France and around 15 ° C for Tunisia. This point represents the temperature at which the effect of the atmospheric correction is reverse .

5. Conclusions

We propose an automatic correction tool for atmospheric and surface emissivity data provided by Landsat remote sensors (L5, L7 and L8). The LANDARTs tool, which stands for Landsat automatic retrieval of surface temperatures, permits the conversion of Landsat mono-window radiances in the infrared band into land surface temperatures (LST) for an entire image, retaining the USGS spatial resolution (30 m for USGS resampled data). LANDARTs was developed in the programming language Python. The tool was built around the radiative transfer model MODTRAN and exploits the model’s outputs in the specific window of Landsat’s thermal spectrum. MODTRAN runs were used to determine the three parameters involved in the radiative transfer Equation (3): the atmospheric transmittance, the upwelling radiance and the downwelling radiance of the atmosphere in the spectral window and for a given vertical atmospheric profile prescribed as the input. The atmospheric profiles were obtained from the ECMWF ERA-Interim full-resolution database. For each Landsat image, a MARS request was sent to the ECMWF server to access the ERA-Interim database. The atmospheric vertical composition of the required variables (air temperature, pressure, relative humidity and ozone concentration) was downloaded on the useful’s footprint of image at the uniform lat/long grids resolution of 0.125 ° and for the heights between the Earth’s surface and 100 km in atmosphere. Calculating the three atmospheric correction parameters with this grid resolution results in a more relevant spatial variability than found in previous studies and thereby improves the spatial precision of LST estimations compared to local scale approaches. To complete the correction process, surface emissivity was obtained through the NDVI empirical relationship (Equation (8)) in accordance with Wittich [37] and adapted by Olioso [39]. The present paper implemented standard values; however, LANDARTs allows users to customize parameters, such as surface soil emissivity and to set limits for maximum and minimum NDVI values according to the user’s knowledge of the region of interest.
A comparison of estimated LST with local surface temperature measurements from two ground networks (the first located in southwestern France near Toulouse, the second located in central Tunisia near Kairouan city) showed good agreement. Validations for the two sites resulted in an R M S E of 2.55 ° C for France and an R M S E of 1.8 ° C for Tunisia with a correlation coefficient of 0.98. These results confirm that, on a local scale, good quality, temporally-stable atmospheric correction was obtained with ERA-Interim atmospheric profiles and MODTRAN model. The LST Landsat data products used with the LANDARTs tool were also compared to two synchronous ASTER tiles for which kinetic temperature products were available. ASTER dates correspond to winter and the L8 sensor for Tunisia and to spring and the L7 (SLC-on) sensor for the French site. As a result, comparisons present spatially-consistent differences that were mostly attributable to emissivity variability, but without apparent atmospheric variability. Pixel-by-pixel comparisons presented a typical RMSE of 2.5 ° C for the Tunisian tile and 3.5 ° C for the French tiles. Both cases showed a systematic bias of 2.3 ° C in favor of ASTER products, which was attributed to the later overpass of ASTER. Finally, we were able to demonstrate that absolute surface temperature estimates, and the absolute range of the atmospheric correction parameters may present significant spatial variability, even at an ASTER-tile scale. At the scale of a Landsat scene, the effect may attain an R M S E of ±2 ° C. Moreover, the variability of the atmospheric correction parameters observed for each scene changes very little over time, at least for the sites considered here. To conclude, we propose to the community a correction tool developed in the Python programming language and based on Barsi’s Atmocorr methodology. This tool implements the commercially-available radiative transfer reference model MODTRAN and the freely-available European Weather reanalysis ERA-Interim. The primary purpose of generating LST outputs with LANDARTs is to provide values that are representative of surface measurements and that offer consistent quality on the scale of an entire Landsat tile.

Acknowledgments

The present study received support from the CNES (Centre National d’Etude Spatiale, France) through funding of student internships and a 12-month engineer vacation. Additional funding was provided by the CNRS (Centre National de Recherche Scientifique, France) in the form of a permanent position for a research engineer.

Author Contributions

The study was initiated and supervised by Vincent Rivalland, Benjamin Tardy completed the technical work supported by Mireille Huc. Sebastien Marcq (CNES) supported the project and defined some specifications. Gilles Boulet and Olivier Hagolle provided general help in the redaction of this technical note.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Anderson, M.C.; Hain, C.; Wardlow, B.; Pimstein, A.; Mecikalski, J.R.; Kustas, W.P. Evaluation of drought indices based on thermal remote sensing of evapotranspiration over the continental United States. J. Clim. 2011, 24, 2025–2044. [Google Scholar] [CrossRef]
  2. Kustas, W.; Anderson, M. Advances in thermal infrared remote sensing for land surface modeling. Agric. For. Meteorol. 2009, 149, 2071–2081. [Google Scholar] [CrossRef]
  3. Anderson, M.C.; Allen, R.G.; Morse, A.; Kustas, W.P. Use of Landsat thermal imagery in monitoring evapotranspiration and managing water resources. Remote Sens. Environ. 2012, 122, 50–65. [Google Scholar] [CrossRef]
  4. Le Traon, P.Y.; Antoine, D.; Bentamy, A.; Bonekamp, H.; Breivik, L.; Chapron, B.; Corlett, G.; Dibarboure, G.; DiGiacomo, P.; Donlon, C.; et al. Use of satellite observations for operational oceanography: Recent achievements and future prospects. J. Oper. Oceanogr. 2015, 8, s12–s27. [Google Scholar] [CrossRef]
  5. Tang, H.; Li, Z.L. Quantitative Remote Sensing in Thermal Infrared: Theory and Applications, 1st ed.; Springer-Verlag: Berlin/Heidelberg, Germany, 2014. [Google Scholar]
  6. Wan, Z.; Zhang, Y.; Zhang, Q.; Li, Z.L. Quality assessment and validation of the MODIS global land surface temperature. Int. J. Remote Sens. 2004, 25, 261–274. [Google Scholar] [CrossRef]
  7. Wan, Z.; Dozier, J. A generalized split-window algorithm for retrieving land-surface temperature from space. IEEE Trans. Geosci. Remote Sens. 1996, 34, 892–905. [Google Scholar]
  8. Masiello, G.; Serio, C.; De Feis, I.; Amoroso, M.; Venafra, S.; Trigo, I.; Watts, P. Kalman filter physical retrieval of surface emissivity and temperature from geostationary infrared radiances. Atmos. Meas. Tech. 2013, 6, 3613–3634. [Google Scholar] [CrossRef]
  9. Masiello, G.; Serio, C.; Venafra, S.; Liuzzi, G.; Göttsche, F.; Trigo, I.; Watts, P. Kalman filter physical retrieval of surface emissivity and temperature from SEVIRI infrared channels: A validation and intercomparison study. Atmos. Meas. Tech. 2015, 8, 2981–2997. [Google Scholar] [CrossRef]
  10. Abrams, M. The Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER): Data products for the high spatial resolution imager on NASA’s Terra platform. Int. J. Remote Sens. 2000, 21, 847–859. [Google Scholar] [CrossRef]
  11. USGS. Pages dedicated to Landsat missions. Calibration Notices of January 6, 2014— Landsat 8 Reprocessing Details. Available online: http://landsat.usgs.gov/calibration_notices.php (accessed on 20 August 2016).
  12. USGS. Pages dedicated to Landsat missions. Calibration Notices of January 29, 2014—Landsat 8 Reprocessing to Begin February 3, 2014. Available online: http://landsat.usgs.gov/calibration_notices.php (accessed on 20 August 2016).
  13. Jiménez-Muñoz, J.C.; Sobrino, J.A. A generalized single-channel method for retrieving land surface temperature from remote sensing data. J. Geophys. Res. 2003, 108. [Google Scholar] [CrossRef]
  14. Jiménez-Muñoz, J.C.; Cristóbal, J.; Sobrino, J.; Soria, G.; Ninyerola, M.; Pons, X.; et al. Revision of the single-channel algorithm for land surface temperature retrieval from Landsat thermal-infrared data. IEEE Trans. Geosci. Remote Sens. 2009, 47, 339–349. [Google Scholar] [CrossRef]
  15. Coll, C.; Galve, J.; Sanchez, J.; Caselles, V. Validation of Landsat-7/ETM+ thermal-band calibration and atmospheric correction with ground-based measurements. IEEE Trans. Geosci. Remote Sens. 2010, 48, 547–555. [Google Scholar] [CrossRef]
  16. Sobrino, J.A.; Jiménez-Muñoz, J.C.; Paolini, L. Land surface temperature retrieval from LANDSAT TM 5. Remote Sens. Environ. 2004, 90, 434–440. [Google Scholar] [CrossRef]
  17. Qin, Z.H.; Karnieli, A.; Berliner, P. A mono-window algorithm for retrieving land surface temperature from Landsat TM data and its application to the Israel-Egypt border region. Int. J. Remote Sens. 2001, 22, 3719–3746. [Google Scholar] [CrossRef]
  18. Barsi, J.A.; Barker, J.L.; Schott, J.R. An atmospheric correction parameter calculator for a single thermal band earth-sensing instrument. Int. Geosci. Remote Sens. Symp. 2003, 5, 3014–3016. [Google Scholar]
  19. Berk, A.; Anderson, G.; Acharya, P.; Chetwynd, J.; Bernstein, L.; Shettle, E.; Matthew, M.; Adler-Golden, S. MODTRAN4 User Manual; Air Force Research Laboratory, Space Vehicles Directorate: Wright-Patterson AFB, OH, USA, 1999. [Google Scholar]
  20. McCarville, D.; Buenemann, M.; Bleiweiss, M.; Barsi, J. Atmospheric correction of Landsat thermal infrared data: A calculator based on North American Regional Reanalysis (NARR) data. In Proceedings of the American Society for Photogrammetry and Remote Sensing Conference, Milwaukee, WI, USA, 1–5 May 2011.
  21. Barsi, J.A.; Schott, J.R.; Palluconi, F.D.; Hook, S.J. Validation of a Web-Based Atmospheric Correction Tool for Single Thermal Band Instruments. Proc. SPIE 2005. [Google Scholar] [CrossRef]
  22. Zhou, J.; Li, J.; Zhang, L.; Hu, D.; Zhan, W. Intercomparison of methods for estimating land surface temperature from a Landsat-5 TM image in an arid region with low water vapour in the atmosphere. Int. J. Remote Sens. 2012, 33, 2582–2602. [Google Scholar] [CrossRef]
  23. Mira, M.; Olioso, A.; Rivalland, V.; Courault, D.; Marloie, O.; Guillevic, P. Quantifying uncertainties in land surface temperature due to atmospheric correction: Application to Landsat-7 data over a Mediterranean agricultural region. In Proceedings of the 2014 IEEE Geoscience and Remote Sensing Symposium, Quebec City, QC, Canada, 13–18 July 2014; pp. 2375–2378.
  24. Chander, G.; Markham, B.L.; Helder, D.L. Summary of current radiometric calibration coefficients for Landsat MSS, TM, ETM+, and EO-1 ALI sensors. Remote Sens. Environ. 2009, 113, 893–903. [Google Scholar] [CrossRef]
  25. Barsi, J.A.; Schott, J.R.; Hook, S.J.; Raqueno, N.G.; Markham, B.L.; Radocinski, R.G. Landsat-8 thermal infrared sensor (TIRS) vicarious radiometric calibration. Remote Sens. 2014, 6, 11607–11626. [Google Scholar] [CrossRef]
  26. Berk, A.; Anderson, G.P.; Acharya, P.K.; Bernstein, L.S.; Muratov, L.; Lee, J.; Fox, M.; Adler-Golden, S.M.; Chetwynd, J.H.; Hoke, M.L.; et al. MODTRAN 5: A reformulated atmospheric band model with auxiliary species and practical multiple scattering options: Update. Proc. SPIE 2005. [Google Scholar] [CrossRef]
  27. Dee, D.P.; Uppala, S.M.; Simmons, A.J.; Berrisford, P.; Poli, P.; Kobayashi, S.; Andrae, U.; Balmaseda, M.A.; Balsamo, G.; Bauer, P.; et al. The ERA-Interim reanalysis: Configuration and performance of the data assimilation system. Q. J. R. Meteorol. Soc. 2011, 137, 553–597. [Google Scholar] [CrossRef]
  28. API ERA-Interim. Available online: https://software.ecmwf.int/wiki/display/WEBAPI/ (accessed on 20 August 2016).
  29. COESA: US Commission/Stand Atmosphere (Compiler); National Oceanic and Atmospheric Administration (Collaborator); National Aeronautics and Space Administration (Collaborator); United States Air Force (Collaborator). U.S. Standard Atmosphere, 1976 (NOAA Document S/T 76-1562), 1st ed.NOAA; NASA; USAF: Silver Spring, MD, USA, 2004.
  30. Lawrence, M.G. The relationship between relative humidity and the dewpoint temperature in moist air—A simple conversion and applications. Bull. Am. Meteorol. Soc. 2005, 86. [Google Scholar] [CrossRef]
  31. Alduchov, O.; Eskridge, R. Improved magnus form approximation of saturation vapor pressure. J. Appl. Meteorol. 1996, 35, 601–609. [Google Scholar] [CrossRef]
  32. Sobrino, J.A.; Jiménez-Muñoz, J.C.; Sòria, G.; Romaguera, M.; Guanter, L.; Moreno, J.; Plaza, A.; Martínez, P. Land surface emissivity retrieval from different VNIR and TIR sensors. IEEE Trans. Geosci. Remote Sens. 2008, 46, 316–327. [Google Scholar] [CrossRef]
  33. Li, Z.L.; Wu, H.; Wang, N.; Qiu, S.; Sobrino, J.A.; Wan, Z.; Tang, B.H.; Yan, G. Land surface emissivity retrieval from satellite data. Int. J. Remote Sens. 2013, 34, 3084–3127. [Google Scholar] [CrossRef]
  34. Olioso, A. Estimating the difference between brigthness and surface temperatures for a vegetal canopy. Agric. For. Meteorol. 1995, 72, 237–242. [Google Scholar] [CrossRef]
  35. Ren, H.; Yan, G.; Chen, L.; Li, Z. Angular effect of MODIS emissivity products and its application to the split-window algorithm. ISPRS J. Photogramm. Remote Sens. 2011, 66, 498–507. [Google Scholar] [CrossRef]
  36. Li, Z.L.; Becker, F.; Stoll, M.; Wan, Z. Evaluation of six methods for extracting relative emissivity spectra from thermal infrared images. Remote Sens. Environ. 1999, 69, 197–214. [Google Scholar] [CrossRef]
  37. Wittich, K.P. Some simple relationships between land-surface emissivity, greenness and the plant cover fraction for use in satellite remote sensing. Int. J. Biometeorol. 1997, 41, 58–64. [Google Scholar] [CrossRef]
  38. Baldridge, A.; Hook, S.; Grove, C.; Rivera, G. The ASTER spectral library version 2.0. Remote Sens. Environ. 2009, 113, 711–715. [Google Scholar] [CrossRef]
  39. Olioso, A.; Mira, M.; Courault, D.; Marloie, O.; Guillevic, P. Impact of surface emissivity and atmospheric conditions on surface temperatures estimated from top of canopy brightness temperatures derived from Landsat 7 data. In Proceedings of the 2013 IEEE International Geoscience and Remote Sensing Symposium, Melbourne, Australia, 21–26 July 2013; pp. 3033–3036.
  40. Zhang, R.H.; Li, Z.L.; Tang, X.Z.; Sun, X.M.; Su, H.B.; Zhu, C.; Zhu, Z.L. Study of emissivity scaling and relativity of homogeneity of surface temperature. Int. J. Remote Sens. 2004, 25, 245–259. [Google Scholar] [CrossRef]
  41. Beziat, P.; Ceschia, E.; Dedieu, G. Carbon balance of a three crop succession over two cropland sites in South West France. Agric. For. Meteorol. 2009, 149, 1628–1645. [Google Scholar] [CrossRef] [Green Version]
  42. Tallec, T.; Beziat, P.; Jarosz, N.; Rivalland, V.; Ceschia, E. Crops’ water use efficiencies in temperate climate: Comparison of stand, ecosystem and agronomical approaches. Agric. For. Meteorol. 2013, 168, 69–81. [Google Scholar] [CrossRef]
  43. Idso, S.B. A set of equations for full spectrum and 8- to 14-μm and 10.5- to 12.5-μm thermal radiation from cloudless skies. Water Resour. Res. 1981, 17, 295–304. [Google Scholar] [CrossRef]
  44. Gillespie, A.R.; Rokugawa, S.; Hook, S.J.; Matsunaga, T.; Kahle, A.B. Temperature/Emissivity Separation Algorithm Theoretical Basis Document; Version 2.4; NASA/GSFC: Greenbelt, MD, USA, 1999. [Google Scholar]
  45. Hagolle, O.; Huc, M.; Pascual, D.V.; Dedieu, G. A multi-temporal method for cloud detection, applied to FORMOSAT-2, VENμS, LANDSAT and SENTINEL-2 images. Remote Sens. Environ. 2010, 114, 1747–1755. [Google Scholar] [CrossRef] [Green Version]
  46. Zhou, J.; Li, M.; Liu, S.; Jia, Z.; Ma, Y. Validation and performance evaluations of methods for estimating land surface temperatures from ASTER data in the middle reach of the Heihe River Basin, Northwest China. Remote Sens. 2015, 7, 7126. [Google Scholar] [CrossRef]
  47. Srivastava, P.; Majumdar, T.; Bhattacharya, A.K. Study of land surface temperature and spectral emissivity using multi-sensor satellite data. J. Earth Syst. Sci. 2010, 119, 67–74. [Google Scholar] [CrossRef]
  48. Gillespie, A.; Rokugawa, S.; Matsunaga, T.; Cothern, J.S.; Hook, S.; Kahle, A.B. A temperature and emissivity separation algorithm for Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) images. IEEE Trans. Geosci. Remote Sens. 1998, 36, 1113–1126. [Google Scholar] [CrossRef]
Figure 1. Localization map of the both network (green triangles), for France (left) and for Tunisia (right). Background are true-color images of L7 (left) and L8 (right) acquired respectively on 29 May 2003 and 29 March 2013. Images are scaled and oriented. The black boxes on images represent the ASTER image footprint used on each location to validate tool (see Section 4.2).
Figure 1. Localization map of the both network (green triangles), for France (left) and for Tunisia (right). Background are true-color images of L7 (left) and L8 (right) acquired respectively on 29 May 2003 and 29 March 2013. Images are scaled and oriented. The black boxes on images represent the ASTER image footprint used on each location to validate tool (see Section 4.2).
Remotesensing 08 00696 g001
Figure 2. Flowchart of the LANDARTs algorithm to generate LST images at a resolution of 30 m.
Figure 2. Flowchart of the LANDARTs algorithm to generate LST images at a resolution of 30 m.
Remotesensing 08 00696 g002
Figure 3. Comparison between Landsat LST and in situ temperatures (in degree Celsius) for the Tunisian site.
Figure 3. Comparison between Landsat LST and in situ temperatures (in degree Celsius) for the Tunisian site.
Remotesensing 08 00696 g003
Figure 4. Comparison between Landsat LST and in situ temperatures (in degree Celsius) for the French site.
Figure 4. Comparison between Landsat LST and in situ temperatures (in degree Celsius) for the French site.
Remotesensing 08 00696 g004
Figure 5. Pixel-by-pixel comparison between the ASTER surface kinetic product and the Landsat LST map at a resolution of 270 m for the Tunisian site.
Figure 5. Pixel-by-pixel comparison between the ASTER surface kinetic product and the Landsat LST map at a resolution of 270 m for the Tunisian site.
Remotesensing 08 00696 g005
Figure 6. Pixel-by-pixel comparison between the ASTER surface kinetic product and the Landsat LST map at a resolution of 270 m for the French site.
Figure 6. Pixel-by-pixel comparison between the ASTER surface kinetic product and the Landsat LST map at a resolution of 270 m for the French site.
Remotesensing 08 00696 g006
Figure 7. Difference in degree Celsius between ASTER surface kinetic product and Landsat LST map at 270 m resolution on Tunisia (29 March 2013).
Figure 7. Difference in degree Celsius between ASTER surface kinetic product and Landsat LST map at 270 m resolution on Tunisia (29 March 2013).
Remotesensing 08 00696 g007
Figure 8. Disparity histogram in degree Celsius between ASTER surface kinetic product and Landsat LST map at 270 m resolution on Tunisia (29 March 2013).
Figure 8. Disparity histogram in degree Celsius between ASTER surface kinetic product and Landsat LST map at 270 m resolution on Tunisia (29 March 2013).
Remotesensing 08 00696 g008
Figure 9. Difference comparison in degree Celsius between ASTER surface kinetic product and Landsat LST map at 270 m resolution on France (29 May 2003).
Figure 9. Difference comparison in degree Celsius between ASTER surface kinetic product and Landsat LST map at 270 m resolution on France (29 May 2003).
Remotesensing 08 00696 g009
Figure 10. Disparity histogram in degree Celsius between ASTER surface kinetic product and Landsat LST map at 270 m resolution on France (29 May 2003).
Figure 10. Disparity histogram in degree Celsius between ASTER surface kinetic product and Landsat LST map at 270 m resolution on France (29 May 2003).
Remotesensing 08 00696 g010
Figure 11. Parameter variability and influence on LST correction for an L7 image taken on 29 May 2003 (199,030). The square black frame represents the ASTER valid data footprint used for the spatial comparison. (a) Transmittance (unitless); (b) upwelling radiance ( W · m 2 · sr 1 · μ m 1 ); (c) map of emissivity estimates derived from NDVI values used in both corrections (unitless); (d) disparity map, in degree Celsius, comparing LST estimated using spatial parameters with LST estimated using the center-scene correction parameters. White pixels are cloudy or snow-covered pixels detected with the level 2A processor MACCS [45].
Figure 11. Parameter variability and influence on LST correction for an L7 image taken on 29 May 2003 (199,030). The square black frame represents the ASTER valid data footprint used for the spatial comparison. (a) Transmittance (unitless); (b) upwelling radiance ( W · m 2 · sr 1 · μ m 1 ); (c) map of emissivity estimates derived from NDVI values used in both corrections (unitless); (d) disparity map, in degree Celsius, comparing LST estimated using spatial parameters with LST estimated using the center-scene correction parameters. White pixels are cloudy or snow-covered pixels detected with the level 2A processor MACCS [45].
Remotesensing 08 00696 g011
Figure 12. Multi-annual temporal evolution of the atmospheric transmittance parameter (unitless) observed for the set of French images (path/row 199/030). The blue line represents the center-scene transmittance, and the red line is the transmittance for the Aurade site extracted from the spatial parameter image at a resolution of 30 m. The green fill represents the amplitude (min/max) of the transmittance for each Landsat scene.
Figure 12. Multi-annual temporal evolution of the atmospheric transmittance parameter (unitless) observed for the set of French images (path/row 199/030). The blue line represents the center-scene transmittance, and the red line is the transmittance for the Aurade site extracted from the spatial parameter image at a resolution of 30 m. The green fill represents the amplitude (min/max) of the transmittance for each Landsat scene.
Remotesensing 08 00696 g012
Figure 13. Multi-annual temporal evolution of the atmospheric upwelling radiance parameter ( W · m 2 · sr 1 · μ m 1 ) observed for the set of French images (path/row 199/030). The blue line represents the center-scene upwelling radiance, and the red line is the upwelling radiance for the Aurade site extracted from the spatial parameter image at a resolution of 30 m. The green fill represents the amplitude (min/max) of the upwelling radiance for each Landsat scene.
Figure 13. Multi-annual temporal evolution of the atmospheric upwelling radiance parameter ( W · m 2 · sr 1 · μ m 1 ) observed for the set of French images (path/row 199/030). The blue line represents the center-scene upwelling radiance, and the red line is the upwelling radiance for the Aurade site extracted from the spatial parameter image at a resolution of 30 m. The green fill represents the amplitude (min/max) of the upwelling radiance for each Landsat scene.
Remotesensing 08 00696 g013
Figure 14. Comparison in degree Celsius between Landsat LST and Tunisia in situ temperature for the differing correction methods.
Figure 14. Comparison in degree Celsius between Landsat LST and Tunisia in situ temperature for the differing correction methods.
Remotesensing 08 00696 g014
Figure 15. Comparison in degree Celsius between Landsat LST and French in situ temperature for the differing correction methods.
Figure 15. Comparison in degree Celsius between Landsat LST and French in situ temperature for the differing correction methods.
Remotesensing 08 00696 g015
Table 1. Landsat thermal bands specifications.
Table 1. Landsat thermal bands specifications.
SensorWavelength ( μ m)Native ResolutionBand
L510.40–12.50120 mB6
L710.40–12.5060 mB6
L810.60–11.19100 mB10
L811.50–12.51100 mB11
Table 2. ASTER thermal band specifications.
Table 2. ASTER thermal band specifications.
Band No.Wavelength ( μ m)Acquisition Resolution
108.125–8.47590 m
118.475–8.82590 m
128.925–9.27590 m
1310.25–10.9590 m
1410.95–11.6590 m
Table 3. Dates available for spatial validation.
Table 3. Dates available for spatial validation.
Site (Lat., Lon.)SensorDateAcquisition Time (hh:mm:ss)
TunisiaL829 March 201309:59:14
( 9 ° 46 13.14 E, 36 ° 5 45.65 N)ASTER29 March 201310:12:52
FranceL729 May 200310:30:47
( 0 ° 11 52.99 E, 43 ° 58 13.04 N)ASTER29 May 200311:00:05
Table 4. Differences between LST retrieved using the different sets of local parameters corresponding to set of maximum transmittance minus and of minimum transmittance. The table values are estimated using equations of regression line that define the green envelope drawn in Figure 14 and Figure 15.
Table 4. Differences between LST retrieved using the different sets of local parameters corresponding to set of maximum transmittance minus and of minimum transmittance. The table values are estimated using equations of regression line that define the green envelope drawn in Figure 14 and Figure 15.
Sensor Temperature ( ° C) δ T Tunisia ( ° C) δ T France ( ° C)
0−1.64−0.97
10−0.520.05
200.591.09
301.712.12
402.823.15
503.944.19
605.055.22

Share and Cite

MDPI and ACS Style

Tardy, B.; Rivalland, V.; Huc, M.; Hagolle, O.; Marcq, S.; Boulet, G. A Software Tool for Atmospheric Correction and Surface Temperature Estimation of Landsat Infrared Thermal Data. Remote Sens. 2016, 8, 696. https://doi.org/10.3390/rs8090696

AMA Style

Tardy B, Rivalland V, Huc M, Hagolle O, Marcq S, Boulet G. A Software Tool for Atmospheric Correction and Surface Temperature Estimation of Landsat Infrared Thermal Data. Remote Sensing. 2016; 8(9):696. https://doi.org/10.3390/rs8090696

Chicago/Turabian Style

Tardy, Benjamin, Vincent Rivalland, Mireille Huc, Olivier Hagolle, Sebastien Marcq, and Gilles Boulet. 2016. "A Software Tool for Atmospheric Correction and Surface Temperature Estimation of Landsat Infrared Thermal Data" Remote Sensing 8, no. 9: 696. https://doi.org/10.3390/rs8090696

APA Style

Tardy, B., Rivalland, V., Huc, M., Hagolle, O., Marcq, S., & Boulet, G. (2016). A Software Tool for Atmospheric Correction and Surface Temperature Estimation of Landsat Infrared Thermal Data. Remote Sensing, 8(9), 696. https://doi.org/10.3390/rs8090696

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop