## Abstract

Land surface temperature (LST) is key parameters in the interaction of land-atmosphere system. This paper proposed a method to inverse LST from multi-temporal thermal infrared remote sensing data based on the theory of split-window algorithm and diurnal temperature cycle model. The new method was validated by a diurnal brightness temperatures data sets corresponding to MSG2-SEVIRI, which was simulated by the atmospheric radiative transfer model MODTRAN 4 with several input parameters under clear sky, including near surface air temperature, atmospheric water, surface temperature and emissivity, and viewing angles, and result showed the root mean square error (RMSE) of LST reaches 1.2K for simulated data and most errors are within ± 2K with accurate parameters input. At the same time, comparison of LST estimated using the proposed method from MSG2-SEVIRI data with that from MOD11B1 V5 product displayed that the RMSE equals to 3.0K and most errors are distributed within ± 6K. But, the method is proposed under no cloudy condition and is tested only in mid-latitude and daytime; more validation should be made in different areas and atmospheric conditions.

©2013 Optical Society of America

## 1. Introduction

Land surface temperature (LST) is an important variable to control fundamental biospheric and geospheric interactions between the Earth’s surface and its atmosphere, and plays many critical roles in estimating soil surface moisture, long-wave radiation, evapotranspiration modeling and so on [1–3]. In many study domain, it is necessary to estimate surface temperature over large spatial and temporal scales, but it is practically impossible to get such information based on ground measurements, while the satellite observations in the Thermal Infra-Red shows to be very attractive since it can give access to global and temporal estimates of LST.

However, the retrieval of the LST from satellite data is a very difficult task due to several reasons as: the tough atmospheric correction due to the complex influence from atmosphere; the hard quantity expression of thermal infra process due to the elaborateness from surface thermal condition; the difficult separation of temperature and emissivity [4]; Seasonal and dynamic states of surfaces and uncertainty of land classification. Up to now, many algorithms for estimating the LST from satellite observations have been proposed and they can be grouped into several categories: the single channel algorithm, the split window algorithm, the multi-channel temperature/emissivity separation method, multi-angle and multi-channel combined algorithm [5].

The single channel method is a simple inversion of the radiative transfer equation providing that the land surface emissivities (LSEs) and the atmospheric profiles are known in advance [6]. The split window method is used to retrieve the LST based on the differential water vapor absorption in two adjacent infrared channels. This method was firstly proposed by McMillin [7] to estimate sea surface temperature from satellite measurements and then was extended to retrieve land surface temperature from thermal infrared data observed by AVHRR, MODIS, Spinning Enhanced Visible and Infrared Imager (SEVIRI) instruments and FengYun geostationary satellite [8–12]. However, coefficients of this method always vary with different study regions, and the method requires a prior knowledge of the pixel emissivity in each TIR channel. The multi-channel temperature/emissivity separation method is proposed by Gillespie et al. to obtain emissivity by using an empirical relationship among five channels’ emissivity and then to estimated LST [13], but its accuracy is low at low emissivity contrast surfaces (e.g. vegetation and water) and also limited by precise atmospheric correction and the empirical relationship between the minimum LSE and the spectral contrast. The multi-angles and multi-channel combined algorithms use the different atmospheric absorption to thermal infrared radiation in different viewing directions. The typical application is realized by Sobrino who proposed a simple double angle model based on Taylor series approximation. The method requires accurate geometric registration and is applied only to homogeneous surfaces [14].

With the launch of geostationary satellite in recent years, such as MSG, FY-2D, the multi-temporal and multi-spectral information can be used to inverse LST [15,16]. Meteosat Second Generation (MSG2) is a new generation of the geostationary satellite developed by European Space Agency (ESA) and EUMETSAT. MSG2’s main payload is the optical imaging radiometer, called Spinning Enhanced Visible and Infrared Imager (SEVIRI). With its 12 spectral channels covering from visible to infrared will provide 48 times more information. The work presented in this paper aims to retrieve diurnal LST cycle from the MSG2 satellite data in two multi-temporal thermal infrared channels (IR1, 10.3-11.3μm and IR2, 11.5-12.5μm) using multi-temporal data and is structured as: at first, a new method was to developed based on split-window algorithm and diurnal LST cycle model in section 2; then, the proposed method is validated using brightness temperature (BT) cycle data simulated from MODTRAN 4 based on field measured LST and near surface air temperature and atmospheric water, at the same time, sensitivity analysis of the method to atmospheric vertical water content (WVC) and noise equivalent temperature difference (NEΔT) are done in section 3; At final, the method is applied to retrieve LST from MSG2/SEVIRI data (denoted hereafter as SEVIRI LST), and compared with LST extracted from the MOD11B1 V5 product (denoted hereafter as MODIS LST) in section 4; the summary and conclusion are given in the final section.

## 2. Theory

#### 2.1 Method

For clear days, the diurnal land surface temperature cycle (DTC) can be described as following [17]:

_{lst}represents approximately the average LST and b represents the half of amplitude, β is the angular frequency and close to π/DD (DD is duration of daytime), td is the time when temperature reaches the maximum T

_{max}, ts is the starting time of temperature decrease and closes to the sunset time, and α is the decay coefficient during the night-time, ${\delta}_{LST}$denotes the error of real and estimated LST with DTC model, and t is time.

Assuming the BT can also be described with DTC model, that is:

_{bt}denotes approximately the average SBT, ${\delta}_{SBT}$denotes the error of real and estimated SBT with DTC model.

The general-split-window (GSW) algorithm proposed by Becker and Li [18] can be used to estimate LST from the two adjacent thermal infrared channels (IR1 and IR2) of MSG2/SEVIRI sensor, which was expressed as:

Where:_{10.8um}and BT

_{12.8um}are the top of atmosphere (TOA) BT measured in channels 9 and 10 of MSG2-SEVIRI, respectively; SBT is the same as that in formula 2, ${\delta}_{split}$denotes the error of real and estimated LST with GSW method.

Combined with Eqs. (1)–(3), formula 4 can be got as following:

Assuming:-P*a_{bt} + c is the function of (a_{bt},θ,wvc) and b-P*b + d is the function of (b,θ,wvc), that is:

Then, formula 6 can be changed to:

According to Eq. (2), β and td can be inversed if more than 6 SBT in a day can be got, therefore, Scale and offset in Eqs. (10)–(12) can be estimated using the least-square method with known variables (time and DBT). Wvc can be estimated using the method developed by Jiang et al. [17] and Li et al. [19] from MSG2-SEVIRI IR1 and IR2 data and VZA can be extracted from sensors. If the specific form of $fun2(b,\theta ,wvc)$can be got, M can be calculated with Scale, further, if the specific form of$fun1({a}_{bt},\theta ,wvc)$ and $fun3(wvc)$are known, a_{lst} can be estimated with M, then, LST in one daytime can be got with Eq. (1). So, the further work is to decide the specific form of all functions and to validate all assumptions.

#### 2.2 Data simulation

Computer simulation is a practical way to simulate diurnal BT in coincidence with the coarse-resolution sensor. The atmospheric radiative transfer model MODTRAN 4 was used to simulate the TOA radiance with the appropriate thermal infrared channel response function of the MSG2-SEVIRI. The model input includes LSTs, LSEs and atmospheric profile. Take into account the angular dependence of the TOA radiance, different view zenith angles (VZA) varying from 0 to 60 with a step of 10 were used in MODTRAN 4 simulations. Diurnal near surface air temperature and atmospheric water of eight cloud free days (28th, 67th, 114th, 170th, 200th, 218th, 243th, and 336th of year) from field measurement in 2006 of YUCHENG field site, ShanDong province as parameters are input to run MODTRAN 4. Considering atmospheric profiles are needed to input as parameters, the other level air temperatures and humidity are interpolated according to two standard atmospheric profiles (mid-latitude summer and mid-latitude winter) while measured near surface air temperature as the first level temperature and measured atmospheric water as total precipitable water. In addition, considering the most land covers, the averaged emissivity varying from 0.9 to 1.0 with a step of 0.01 and the emissivity difference from −0.02 to 0.02 with a step of 0.01, were used in our simulation.

Then for a given LST, TOA radiance of different channel can be simulated by MODTRAN 4, corresponding BT can be got with formula 13. In total, diurnal LST and related TOA measured brightness temperature were obtained for different situation.

^{−5}mW (cm

^{−1})

^{−4}m

^{−2}sr

^{−1},${C}_{2}^{\text{'}}$ = 1.43877 K (cm

^{−1})

^{−1}, L is the TOA radiance, i is the channel, Vc is the central wavenumber of the SEVIRI channel, and A and B are coefficients listed in Table 1.

## 3. Results and analysis

#### 3.1 DTC model for LST and BT

Simulated BT of MSG2-SEVIRI and field measurements were chosen to validate the DTC model described in Eqs. (1) and (2). Figure 1 shows the predicted and measured LST of eight days at YuCheng site in 2006, with root mean square errors (RMSE) ranging from 0.3K to 1K. The result indicated that DTC model performed well in real situation. Figure 2 displays the histogram of errors between fitted and simulated average of channels BT1 and BT2 for 168800 situations. The RMSE of all BTs is 0.4K and the errors between estimated and real BT are within ± 2.5K, and most errors are within ± 0.5K, suggesting that the DTC model can describe the diurnal BT evolution for a clear day with a high accuracy.

#### 3.2 LST retrieval with split-window algorithm

As is shown in Eq. (3), LST can be well estimated
with known parameters A, P and M. For a given LST, channel T_{i} can be determined
according to Eq. (13) and TOA radiance got by
running MODTRAN 4, the coefficients A, P and M in Eq.
(3) can be obtained through statistical regression method. But, A, P and M is the
function of variables including emissivity, view angle and atmospheric water. In the process of
estimating the coefficients, we found the parameter P varies from 0.98 to 1.1, and the average
of P is 1.016, so, P is fixed as 1.016 in the process of regression. Figure 3 gives the histograms of errors of estimated and real LSTs, the
figure indicates that all differences of actual and estimated LST are within ± 4K and
most errors are distributed within ± 1K, and the RMSE is 1.05K. So, we think the GSW algorithm can be used to calculate the LST when P is fixed as
1.016.

**3.3 DTC model for ${\delta}_{LST}-P*{\delta}_{SBT}\text{-}{\delta}_{splitw}$**

In the work, we assumed daily ${\delta}_{LST}-P*{\delta}_{SBT}\text{-}{\delta}_{splitw}$ (SE) evolution can be expressed as cosine function just described as Eq. (5). The histogram of errors between real SE and that fitted with Eq. (5) is displayed in Fig. 4. From Fig. 4, one can find that the RMSE of estimated and actual SE is 0.09K and most errors are within ± 0.1K (the range of SE is from −2.994K to 2.71K). So, a conclusion that daily SE evolution can be described with cosine function is drawn.

#### 3.4 The specific form of function (1-3)

Keep in mind that the specific form of Eqs. (7) and
(8) should be given out if LSTs can be
estimated. Just as described in 2.2, for given LST, BT of TOA for channel 9 and 10 can be
calculated by running MODTRAN 4, so, coefficient a_{bt} in Eq. (2) can be estimated using more than six BTs. In
addition, parameters of c and d in Eq. (5) can
also inversed using more than two SEs. In combination with VZA and WVC, the specific form of
fun1 and fun2 can be got using the least-square method and the forms are just listed as Eqs. (14) and (15), fun3 (Eq. (16)) can be obtained
when deciding the coefficients of Eq. (3) with
GSW algorithm. Figures 5 and 6 give histogram distribution of all the errors for Eqs. (14) and (15) respectively. From Fig. 5, one can see that RMSE of errors
between real and estimated $\text{-}P*{a}_{bt}+c$ is 0.95K and most errors are distributed ± 2K; From Fig. 6, one can see that RMSE of errors between real and
estimated b-P*b + d is 0.17K and most errors are distributed ± 0.2K. Because
$\text{-}P*{a}_{bt}+c$ ranges from 265.45K to 305.103K, the RMSE from Eq. (14) can be omitted and Eq. (14) can be used to substitute
$\text{-}P*{a}_{bt}+c$, while b-P*b + d ranges from −0.057K to 1.158K, RMSE from
Eq. (15) is relatively bigger than that from
Eq. (14), the result perhaps affect the
accuracy of estimation for LST.

#### 3.5 Determination of WVC

According to Li et al. [19], the atmospheric WVC can be derived by the use of the transmittance ratio of split-window channels (assuming emissivity of two channels is identical),

Where:_{1}and c

_{2}can be derived as functions of VZA just listed in Eq. (18), i and j are the split-window channels, the subscript k denotes pixel k, and the ${\overline{T}}_{i}$/ ${\overline{T}}_{j}$are the TOA mean (or the median) channel brightness temperatures of the N neighboring pixels considered for channels i and j, respectively.

#### 3.6 Results from multi-temporal inverse method for simulated data

According to the proposed method, LST are computed and compared to the real LST for simulated data, the results are displayed in Fig. 7. It is noted that the RMSE is 1.18K and most errors concentrated on ± 2K and maximum of absolute errors is less 5K from Fig. 7. Then, we can draw that LST can be estimated with the proposed method.

#### 3.7 Sensitivity analysis to the atmospheric WVC

It is well known that the WVC in the atmosphere is not easily determined from satellite data. In order to see how significant the effect of the uncertainty of the WVC on the retrieval of LST with proposed algorithm, the $\pm 20\%$errors of the WVC are added to the real WVC, then LSTs are estimated with the wrong WVC. Figure 8 gives an example of the uncertainty of the WVC. From Fig. 8, one can see, when we estimate LST using adding 10% errors of WVC, RMSE of LST is 1.56K and most errors are within ± 3K; while we estimate LST using adding 20% errors of WVC, RMSE of LST reach 2.18K, but most errors are still concentrated on within ± 3K and maximum of absolute errors is less than 5K.

#### 3.8 Sensitivity analysis to noise equivalent temperature difference (NEΔT)

In order to see how significant the effect of the instrumental NEΔT on the retrieval of LST using proposed method, 0.1 K, 0.2 K and 0.5 K are, respectively added to the TOA BT in Eq. (2). Then we estimate the LST using proposed algorithm with the noised TOA BT.

Figure 9 gives an example of the uncertainty of NEΔT. Comparing the actual LST with the estimated LST, the RMSE = 1.29K for NEΔT = 0.1k, the RMSE = 1.31K for NEΔT = 0.2K and the RMSE = 1.33K for NEΔT = 0.5K, and most errors are concentrated on within ± 3K.

## 4. Application to actual MSG-SEVIRI satellite data

The objective of the present work is to estimate the LST from MSG2-SEVIRI operational geostationary meteorological satellite data for cloud-free skies. Figure 10 gives the comparison of LST retrieved using proposed method from MSG2-SEVIRI data with MOD11B1 V5 product for part of southeast Europe and west Asia with geospatial coverage of latitude 30°N to 50°N and longitude 25°E to 45°E on August 5, 2006 at 10:30 local time for no cloudy condition. The model inputs are the TOA BTs, VZA and WVC. The TOA BTs and VZA are directly extracted from MSG-SEVIRI satellite data and the WVC are obtained according to Eq. (17).

Because LST varies strongly in space and time, the comparison of LST between two different sensors must be at the same place, and the time is identical. Take into account that these two satellite-derived products have different spatial resolutions, they must be aggregated to the same spatial resolution for comparison. In this study, the SEVIRI LST are aggregated to the same spatial resolution of the MOD11B1 V5 LST products using nearest neighbor method and SEVIRI LST are calculated using DTC model according to the MODIS overpass time from MOD11B1 V5 LST products.

From Fig. 10, one can see that the RMSE reaches 3.0K and the scatters are distributed near 1:1 line. The errors may be caused by the inaccuracy of proposed method and the errors of MODIS LST products and misregistration of the two products.

## 5. Conclusions

In this paper, we have proposed a new method to estimate land surface temperature (LST) based on diurnal temperature cycle model and split-window method from Meteosat Second Generation (MSG2) Data in two thermal infrared channels IR1 (10.3-11.3μm) and IR2 (11.5-12.5μm).

Taking account of the fact that there are few databases of in situ LST measurements in coincidence with the coarse-resolution sensor, using air temperature, atmospheric water and land surface temperature measured at Yuncheng, Shandong in 2006 over a wide range of surface conditions as input parameters, a databases responding to MSG2-SEVIRI were built with an accurate atmospheric radiative transfer model MODTRAN 4. The simulation analysis showed that the LST could be estimated by the proposed algorithm with RMSE less than 1.5 K and most errors concentrated on ± 2K and maximum of absolute errors is less 5K for mid-latitude region and for daytime under clear skies

As the proposed algorithm requires atmospheric vertical water content (WVC) as model input, the atmospheric WVC can be determined using the method developed by Li et al. [19].

In addition, the sensitivity and error analyses in term of the uncertainty of WVC as well as the instrumental noise were also performed in this work. The results show that the accuracy of retrieval LST can be affected by 0.01K for NEΔT = 0.1 K, by 0.13K for NEΔT = 0.2 K, and by 0.15K for NEΔT = 0.5 K for RMSE of all data; and the effect of the uncertainty of the WVC on the retrieval LST could be around 0.4 K and 1.0K respectively for adding 10% and 20% error to real WVC.

Last, LST were calculated using MSG2-SEVIRI data with new method and compared with MOD11 temperature products. The results showed that the trend is similar between LST estimated with proposed method for MSG2-SEVIRI and LST from MODIS products, and the RMSE reaches 3K and most errors concentrated on ± 6K. It also illustrates the method can be used to estimated LST for MSG2-SEVIRI data.

Need to be noticed, the method was validated only in mid-latitude region and daytime under clear skies. In other region or other condition, the method still needed to be checked. In addition, accuracy of the method is limited by initial parameters values, prior knowledge and inverse algorithm, these questions will be researched in next step.

## Acknowledgments

This work was supported by the National Natural Science Foundation of China under Grant 41271381/40971199. We also thank EUMETSAT for providing the MSG-2/SEVIRI data.

## References and links

**1. **J. C. Price, “The potential of remotely sensed thermal infrared data to infer surface soil moisture and evaporation,” Water Resour. **16**, 787–795 (1990).

**2. **B.-H. Tang and Z.-L. Li, “Estimation of instantaneous net surface longwave radiation from MODIS cloud-free data,” Remote Sens. Environ. **112**(9), 3482–3492 (2008). [CrossRef]

**3. **Z.-L. Li, R. L. Tang, Z. M. Wan, Y.-Y. Bi, C. H. Zhou, B.-H. Tang, G. J. Yan, and X. Y. Zhang, “A review of current methodologies for regional evapotranspiration estimation from remotely sensed data,” Sensors (Basel) **9**(5), 3801–3853 (2009). [CrossRef] [PubMed]

**4. **Z.-L. Li, H. Wu, N. Wang, S. Qiu, J. A. Sobrino, Z. Wan, B.-H. Tang, and G. J. Yan, “Land surface emissivity retrieval from satellite data,” Int. J. Remote Sens. **34**(9–10), 3084–3127 (2013). [CrossRef]

**5. **Z.-L. Li, B.-H. Tang, H. Wu, H. Z. Ren, G. J. Yan, Z. M. Wan, I. F. Trigo, and J. A. Sobrino, “Satellite-derived land surface temperature: Current status and perspectives,” Remote Sens. Environ. **131**, 14–37 (2013). [CrossRef]

**6. **J. Cristóbal, J. C. Jiménez-Muñoz, J. A. Sobrino, M. Ninyerola, and X. Pons, “Improvements in land surface temperature retrieval from the Landsat series thermal band using water vapor and air temperature,” J. Geophys. Res. **114**(D8), D08103 (2009). [CrossRef]

**7. **L. M. McMillin, “Estimation of sea surface temperature from two infrared window measurements with different absorption,” J. Geophys. Res. **80**(36), 5113–5117 (1975). [CrossRef]

**8. **Z. Wan and J. Dozier, “A Generalized split-window algorithm for retrieving land-surface temperature from space,” IEEE Trans. Geosci. Rem. Sens. **34**(4), 892–905 (1996). [CrossRef]

**9. **B.-H. Tang, Y. Y. Bi, Z.-L. Li, and J. Xia, “Generalized Split-Window algorithm for estimate of land surface temperature from Chinese geostationary FengYun meteorological satellite (FY-2C) data,” Sensors **8**(2), 933–951 (2008). [CrossRef]

**10. **J. A. Sobrino and M. Romaguera, “Land surface temperature retrieval from MSG1-SEVIRI data,” Remote Sens. Environ. **92**(2), 247–254 (2004). [CrossRef]

**11. **T. Schmugge, A. French, J. C. Ritchie, A. Rango, and H. Pelgrum, “Temperature and emissivity separation from multispectral thermal infrared observations,” Remote Sens. Environ. **79**(2–3), 189–198 (2002). [CrossRef]

**12. **A. J. Prata and C. M. R. Platt, “Land surface temperature measurements from the AVHRR,” in *Proceedings of the 5th AVHRR Data Users Conference*, Tromso, Norway, June 25–28, EUM P09 (1991), pp. 433–438.

**13. **A. Gillespie, S. Rokugawa, T. Matsunaga, T. Matsunaga, and J. S. Cothern, “A temperature and emissivity separation algorithm for Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) images,” IEEE Trans. Geosci. Rem. Sens. **36**(4), 1113–1126 (1998). [CrossRef]

**14. **J. A. Sobrino, Z.-L. Li, M. P. Stoll, and F. Becker, “Multi-channel and multi-angle algorithms for estimating sea and land surface temperature with ATSR data,” Int. J. Remote Sens. **17**(11), 2089–2114 (1996). [CrossRef]

**15. **C. Gao, X. Jiang, H. Wu, B.-H. Tang, Z. Li, and Z.-L. Li, “Comparison of land surface temperatures from MSG-2/SEVIRI and Terra/MODIS,” Int. J. Remote Sens. **6**(1), 063606 (2012).

**16. **G. M. Jiang and Z.-L. Li, “Split-window algorithm for land surface temperature estimation from MSGI-SEVIRI data,” Int. J. Remote Sens. **29**(20), 6067–6074 (2008). [CrossRef]

**17. **F. M. Göttsche and F. S. Olesen, “Modeling of diurnal cycles of brightness temperature extracted from METEOSAT data,” Remote Sens. Environ. **76**(3), 337–348 (2001). [CrossRef]

**18. **F. Becker and Z.-L. Li, “Toward a local split window method over land surface,” Int. J. Remote Sens. **11**(3), 369–393 (1990). [CrossRef]

**19. **Z.-L. Li, L. Jia, Z. B. Su, Z. M. Wan, and R. Zhang, “A new approach for retrieving precipitable water from ATSR2 split-window channel data over land area,” Int. J. Remote Sens. **24**(24), 5095–5117 (2003). [CrossRef]