Derivation of Land Surface Temperature for Landsat-8 TIRS Using a Split Window Algorithm

Land surface temperature (LST) is one of the most important variables measured by satellite remote sensing. Public domain data are available from the newly operational Landsat-8 Thermal Infrared Sensor (TIRS). This paper presents an adjustment of the split window algorithm (SWA) for TIRS that uses atmospheric transmittance and land surface emissivity (LSE) as inputs. Various alternatives for estimating these SWA inputs are reviewed, and a sensitivity analysis of the SWA to misestimating the input parameters is performed. The accuracy of the current development was assessed using simulated Modtran data. The root mean square error (RMSE) of the simulated LST was calculated as 0.93 °C. This SWA development is leading to progress in the determination of LST by Landsat-8 TIRS.


Introduction
Land surface temperature (LST) is related to surface energy and water balance, at local through global scales, with principal significance for a wide variety of applications, such as climate change, urban climate, the hydrological cycle, and vegetation monitoring [1][2][3][4]. LST variations in space and time, measured by satellite remote sensing, are used for the estimation of a multitude of geophysical variables, such as evapotranspiration, vegetation water stress, soil moisture, and thermal inertia [5][6][7]. With the increasing recognition of the importance of LST, methods for its estimation from space have continuously been developed [8]. In recent decades, sensors, such as the Moderate-resolution Imaging Spectroradiometer (MODIS) and the Advanced Very High Resolution Radiometer (AVHRR), have provided public domain global thermal data twice daily, using two long-wave infrared (LWIR) bands. Landsat-5 Thematic Mapper (TM) and Landsat-7 Enhanced Thematic Mapper Plus (ETM+) provide thermal data using just one long-wave infrared (LWIR) band, with a higher spatial resolution but with a 16-day temporal resolution. Since satellite remote sensing provides a repetitive synoptic view in short intervals of the Earth's surface, it is a vital tool for monitoring LST.
Landsat-8 was successfully launched on 11 February 2013 and deployed into orbit with two instruments on-board: (1) the Operational Land Imager (OLI) with nine spectral bands in the visual (VIS), near infrared (NIR), and the shortwave infrared (SWIR) spectral regions; and (2) the Thermal Infrared Sensor (TIRS) with two spectral bands in the LWIR. The relative spectral response of the TIRS bands is presented in Figure 1. The two TIRS bands were selected to enable the atmospheric correction of the thermal data using a split window algorithm (SWA) [9,10]. The use of two separate, relatively narrow, thermal bands has been shown to minimize the error in the retrieval of LST [11]. The spatial resolution of TIRS data is 100 m with a revisit time of 16 days, and as a result, applications are different than those of other sensors with coarser spatial resolutions and shorter revisit times. While Landsat-8 images are already freely distributed through the U.S. Geological Survey (USGS), to the best of our knowledge, no SWA for LST retrieval from TIRS has been published. Although several SWAs have been developed for use with other sensors [12][13][14][15], some adaptations are required in order to implement them for the TIRS spectral bands. Therefore, the objective of this letter is to develop a SWA, adapted for use with Landsat-8 TIRS data, along with its accuracy assessment.

Split Window Algorithm (SWA)
The SWA was first proposed by McMillin [17] who suggested using the differences in the atmospheric absorbance of two adjacent LWIR bands in order to accurately retrieve the sea surface temperature (SST). In order to make the transition from SST to LST retrieval, one has to assume the land surface emissivity (LSE) in both bands a priori [12]. Qin et al. [18] have presented a SWA for AVHRR that requires only two essential variables: LSE and atmospheric transmittance. They tested their algorithm and found its accuracy to be 1.75 °C for real world data. Additionally, they found that it was preferable to the other SWAs that also performed well but required some parameters that are difficult to estimate [19,20]. Therefore, the SWA suggested by Qin et al. [18] was chosen to be adapted for TIRS, not only because it was tested and proved to be accurate, but also because it is reasonable to estimate its input parameters, as will be discussed in the next section of the paper. Further, the current work complements that of Qin et al. [18] to which we refer the reader for theoretical background and algorithm development.
The adaptation of the SWA for TIRS bands relies on the determination of parameter L i for the TIRS-specific spectral bands. L i is defined as: (1) in which L i has the dimension of temperature in Kelvin, B i (T) is the Planck function radiance, spectrally integrated over each of the TIRS bands at temperature T, and ∂B i (T) is the derivative of the Planck function for band i at temperature T. Therefore, ∂B i (T)/∂T can be calculated as: (2) We computed L i numerically by using ∂B i (T)/∂T and accurately fit it using a linear regression  Table 1. From this table, we can see that as the range of T decreases, a better accuracy, or a lower SEE, is achieved. In this paper, the SWA accuracy assessment was carried out for the strict case of the extreme temperature range of 0-60 °C. However, in order to obtain the LST as accurately as possible, it is advisable to select coefficients according to the temperature range in the image.  The SWA is derived from a first-order Taylor-series linearization of the radiative transfer equation [21], and its formulation by Qin et al. [18] takes the general form of: where T s is the LST, T 10 and T 11 are the brightness temperatures of TIRS bands 10 and 11, respectively, and A 0 , A 1 , and A 2 are coefficients determined by the atmospheric transmittance and LSE in both TIRS bands: Based on the algorithm developed by Qin et al. [18], we define: where ε i is the LSE of band i and τ i (θ) is the atmospheric transmittance for a given zenith view angle θ in band i . Accordingly, the parameters in Equation (4) are defined: The algorithm in Equation (4) uses the estimation of two geophysical parameters, namely atmospheric transmittance and LSE, in order to estimate the LST. The next section will briefly discuss the estimation of these parameters.

Determination of Atmospheric Transmittance
Some SWAs were developed for wide swath sensors (e.g., NOAA-AVHRR, MODIS), and consequently emphasize the correction for zenith view angle effects on the atmospheric transmittance. However, in the case of the Landsat-8 TIRS at an altitude of 705 km with a swath of 185 km, the maximum zenith view angle is about 7.5°. At that angle, the effect on the atmospheric transmittance in both LWIR bands is negligible [18]. Thus, the term θ may be removed from Equations (5) and (6) for the purpose of implementing this SWA for TIRS.
In addition, transmittance is wavelength dependent and, therefore, different for each of the TIRS bands. Absorption in the 10.5-12.5 µm atmospheric window is mainly affected by water vapor that has a high spatial variability, since other atmospheric gases, such as 2 CO , 2 N , and 3 O , are well-mixed, and their effect can be considered constant throughout an image for the purpose of this analysis [22,23]. The atmospheric water vapor content, at the time of image acquisition, can be obtained from local measurements in-situ, or at nearby meteorological stations. The sun photometer measurements of the AErosol RObotic NETwork (AERONET), operated by NASA/GSFC [24], are suggested as a good source of water vapor data by Qin et al. [18], but they only provide data when a direct line of sight can be established between the station and the sun. Thus, sun photometer readings are a good source of information only for the processing of daytime and clear sky images. Unlike LSE (discussed in the next section), it is not practical to estimate the atmospheric transmittance per pixel. Consequently, when using sun photometer point measurements, we make a latent assumption that the atmosphere is constant throughout the scene. While this is not always the case, this first-order estimation of the atmospheric transmittance for the whole scene, which relies on measured water vapor, is preferable to not accounting for it at all [25]. To roughly meet this assumption, cloud pixels must be masked. The atmospheric transmittance can then be simulated for the entire scene, based on the point measurement of the total content of water vapor in the column and standard atmospheric profiles using radiative transfer models, such as Modtran.
The use of sun photometers for estimating water vapor has some limitations. As mentioned earlier, they do not provide information when the sky is cloudy or at night. Furthermore, due to the lack of global coverage by ground-based instruments, AERONET data might not be available for all TIRS users for estimating the atmospheric water vapor content, especially in the northern latitudes, Asia, Australia, and central Africa. In these areas, other means are suggested, such as the use of a total column water vapor product provided by European Centre for Medium-Range Weather Forecasts (ECMWF) or the Moderate Resolution Imaging Spectroradiometer (MODIS) Precipitable Water product (MOD05_L2, MYD05_L2). While these products' spatial resolutions are significantly coarser than TIRS, they capture some spatial variability in water vapor distribution. Accordingly, ECMWF and MODIS data can be used when the users wish to refrain from assuming constant water vapor throughout the scene. The disadvantage of these data products is that in many cases, they do not accurately represent the state of the atmosphere at the time of TIRS overpass. Unlike the sun photometer that can perform measurements instantaneously or at a close temporal proximity to TIRS, the water vapor estimations derived from these products may be produced within several hours of the TIRS acquisition, and thus do not capture the exact atmospheric conditions at the acquisition time. This is indeed a limitation, and for this reason, these products should be used with caution. As a rule, users should select the data source according to its accessibility. When multiple sources of data are available, temporal proximity to the TIRS image acquisition and spatial proximity of the water vapor measurement site to the study area should be preferred. When facing a dilemma between temporal and spatial precision, the users have to consider this trade-off and make an optimal decision based on their experience.
The results of Modtran 4.0 simulations, conducted for a mid-latitude summer and for a 1976 standard US atmospheric profile to determine the relation between water vapor and atmospheric transmittance, are presented in Table 2. Throughout this paper, the model developed for mid-latitude summer is used as an example. Users may opt to use the coefficients for the 1976 standard US atmospheric profile where the atmosphere is modeled better. Based on AERONET measurements, the water vapor content in the current research area, in the northern Negev Desert, Israel, ranges from 0.5 to 3 g·cm −2 . As can be seen in Table 2, the relation between water vapor and transmittance is close to linear. Qin et al. [18] showed that when this relation is evaluated for a larger range of water vapor values, it is better to divide the range into several sections and evaluate each of them separately in order to achieve a better accuracy. Since the atmosphere in our research area is relatively dry, and consequently, the range of water vapor values is relatively small, we are satisfied with the accuracy achieved by treating the entire range as a whole. Taking into account a plausible error in the water vapor content estimation of 0.2 g·cm −2 , as suggested by Qin et al. (2001) [18] and according to the regression coefficients in Table 2, we can determine that the error in the atmospheric transmittance estimation is less than 0.031, which is slightly lower than the value obtained by Qin et al. [18].

Determination of LSE
The emissivity of land, in contrast to that of the ocean, is significantly different than unity, and varies with the heterogeneity of vegetation, surface moisture, roughness, and viewing angle [26]. Since LSE can change substantially over short distances, it is important to estimate its value for every pixel prior to applying the SWA. Several methods have been suggested to estimate the emissivity for other sensors and can also be applied to TIRS. Techniques for emissivity estimations from infrared and visible data are reviewed and discussed in detail elsewhere [27,28]. Adapting some of these techniques for use with the Landsat-8 requires the use of OLI bands to indirectly estimate the LSE in the TIRS bands. For instance, the LSE could be obtained from a land-cover classification, in which the emissivity values for each class are assumed [29]. This type of approach is exercised for MODIS LST and emissivity products. However, the estimated emissivity in arid and semi-arid areas is potentially uncertain, and users are advised to exercise caution in their applications. Of course, since vegetative cover tends to change with time [30], good knowledge of the study site and in situ LSE measurements of representative ground covers of the different classes that coincide with the satellite overpass are desired. However, if the required conditions for the classification approach are not met, it is possible to use NDVI by retrieving the proportions of soil and vegetation in order to estimate the LSE [31][32][33].

Sensitivity Analysis
Several scenarios were considered in order to estimate the possible LST estimation error due to misestimating the SWA input parameters: atmospheric transmittance and LSE. These scenarios included an LST range of 0-60 °C and a T 10 -T 11 range from −3 to 3 °C.

Sensitivity Analysis to Water Vapor Content
Since transmittance is derived from atmospheric water vapor content, it is expected that transmittance estimation errors will occur simultaneously in TIRS bands 10 and 11. Therefore, the sensitivity analysis was conducted for water vapor content, which serves as the input to the model, and by which the simultaneously occurring transmittance errors can be estimated according to the regression coefficients in Table 2.
Estimation error of LST is almost independent of temperature change. It changes less than ±0.02 °C over the temperature range of 0-60 °C, assuming T 10 -T 11 = −2.3 (the average case for the Sinai-Negev dune field, as seen in Figure 3A), e 10 = 0.967, e 11 = 0.971, and underestimating the water vapor content in the atmospheric column by 0.2 g·cm −2 for the water vapor content range of 0.5-3 g·cm −2 . This minute change is practically negligible.  Figure 3B.
The LST estimation error increases when the atmospheric water vapor content decreases (and, thus, the atmospheric transmittance increases). This effect increases rapidly as the brightness temperature difference between TIRS bands 10 and 11 increases, as seen in Figure 3B,C. However, the contribution of the water vapor estimation error to the LST estimation is complex and also depends on the emissivity in both channels, as seen in Figure 4.

Sensitivity Analysis to LSE
An error in LSE estimation can occur simultaneously for both of the TIRS bands, but a separate error for each of the bands is possible. In the following analysis, we present the example of a simultaneous error in both bands. Unlike the relative indifference of the LST estimation error to temperature when the water vapor content is misestimated, the LST error is sensitive to temperature when the LSE is misestimated. Figure 5A,B presents this sensitivity. The LST error increases linearly, with the LST, and decreases linearly when the LSE increases. Figure 6A,B presents a similar linear dependence of the LST error by T 10 -T 11 ; as T 10 -T 11 increases, the LST error decreases. In both Figures 5 and 6, the change of the LST error increases at the same rate as the LSE estimation error. In addition, the LST estimation error is sensitive to the water vapor content when the LSE is misestimated, as depicted in Figure 7. When water vapor content is higher, and the atmospheric transmittance is lower, the LST error decreases. Figures 4-7 show that the LST error is oppositely related to the LSE.

Accuracy Assessment of the Proposed SWA
In order to assess the accuracy of our SWA, we used Modtran 4.0 to simulate the thermal radiance reaching the sensor for an input Mid-Latitude Summer atmospheric profile with known LST and LSE. The simulated radiance was then converted into brightness temperatures for the TIRS bands, and used as inputs to the SWA. The LST estimation errors for 60 different scenarios, featuring different combinations of LST, LSE, and atmospheric water vapor content, are presented in Table 3. The RMSE for the LST estimation errors in Table 3 is 0.93 °C.

Summary
We presented a SWA for Landsat-8 TIRS data. The ability to derive LST temperatures accurately has immediate impacts for potential Landsat-8 users and applications. As data from the new Landsat-8 are distributed freely, its high resolution thermal abilities can allow for new scientific advances in earth science. Although the spatial resolution of TIRS is degraded to 100 m, in comparison with the 60 m resolution of the Landsat-7 ETM+ thermal band, the added value of having two bands is more accurate LST estimations with TIRS than with its predecessor. In addition, the 100 m resolution is sufficient for water consumption measurements over fields irrigated by center pivot systems [9], as well as for other uses over relatively homogeneous areas.
In this paper, we have outlined possible strategies for evaluating the atmospheric water vapor content and LSE. When selecting a data source for estimating the water vapor content, the tradeoff between temporal and spatial precision has to be considered by the users on a case-by-case basis. Strategies to estimate LSE, based on techniques employed by previous sensors, have been suggested. These techniques can still be refined and adapted for OLI's new spectral band configuration, which is different than its predecessors. Therefore, future work should focus on evaluating the best methods of estimating the SWA input parameters (e.g., atmospheric transmittance, LSE), the cross-sensor validation of the input parameters and the resulting LST.
Yevgeny Derimian contributed to the algorithm adaptation, the atmospheric modeling, and interpretation of the modeling results. Arnon Karnieli conceived the idea for the project, and supervised and coordinated the research activity. All authors participated in the editing of the paper.