Sharpening method of satellite thermal image based on the geographical statistical model

Pengcheng Qi; Shixiong Hu; Haijun Zhang; Guangmeng Guo

doi:10.1117/1.JRS.10.025013

12 May 2016 Sharpening method of satellite thermal image based on the geographical statistical model

Pengcheng Qi, Shixiong Hu, Haijun Zhang, Guangmeng Guo

Author Affiliations +

Journal of Applied Remote Sensing, Vol. 10, Issue 2, 025013 (May 2016). https://doi.org/10.1117/1.JRS.10.025013

Abstract

To improve the effectiveness of thermal sharpening in mountainous regions, paying more attention to the laws of land surface energy balance, a thermal sharpening method based on the geographical statistical model (GSM) is proposed. Explanatory variables were selected from the processes of land surface energy budget and thermal infrared electromagnetic radiation transmission, then high spatial resolution (57 m) raster layers were generated for these variables through spatially simulating or using other raster data as proxies. Based on this, the local adaptation statistical relationship between brightness temperature (BT) and the explanatory variables, i.e., the GSM, was built at 1026-m resolution using the method of multivariate adaptive regression splines. Finally, the GSM was applied to the high-resolution (57-m) explanatory variables; thus, the high-resolution (57-m) BT image was obtained. This method produced a sharpening result with low error and good visual effect. The method can avoid the blind choice of explanatory variables and remove the dependence on synchronous imagery at visible and near-infrared bands. The influences of the explanatory variable combination, sampling method, and the residual error correction on sharpening results were analyzed deliberately, and their influence mechanisms are reported herein.

1. Introduction

Satellite-based thermal infrared remote sensing is the main method used to acquire distributed land surface thermal information. However, the spatial resolution of existing satellite-based remote sensing thermal infrared images is generally low. For example, the spatial resolution of thermal infrared imageries collected by the moderate resolution imaging spectroradiometer (MODIS) and the advanced very high resolution radiometer (AVHRR) is $\sim 1000 m$ , and the resolution of thermal infrared images collected by the Geostationary Operational Environmental Satellite Imager is 4000 m. Such low spatial resolution significantly restricts the usefulness of these images in drought monitoring and other agricultural and forestry practices,¹^–⁴ particularly in mountainous regions with complex landscape. High-resolution images can be obtained through sharpening of coarse resolution images under the premise that most physical information is maintained. Sharpening is an effective method used to solve the problem of low spatial resolution of thermal images.³ Thermal image sharpening, which is also called thermal information downscaling or subpixel temperature estimation, has been an important topic in the international remote sensing field in recent years.¹^–¹⁵ The two following paragraphs will provide an analysis of previous studies from the study area types and explanatory variables.

The studied areas in previous studies can be classified into agricultural areas (e.g., see Refs. 2, 3, and 5), urban areas (e.g., see Refs. 11, 16 17.18.–19), and areas with complex landscapes (e.g., see Refs. 8 and 20). The performance of a certain method in different landscapes was reported.¹⁴^,²¹ Most of these areas featured flat terrain. Thus far, studies on thermal image sharpening in mountainous areas have been relatively rare.¹ Functioning as complicated energy budget systems, mountainous regions refer to regions that are characterized by a certain altitude, relative height, and gradient.²² Different from other geographical units, mountainous regions have their own laws in weather and climate, including: (1) the temperature, water vapor pressure, solar radiation, rainfall, and other factors vary with the change in altitude.²³^,²⁴ (2) Solar radiation, wind speed, and other factors vary with the difference in slope aspect and gradient.²³^–²⁶ (3) Different terrains exert differentiated aggregating effects on water, mineral substances, plant litter, and other substances.²² (4) Terrain has a significant effect on plant species and growth status.²² The aforementioned influences have direct or indirect effects on the surface energy balance (SEB) in mountainous regions, and incoming and outgoing items of SEB vary with topographical change. As a state variable in the surface energy budgeting process, land surface temperature (LST) in mountainous regions results in more obvious spatial variation than that in flat areas and has intrinsic relationships with more environmental factors.

At the beginning of this research field, disTrad (later known as TsHARP) was the typical thermal sharpening method.²⁷ It used the normalized difference vegetation index (NDVI) as the explanatory variable (the term of “explanatory variable” mentioned in the paper is basically tantamount to “kernal” and “scale factor” involved in similar literature, which have been summarized in Ref. 1). Subsequent scholars improved TsHARP and used the transformed NDVI as explanatory variables [e.g., Refs. 5 and 28 used vegetation coverage (VC) as explanatory variable] or introduced new explanatory variables (e.g., albedo¹¹ and normalized multiband drought index²⁰). Some scholars still used NDVI as explanatory variables but introduced new regression approaches or interpolation methods (such as least-median square regression downscaling,² regression-kriging,³ the combination of TsHARP and thin plate spline,²¹ and so on). An emerging tendency for this research field is to use multiple factors as explanatory variables and adopt machine learning algorithm for modeling.⁸^,¹⁴ Most previous studies in this field used merely the optical image, without further exploring the role that external data [such as digital elevation model (DEM), meteorological data, microwave data, and so on] play on thermal image sharpening (see the discussion in Ref. 1). As stated earlier, LST in mountainous regions features special laws in spatial variation and is subject to many environmental factors. Viewed from this perspective, it is impossible to satisfyingly sharpen thermal images in mountainous regions if only explanatory factors from previous studies are used. Still further studies should be conducted on issues like the selection of explanatory variables suited to thermal image sharpening in mountainous regions and the ways for spatially depicting these explanatory variables at high resolution.

To improve the effectiveness and accuracy of the thermal sharpening technique in mountainous regions, this study proposed a method based on the geographical statistical model (GSM). Given that LST is the result of land SEB,²⁹^,³⁰ it is determined by land surface properties and energy balance together. This study selected various parameters in the SEB process as the explanatory variables. These variables were spatially simulated to generate their raster layers. On this basis, the local adaptation statistical relationship between low-resolution top-of-atmosphere brightness temperature (BT) and explanatory variables, i.e., the GSM, was established through multivariate adaptive regression splines (MARS).³¹ Then, the GSM was applied to high-resolution explanatory variables, thereby obtaining high-resolution BT image. In this paper, the GSM method was assessed by using it in sharpening a simulated coarse resolution (1026 m) BT image to fine resolution (57 m). A series of experimental designs was made in this study to analyze the influences of explanatory variable combination, sampling method, and residual error correction on sharpening results, as well as their influence mechanisms.

2. Method

2.1.

Study Area and Data Source

The study area, which is located west of Oregon in the northwestern United States, is an approximate rectangular region between the west longitude of 123 to 124 deg and the northern latitude of 43 to 44 deg. The elevation of the region is 0 to 1299 m (335 m on average), and the slope degree is 0 to 54.5 deg (15.9 deg on average). The main land cover types (LCTs) in this region include forest, cutover land, grassland, farmland, and water bodies. This region was selected as the study area because (1) it is a typical mountainous region, (2) it has rich LCTs and rather complicated landscape, and (3) sufficient data were available in this region, especially DEM with high spatial resolution and accuracy.

The Landsat enhanced thematic mapper plus (ETM+) remote sensing images used, whose path and row number was P46R30, were generated on September 7, 1999. These images were L1T terrain correction data products downloaded from United States Geological Survey (USGS) Earth Explorer. The resolution of these images is 57 m for thermal band and 28.5 m for visible and near-infrared (Vis-NIR) bands. The thermal band was used to simulate a coarse resolution (1026 m) thermal image (see Sec. 2.2 for details). In the following process, the coarse resolution thermal image would be sharpened back to 57-m resolution, and then the sharpening result would be verified based on the initial 57-m resolution thermal image. Using ETM+ data in the study is justified by the following considerations: first, it provides a resolution (57 m) in the thermal infrared band much higher than most legacy and current satellite-borne sensors do. As to regions having a strong spatial heterogeneity in LCT, a more practical value may be provided if the study effort focuses on how to sharpen the thermal images from a 1000- to 60-m resolution level. Second, from May 2003 to the present, about a quarter of the data in ETM+ images are missing because of a scan line corrector (SLC) failure, but fortunately, the remaining three quarters are accurate. Using these data to evaluate the thermal image sharpening effect creates a new application of these data. (As a basic study, this paper used the ETM+ data before the SLC failure; however, the approach developed in such a way is likewise suitable for post-SLC-failure data).

A DEM, with a spatial resolution of 10 m and a mean square error of $< 5 m$ , was obtained from the National Elevation Dataset of USGS. Precipitation data 2 weeks before the imaging date and air temperature (AT) and atmospheric precipitable water data on the imaging date were collected from the U.S. National Climatic Data Center of National Oceanic and Atmospheric Administration. The solar radiation data of the stations on the imaging date and the atmospheric data related with solar radiation were obtained from the solar radiation database of the American National Renewable Energy Laboratory (ANREL). The National Land Cover Database (NLCD; Edition 2) of USGS was used as a reference in the classification of remote sensing images.

2.2.

Preparation of Images for Sharpening

The sharpening object in previous studies can be divided into three classes: digital number (DN) images (e.g., Ref. 6), top-of-atmosphere BT images (e.g., Refs. 11 and 14), and land surface real temperature images (e.g., Refs. 2, 5, and 10). All these classes can be called a “thermal image.” Because there are various retrieval methods from DN or top-of-atmosphere BT to land surface real temperature (e.g., Refs. 31 32.33.34.35.–36) and different methods will get different results,³⁷^,³⁸ the relatively original data, top-of-atmosphere BT images, are better sharpening objects compared to the land surface real-temperature images. This study used the top-of-atmosphere BT image as a sharpening object.

First, the original DN image of ETM+ Band6 with 57-m resolution was converted into the radiance image at the entrance pupil ( $L_{6}$ ) according to the gain and bias parameters in its header file. Then, $L_{6}$ was converted into a top-of-atmosphere BT image (abbreviated as BT57) according to the Planck formula. Finally, the BT image with 1026-m resolution (abbreviated as BT1026) was simulated based on BT57. In the following procedure (see Sec. 2.7), BT1026 would be used as sharpening object. The reasons that we selected 1026 m as the source resolution were as follows: (1) 1026 m is integral multiple of the resolution of BT57, which made the aggregating process relatively easy; (2) 1026 m is approximately equal to the spatial resolution of thermal infrared images collected by the MODIS or AVHRR; via the GSM’s performance in sharpening the simulated 1026-m resolution image, its performance in sharpening MODIS or AVHRR thermal image can be evaluated. The procedure of generating BT1026 was as follows:¹⁴ every $18 \times 18$ BT57 pixel aggregates into one pixel with a resolution of 1026 m, and the pixel value was calculated from the values of the corresponding $18 \times 18$ BT57 pixels with the equation:

Eq. (1)

BT 1026 (m, n) = {[1 / (18 \times 18)] \sum_{j = 1}^{18} \sum_{i = 1}^{18} {[B T 57 (i, j)]}^{4}}^{0.25} .

2.3.

Selection of Explanatory Variables and Establishment of Their Raster Layers

LST is related to many natural factors, and a good understanding of the intrinsic relationships between LST and its determinants is the basis of downscaling. Since the basic idea of this method was to simulate LST as the result of SEB, the main environmental factors of the SEB were selected as explanatory variables. These environmental factors were judged on the basis of previous studies focusing on SEB.²⁹^,³⁰^,³⁹^,⁴⁰ In addition, some other explanatory variables were proxy variables, which were used to replace those environmental factors proved by previous studies as directly related to the BT, whereas their raster layers were difficult to obtain. This process could guarantee the simplicity and feasibility of this method. In this study, 13 kinds of environmental factors were selected according to the aforementioned principles. The specific selection basis and the raster layer establishment method of these factors are detailed in Secs. 2.3.1–2.3.6

2.3.1.

Total solar radiation

In mountainous regions, the spatial variation of total solar radiation ( $I$ ) received by the land surface caused by topographic relief²⁵^,²⁶ is a significant cause of the small-scale spatial variation of surface temperature in mountainous regions. For instance, a simulation study based on an ecohydrological model shows that for the climate conditions considered in that study, the incoming solar radiation is one of the major factors controlling LST spatial distribution.⁴¹ In addition, because there is a time lag between the surface temperature intraday change and $I$ intraday change, the effect of $I$ in some period before imaging on thermal image is possibly superior to the effect of $I$ simultaneous with imaging on thermal image;⁴² Ref. 42 shows that under the conditions of their experiment, a time lag of 1 h considerably increases the linear relation between vegetation canopy temperature and local insolation angle. Therefore, $I$ in several tiny periods before imaging was selected as an alternative explanatory variable and the simulation method is as follows:

1. Tiny period definition: to facilitate calculation, we divided the period from 8:06:00 to 10:36:00 (i.e., 2.5 h before ETM+ imaging) into 30 tiny periods, i.e., each tiny period is 5 min. The tiny periods were named 5-min-periods. In each 5-min-period, all factors related to solar radiation were assumed to remain unchanged and be the status at the 5-min-period’s midpoint.
2. Astronomical solar radiation ( $I_{0}$ ): it was determined whether each DEM grid cell is irradiated directly in a specific 5-min-period(s) with the method indicated in Ref. 43. For each DEM grid cell that can be irradiated directly, $I_{0}$ was calculated as follows:
Eq. (2)
$I_{0} = T \times S_{Const} \times E_{0} \times \cos θ,$
where $T$ is the length of the 5-min-period (300 s); $S_{const}$ is solar constant with the value of $1367 W / m^{2}$ ; $E_{0}$ is the correction parameter of earth–sun distance (unitless) and can be calculated from the day number of the year;²⁵ and $\cos θ$ is the cosine of the insolation angle on a tilted surface, which can be calculated as follows:⁴³
Eq. (3)
$\cos θ = (\sin φ \cos α - \cos φ \sin α \cos β) \sin δ + (\sin φ \sin α \cos β + \cos φ \cos α) \cos δ \cos ω + \sin α \sin β \cos δ \sin ω,$
where $φ$ is latitude (rad), $α$ slope degree (rad), and $β$ slope aspect (rad); $δ$ is solar declination (rad) and can be calculated from the day number of the year;²⁵ and $ω$ is the solar hour angle (rad).
3. Direct solar radiation ( $I_{dir}$ ): $I_{dir}$ of each grid cell that is irradiated directly in clear sky condition was calculated according to Iqbal C model²⁵ as follows:
Eq. (4)
$I_{dir} = 0.9751 I_{0} τ_{r} τ_{a} τ_{o_{3}} τ_{w} τ_{g},$
where the factor 0.9751 is included because the spectral interval considered is 0.3 to $3 μ m$ ; $τ_{r}$ , $τ_{a}$ , $τ_{o}$ , $τ_{w}$ , and $τ_{g}$ are Rayleigh scatting transmittance, aerosol scatting transmittance, ozone absorption transmittance, water vapor absorption transmittance, and mixed gas absorption transmittance, respectively. Calculation of these five transmittance parameters was based on the relevant formula in Ref. 25. Among critical parameters needed in the calculation, atmosphere moisture content and aerosol optical thickness came from the solar radiation database of ANREL, and vertical thickness of the ozone layer was calculated with latitude and date according to Ref. 44.
4. Scatting radiation: First, solar scattering radiation in the clear sky ( $I_{dif}^{H}$ ) was calculated on the assumption of “level ground surface without obstruction around” on basis of Iqbal C model.²⁵ Then, solar scattering radiation in clear sky ( $I_{dif}$ ) was calculated on the condition of “possibly inclined ground surface with potential obstruction around toward sky” according to Ref. 43.
5. $I$ : $I$ is the sum of $I_{dir}$ and $I_{dif}$ . $I$ of six 5-min-periods before imaging was used as alternative explanatory variables, which are denoted as $I_{dur_d}$ ( $d = 1,2 \dots 6$ ). In addition, every three adjacent 5-min-periods were taken as one group, and their $I$ was summed up generating eight $I_{DUR_D}$ ( $D = 1,2 \dots 8$ ), which were used as alternative explanatory variables also. In other words, two kinds of $I$ with different duration were generated including $I_{dur_d}$ with the period of 5 min and $I_{DUR_D}$ with the period of 15 min. The initial resolution of $I_{dur_d}$ and $I_{DUR_D}$ layers was 10 m depending on DEM used. We converted these layers to 57-m resolution layers via the following process: at first, via a bilinear resample method, we converted these layers into 11.4-m ones (11.4, a number close to 10, was obtained by dividing 57 by integer 5). Then, we converted layers of 11.4-m resolution into 57-m resolution layers using aggregation method, which meant that each $5 \times 5 pixel$ set was grouped into one pixel, whose average formed the value in their corresponding pixel in the new layers.

2.3.2.

Precipitation before imaging and topographic wetness index

Soil moisture affects the plant canopy temperature by influencing the transpiration intensity of vegetation and has a profound effect on the energy distribution ratio to sensible heat flux and latent heat flux.³⁰^,³⁹^,⁴⁵^,⁴⁶ Generally, the more moist the soil, the more evaporation removes heat. Therefore, soil moisture is an indispensable factor of thermal image modeling.⁴¹ Because it is difficult to simulate soil moisture accurately, precipitation data 2 weeks before the imaging date (PBI) and topographic wetness index (TWI) were used as proxy variables of soil moisture.

Spatial interpolation was conducted to PBI of 37 meteorological stations in the study area and its surrounding areas by using the common Kriging method, thereby obtaining 57-m resolution raster data of PBI.

The 57-m resolution layer of TWI was generated via the following process: first, using the bilinear resample method, a 10-m-resolution DEM was converted to a 11.4-m-resolution DEM. Second, the latter was depression filled. Third, flow direction was calculated using the single-flow direction algorithm⁴⁷ on the basis of the depression-filled DEM, thereby deriving the flow accumulation ( $A$ ). Thus, TWI could be calculated as follows:

Eq. (5)

TWI = \ln (A / \tan α),

where

α

is the slope degree. Finally, the TWI layer of 11.4-m resolution was converted to a 57-m-resolution layer using the aggregation method described in Sec. 2.3.1.

2.3.3.

Vegetation index and vegetation coverage

Transpiration, which is the significant expenditure item in the energy balance of vegetation canopy, controls vegetation canopy temperature.²⁹ Transpiration is largely determined by vegetation quantity and flourishing degree. Therefore, the terrain-corrected NDVI (NDVI′) and VC were selected as the explanatory variables. Their computation methods are described in the following.

The original DN image of the third and fourth bands was converted into the radiance image at the entrance pupil. On this basis, the top-of-atmosphere reflectance was calculated. Then, the $C$ calibration model was used for topographical correction, as follows:⁴⁸^,⁴⁹

Eq. (6)

ρ_{k}^{'} = ρ_{k} \times (\cos θ_{z} + C_{k}) / (\cos θ + C_{k}),

where

ρ_{k}^{'}

is the reflectance of the

k

’th band after topographical correction and

ρ_{k}

is the reflectance of the

k

’th band before topographical correction, in which

k

equals 3 or 4;

θ

is the incidence angle and

θ_{z}

is the solar zenith angle. With respect to the

k

’th band,

C_{k}

was a constant, which was determined by the method of exhaustion, as follows. (1)

C_{k}

was supposedly an arbitrary number from 0 to 3 (

step length = 0.1

); substitute

C_{k}

into Eq. (6) for automatic batch computation, thereby obtaining 30

ρ_{k}^{'}

layers; (2) Based on the historic high-resolution remote sensing images of Google Earth, a total of 10 “observation blocks” were recognized in the study area through visual interpretation. In each “observation block,” land-cover types and vegetation growth at different slope directions were basically the same (the forest managers mainly adopted even-aged regeneration method to manage the forests; therefore, the forests on the same stand belonged to the same age class and were managed with the same measure; this facilitated the selection of observation blocks). (3) In these observation blocks, the topographical correction results (elimination degree of the illumination effect) of

ρ_{k}^{'}

images that were produced by different

C_{k}

were estimated visually;

C_{k}

, which could produce best visual effect, was considered the proper

C_{k}

. Based on these

C_{k}

,

ρ_{k}^{'}

values were calculated. Then, the raster layer of NDVI′ was calculated from

ρ_{k}^{'}

. The NDVI′ layers’ resolution was converted from 28.5 to 57 m using the aggregation method.

The raster layer of VC was calculated by using the 57-m resolution NDVI′, as follows:

Eq. (7)

VC = (NDV I^{'} - {NDVI}_{S}^{'}) / ({NDVI}_{V}^{'} - {NDVI}_{S}^{'}),

where

{NDVI}_{S}^{'}

and

{NDVI}_{V}^{'}

are NDVI′ values of the bare soil pixel and the pixel that is completely covered by vegetation, respectively.

2.3.4.

Emissivity, water vapor content of air column at imaging, and air temperature at imaging

Given that the sharpening object is a top-of-atmosphere BT image, according to the radiative transfer equation,³²^,³⁸ atmospheric transmissivity, atmospheric downward and upward long-wave radiations, and emissivity ( $ϵ$ ) should be considered. Their simulation or proxy methods are as follows:

1. Although atmospheric transmissivity of the thermal infrared band is influenced by many factors (e.g., air pressure, aerosol content, and content of various atmospheric components), its spatial variation mainly depends on the atmospheric water vapor content (AWVC).⁵⁰ Therefore, the AWVC during imaging was used as the proxy factor of atmospheric transmissivity. The atmospheric precipitable water data of 19 meteorological stations in the study area and its surrounding areas were used to generate the 57-m-resolution AWVC raster data on the basis of the macrofactor regression and residual interpolation method.⁵¹
2. The AT at imaging was used to represent the atmospheric downward and upward long-wave radiations indirectly. Based on the AT data of 12 meteorological stations in the study area and its surrounding areas, 57-m resolution AT raster data were generated by using the combined method.⁵²
3. Based on VC, 57-m resolution $ϵ$ raster data were generated by using the method introduced in Ref. 37.

2.3.5.

Profile curvature, slope degree, and elevation

Profile curvature (PC) affects the acceleration or deceleration of flow across the surface. Many environmental factors, such as soil thickness, litter thickness, and soil humidity, have a certain degree of correspondence with PC.⁵³^,⁵⁴ Slope degree ( $α$ ) and elevation (ELEV) were employed as proxy variables of factors, which are related to terrain but unknown currently and their raster data were generated. The 57-m-resolution layers of these three variables were generated via the following process: first, via the bilinear resample method, we converted the 10-m-resolution DEM to a 11.4-m-resolution DEM. Then, we generated 11.4-m-resolution layers for PC and $α$ from the 11.4-m-resolution DEM using corresponding tools in ArcGIS. Finally, we converted these layers of 11.4-m resolution to layers of 57-m resolution using aggregation method mentioned in Sec. 2.3.1.

2.3.6.

Land cover type

Aside from VC, vegetation index, soil moisture, and emissivity, several other thermodynamic factors of the land surface (e.g., specific heat, soil density, thermal conductivity, surface roughness, and albedo) also affect surface temperature. Unfortunately, the spatial distributions of these factors are difficult to describe accurately. Previous studies proved that there is a significant relationship between these factors and LCT.⁵⁵^–⁵⁷ Therefore, the spatial variation of the LCT is one of the major factors causing the spatial variation of LST and SEB.⁴¹^,⁵⁸^,⁵⁹ Reference 41 shows that in addition to solar radiation, the land cover variability is an another major factor causing LST spatial variation. The procedure of generating the LCT raster layer is as follows: first, clustering was conducted with unsupervised classification (iterative self-organization analysis algorithm). Second, the attribute of the unsupervised classification results was judged according to the high-resolution remote sensing images (from Google Earth) and NLCD. As a result, eight LCTs (i.e., evergreen forest, mixed forest, open forest, shrubbery, grassland, farmland, water body, and construction land) were obtained. LCT is a categorical variable, which is unsuitable for scale conversion. Therefore, LCT data were converted to quantitative variables, that is, the area percentage of each LCT in each 57-m resolution pixel; these variables were recorded as $P_{lct}$ , ( $lct = 1,2, \dots, 8$ ).

2.3.7.

Radiance of visible and near infrared bands

Radiance at the entrance pupil of the first to fifth and seventh bands of ETM+ [ $L_{k}$ , ( $k = 1,2, \dots, 5,7)$ ] could also be used as explanatory variables to determine whether adding $L_{k}$ into the aforementioned variables could further improve the sharpening effects. So, at-entrance-pupil radiance layers of these six bands were calculated using the parameters of gain and bias in metadata. The layers’ resolution was converted from 28.5 to 57 m using the aggregation method.

Through the aforementioned steps, thirty-eight 57-m-resolution raster layers (including six $I_{dur_d}$ , eight $I_{DUR_D}$ , PBI, TWI, NDVI′, VC, AWVC, AT, $ϵ$ , PC, $α$ , ELEV, eight $P_{lct}$ , and six $L_{k}$ ) of the explanatory variables were acquired. For every 57-m-resolution raster layer of the explanatory variables, every $18 \times 18 pixel$ set was divided into one group (corresponding to a BT1026 pixel) for pixel value averaging, thereby obtaining thirty-eight 1026-m-resolution raster layers of the explanatory variables.

2.4.

Combining Alternative Explanatory Variables

Considering that testing all combinations of 38 variables was laborious and unnecessary, only specific variable groups with particular purposes were tested (Table 1). Among these groups, group 1 was the default explanatory variable group by GSM, and groups 2 to 9 were used to test the influence of cutting or adding explanatory variables on sharpening.

Table 1

Different combinations of alternative explanatory variables.

Group No.	Alternative explanatory variables	Description
1	$I_{DUR_D} (D = 1,2 \dots 8)$ , $I_{dur_d} (d = 1,2 \dots 6)$ , NDVI′, VC, TWI, PC, ELEV, $α$ , $ϵ$ , PBI, AT, $P_{lct} (lct = 1,2, \dots, 8)$ , AWVC	Basic explanatory variable group
2	NDVI′, VC, TWI, PC, ELEV, $α$ , $ϵ$ , PBI, AT, $P_{lct} (lct = 1,2, \dots 8)$ , AWVC	Subtract solar radiation based on group 1
3	$I_{DUR_D} (D = 1,2 \dots 8)$ , $I_{dur_d} (d = 1,2 \dots 6)$ , TWI, PC, ELEV, $α$ , $ϵ$ , PBI, AT, $P_{lct} (lct = 1,2, \dots 8)$ , AWVC	Subtract NDVI′ and VC based on group 1
4	$I_{DUR_D} (D = 1,2 \dots 8)$ , $I_{dur_d} (d = 1,2 \dots 6)$ , NDVI′, VC, PC, ELEV, $α$ , $ϵ$ , PBI, AT, $P_{lct} (lct = 1,2, \dots 8)$ , AWVC	Subtract TWI based on group 1
5	$I_{DUR_D} (D = 1,2 \dots 8)$ , $I_{dur_d} (d = 1,2 \dots 6)$ , NDVI′, VC, TWI, PC, ELEV, $α$ , $ϵ$ , PBI, AT, AWVC	Subtract $P_{lct}$ based on group 1
6	$I_{DUR_D} (D = 1,2 \dots 8)$ , $I_{dur_d} (d = 1,2 \dots 6)$ , TWI, PC, ELEV, $α$ , PBI, AT, AWVC	Subtract all variables derived from the remote sensing image based on group 1
7	PBI, AT, AWVC, NDVI, $P_{lct} (lct = 1,2, \dots 8)$	Subtract the variables derived from DEM based on group 1
8	$I_{DUR_D} (D = 1,2 \dots 8)$ , $I_{dur_d} (d = 1,2 \dots 6)$ , NDVI′, VC, TWI, PC, ELEV, $α$ , $ϵ$ , $P_{lct} (lct = 1,2, \dots 8)$	Subtract the macro variables based on group 1
9	$I_{DUR D} (D = 1,2 \dots 8)$ , $I_{dur d} (d = 1,2 \dots 6)$ , NDVI′, VC, TWI, PC, ELEV, $α$ , $ϵ$ , PBI, AT, $P_{lct} (lct = 1,2, \dots 8)$ , AWVC, $L_{k} (k = 1,2, 3,4, 5,7)$	Add radiance at the Vis-NIR bands based on group 1

2.5.

Selecting the Sample Pixels

After this process, images of BT and explanatory variables at 1026-m resolution were prepared. All or a part of the pixels in these images could be selected as sample pixels, values of which (i.e., values of BT and explanatory variables) would be used for building the MARS model. The internal homogeneity degrees varied for different sample pixels. Previous research posited that the statistical model built based on samples with high homogeneity degrees is less scale dependent than one built based on samples with low homogeneity degrees.¹⁴ To observe the influence of the homogeneity degree of sample pixels on the model for the mountainous region, two groups of pixel samples were selected. The first group was the homogeneous pixel sample (HoPS) group. To select HoPSs, the variation coefficients of $18 \times 18 pixel$ values at 57-m resolution corresponding to each pixel at 1026-m resolution were calculated and ordered from high to low for each kind of explanatory variable; if the variation coefficients of all explanatory variables of a pixel at 1026-m resolution were listed out of top 30%, then this pixel was regarded as a HoPS. The second group was the heterogeneous pixel sample (HePS) group, in which all pixels at 1026-m resolution were selected as the pixel samples without considering their homogeneity degrees.

2.6.

Multivariate Adaptive Regression Spline

MARS is an automatic regression modeling method that can predict an output variable according to a group of explanatory variables. This method excels at finding the complex data structure that often hides in high-dimensional data. The method divides the entire high-dimensional data space into several regions. In each specific region, an independent function is used for fitting. The core formula of MARS model is

Eq. (8)

y = β_{0} + \sum_{m = 1}^{M} β_{m} B F_{m},

where

{B F}_{m}

is the basis function,

β_{0}

and

β_{m}

are coefficients to be estimated,

m

is the series number of BF and its coefficient

B F_{m}

can be expressed in the following forms:

Eq. (9)

B F_{m} = \max (0, x_{f} - t),

Eq. (10)

B F_{m + 1} = \max (0, t - x_{f}),

where

t

is the knot location and

x_{f}

is the

f

’th explanatory variable. BF is important in specifying the action area of the explanatory variable. In modeling, BFs were first constructed (by selecting the explanatory variable and its knot location) and were then gradually added to the model until the preset maximum number of BFs was achieved (the preset number in this study was 300). This process generated a group of candidate MARS models. The optimal model was selected from the candidate MARS models according to the mean-squared error (MSE). The reserved 20% of samples was adopted as the verification set to calculate MSE, and the model with the minimum MSE was regarded as the optimal MARS model.

In this study, nine groups of the aforementioned explanatory variables (see Table 1) were taken as the dependent variables ( $x_{f}$ ) to build the MARS models one by one. Given that two sampling methods (i.e., HoPS and HePS methods) were used, 18 MARS models were developed.

2.7.

Simulation of Thermal Images at 57-m Resolution and Residual Redistributing

The 57-m-resolution explanatory variable raster layers were input into the 18 developed MARS models using the raster calculation tool in GIS software, and then 18 simulation results of BT at 57-m resolution (i.e., sharpening results, written as BT57′ images) were obtained. BT1026 image was used to modify each BT57′ image through the following procedures. First, a BT57′ image was reaggregated into an image at 1026-m resolution (each set of $18 \times 18$ BT57′ pixels were aggregated into a pixel at 1026-m resolution, and its pixel value was the average of the corresponding $18 \times 18$ BT57′ pixels) written as BT1026′. Second, the BT1026′ image was subtracted from the BT1026 image, resulting in $⊿ BT 1026$ , which was then added to the BT57′ image to obtain the result after modification, written as BT57″; this process is called residual redistributing (RR) below. The BT57′ and BT57″ images are the final sharpening results.

2.8.

Validation of the Sharpening Results

The mean absolute error (MAE) and root mean-square error (RMSE) of the sharpening results (BT57′ and BT57″) were calculated with the BT57 image as a reference to measure accuracy. The MAE and RMSE were calculated as follows:

Eq. (11)

MAE = [\sum_{m = 1}^{M} \sum_{n = 1}^{N} | SR (m, n) - BT 57 (m, n) |] / (M \times N),

Eq. (12)

RMSE = \sqrt{{\sum_{m = 1}^{M} \sum_{n = 1}^{N} {[SR (m, n) - BT 57 (m, n)]}^{2}} / (M \times N)},

where

SR (m, n)

represents the value of the cell at the

m

’th column and

n

’th row in the sharpening results (BT57′ or BT57″),

M

is the amount of columns in BT57, and

N

is the amount of rows in BT57.

3. Results

3.1.

Evaluation of the Sharpening Effect

As previously mentioned, 18 MARS models and 36 sharpening results (18 BT57′ images and 18 BT57″ images) were obtained in this study. Table 2 lists the basic parameters of these MARS models and the MAEs and RMSEs of the corresponding sharpening results. In Table 2, the sharpening results of group 1, group 3, group 4, and group 9 based on the HePS method have the lowest MAEs and RMSEs. These MAEs before and after RR are as low as about 1.3 and 1.2 K respectively, showing the sharpening results are of high accuracy.

Table 2

Basic parameters of the MARS models and the sharpening errors. The data before and after “\” are the errors of results of HePS and HoPS method, respectively.

Group No.	Basic parameter of MARS model			Sharpening error (K)
Group No.	Number of alternative explanatory variables	Number of BFs	Self error of the model	BT57′ image MAE	BT57″ image MAE	BT57′ image RMSE	BT57″ image RMSE
1	32	42\19	0.359\0.179	1.258\2.115	1.134\1.655	1.985\3.768	1.871\2.671
2	18	61\31	0.393\0.210	1.361\3.613	1.227\3.212	2.113\7.508	2.000\5.862
3	30	48\28	0.356\0.162	1.263\1.774	1.112\1.444	1.767\2.806	1.624\2.203
4	31	50\22	0.367\0.188	1.241\2.341	1.106\2.119	1.748\3.464	1.611\2.929
5	24	55\23	0.512\0.280	1.555\2.194	1.407\1.929	2.330\4.580	2.168\4.269
6	21	38\19	2.734\0.719	2.380\3.600	2.330\2.944	3.224\4.851	3.127\3.881
7	12	99\18	0.540\0.275	1.471\2.455	1.310\1.894	2.160\4.379	1.991\3.007
8	29	73\12	0.735\0.550	2.412\2.487	1.756\1.846	3.584\4.503	2.949\3.474
9	38	92\28	0.265\0.158	1.242\1.436	1.166\1.275	1.705\2.481	1.621\2.280

To measure the sharpening effects of the proposed method, we compared them with the effects of two simple processing methods (bilinear resampling and cubic convolution resampling) that can be regarded as benchmarks for measuring the effect of any thermal sharpening method. The results show that the MAE of both simple processing methods is 1.776 K. The MAE of the best results obtained in this study is smaller than that of the simple processing methods by $\sim 0.6 K$ .

We qualitatively evaluate the sharpening effects by visual inspection. The visual presentation of the various sharpening results is closely correlated to their error. The smaller the error is, the better the visual presentation. To save space, we only present part of the sharpening results of the models constructed based on the first group of explanatory variables (Fig. 1). Comparing Fig. 1(e) with Fig. 1(b) reveals that the sharpening result based on the default explanatory variable group, and the HePS method shows significant improvement in its ability to reflect details of spatial variation of thermal information. Comparing Fig. 1(e) with Fig. 1(c) reveals that this sharpening result is very similar to the reference image. Figure 1(e) is much better than the result of cubic convolution resampling [see Fig. 1(d)]. Figures 1(g) and 1(h) are less similar than Figs. 1(e) or 1(f) to the reference image, indicating that the HoPS method is worse than the HePS method in a mountainous region. Very slight artificial box-like traces are found in specific areas in Fig. 1(f), and rather serious artificial box-like traces are found in Fig. 1(h).

Fig. 1

(a) True color composite, (b) BT image to be sharpened (BT1026), (c) reference BT image (BT57), (d) sharpening result based on cubic convolution resampling, (e) sharpening result based on first variable group and HePS and before RR, (f) sharpening result based on first variable group and HePS and after RR, (g) sharpening result based on first variable group and HoPS and before RR, and (h) sharpening result based on first variable group and HoPS and after RR.

3.2.

Effect of Sampling on Sharpening

Table 2 shows that in terms of the self errors, the errors of the models built by the HoPS method are smaller than those of the models built by the HePS method. However, in terms of the errors of the sharpening results (MAE and RMSE), the errors of the model built by the HoPS method are larger than those of the models built by the HePS method (before and after RR). The reason behind this condition was the fact that the HoPS method pursued the consistency of $18 \times 18$ 57-m resolution pixels covered by each 1026-m-resolution sample pixel on explanatory variables, and thereby the majority of the extracted sample pixels belonged to the flat region with weak regional representativeness of the large area of rugged terrain. The values of the HoPSs were more concentrated than the overall value distribution. The MARS models built by the HoPS method had a simple model structure (i.e., small number of BFs) and small self error. Nonetheless, given the obvious difference between the value distribution of the HoPSs and the overall value distribution, large errors generated when the HoPS-based models were applied to the entire research region. Thus, the HoPS method cannot adapt to the thermal image sharpening of the mountainous region. The analysis in the subsequent paragraph analyzes only the HePS-based models without considering HoPS-based models.

3.3.

Effect of Explanatory Variables on Sharpening

Table 2 indicates that the errors of the sharpening effect obtained by different explanatory variable groups significantly varied. The variables in group 1 were all alternative ones [except $L_{k} (k = 1,2, 3,4, 5,7)$ ], which generated a sharpening result with very low MAE and RSME. Groups 2 to 8 eliminated one or several variables based on group 1. The MAE and RSME change after specific variable(s) was eliminated that could reflect the importance of the eliminated variable(s) to sharpening. The elimination of solar radiation increased MAE by $\sim 0.1 K$ (see group 2), indicating that this variable is valuable for thermal sharpening. Eliminating NDVI′ and VC only slightly influenced MAE (see group 3). The reason behind the case was the condition that another explanatory variable $ϵ$ was highly correlated with NDVI′; $ϵ$ replaced the explanatory effect of NDVI′ and VC. After the elimination of TWI, MAE did not rise but reduced slightly (see group 4). This finding may be due to the fact that TWI inaccurately reflected the spatial differentiation of soil moisture decided by the terrain and introduced numerous new noises. After $P_{lct}$ was eliminated, MAE rose by $\sim 0.3 K$ (see group 5), implying that these variables are very valuable for thermal sharpening. If all variables derived from images at Vis-NIR bands were eliminated [i.e., NDVI′, VC, $ϵ$ , and $P_{lct}$ (see group 6)], then the MAE of the sharpening result would reach over 2 K. This phenomenon suggests that thermal sharpening could not be separated from these variables completely. If all variables reflecting topographic features (see group 7) were eliminated, then MAE would increase by nearly 0.2 K. If three macrovariables (i.e., PBI, AT, and AWVC) were eliminated simultaneously, then the MAEs of BT57′ and BT57″ images would increase by $\sim 1.2$ and 0.6 K (see group 8), respectively. The common characteristic of the macrovariables is that they have the spatial variation only on the macroscale and that they are nearly uniform inside each 1026 pixel. As such, the BT variation inside the $1026 \times 1026 m$ area should bear no relation to these macrovariables. However, the results of group 8 show the importance of macrovariables to sharpening. The reason behind this event is the case that the loss of macrovariables sharply increased the self errors of the MARS model (see Table 2), which in turn caused large errors in the sharpening result. Group 9 added the Vis-NIR radiance variables based on the variables in group 1, and the MAE of the sharpening results reduced slightly. Nevertheless, group 9 did not get rid of the dependence on synchronously acquired Landsat Vis-NIR images. Vis-NIR radiance changes rapidly with date; hence, the Vis-NIR radiance of the target date cannot be accurately calculated with Vis-NIR images on other dates. This makes such explanatory variable group less feasible in practical applications.

The explanatory variable group built based on the prior knowledge about the energy budge and thermal infrared electromagnetic radiation transmission (i.e., group 1 in this study) could accurately conduct sharpening on thermal image. In practice, choosing the explanatory variables is difficult, especially in regions with complex environmental conditions. This study suggests that in this case, the explanatory variable group should be determined based on the laws of land surface energy budget and thermal infrared electromagnetic radiation transmission, and thereby the blind choice of explanatory variables can be avoided.

4. Discussion

4.1.

Comparison with Previous Methods

In previous studies of thermal imagery sharpening based on statistical models, the NDVI or its transformed form was most commonly used as a scale factor.²^,³^,⁵^,²¹^,²⁷ Reference 27 performed sharpening of multiple resolutions using the NDVI as the scale factor. As in the current study, that study downscaled the resolution from 1536 to 192 m, and the RMSE of the results was 1.61 K. Reference 5 used four forms of NDVI-based function to sharpen the resolution from 960 to 60 m, and the RMSEs ranged from 1.20 to 3.27 K. Reference 10 used a global model, a resolution-adjusted global model, a piecewise regression model, a stratified model, and a local model to sharpen the resolution from 990 to 90 m, and the RMSEs ranged from 3.24 to 3.68 K. Using an artificial neural network method and 11 explanatory variables including the NDVI, difference vegetation index, and leaf area index, Ref. 8 sharpens LST images; the RMSEs of the results sharpened from 1080 to 90 m were 2.443 K (vegetated area) and 1.573 K (bare area). Using a regression tree method and six Vis-NIR bands as the explanatory variables, Ref. 14 sharpens BT images; the MAEs of the sharpening results in flat areas were generally smaller than 1 K, while large errors were obtained in mountainous regions (960 m to 120 m, MAE 1.67 K).

Previous research has shown that within mountainous regions, it is more difficult to sharpen the thermal imagery.¹⁴^,¹⁵ In the current study, we performed sharpening of the thermal imagery in a typical mountainous region. The errors in the current results can be controlled at a relatively low level. For the model built with the default explanatory variable group of GSM and based on HePS method, the MAEs were 1.258 K (before RR) and 1.134 K (after RR) and the RMSE were 1.985 K (before RR) and 1.871 K (after RR). Thus, we tentatively put forward that the GSM-based method is superior to previous methods for mountainous regions.

For a more accurate comparison with TsHARP, we performed sharpening using five TsHARP-based methods proposed by previous study (an NDVI-based linear function model,⁵^,¹⁰^,²⁷ an NDVI-based quadratic function model,⁵ a vegetation cover model,⁵ a simplified vegetation cover model,⁵ and a local NDVI-based model¹⁰), and the thermal image to be sharpened was BT1026. The MAEs of the sharpening results were 1.631, 1.627, 1.668, 1.668, and 1.591 K, and the RMSE values were 2.223, 2.425, 2.226, 2.226, and 2.214 K, respectively. These MAEs and RMSEs were larger than corresponding MAEs and RMSEs of the sharpening results in the GSM-based method.

4.2.

Localization of the Model

Statistical sharpening models can be divided into two kinds, namely, global and local models. The global model refers to a single model built by the sample pixels from the entire study area. Thermal image sharpening with a global model assumes that the target variable–explanatory variables relation is uniform across the entire area, but this assumption often contradicts the reality. In the real world, the target variable–explanatory variables relation cannot remain the same globally. If samples are collected in each local region in the form of a sliding window and the special model for each local region is established based on these samples, then the model is a local model. Several studies have shown that for thermal image sharpening, a local model performs better than a global one.¹⁴^,¹⁸

In this study, the statistical models between BT and explanatory variables were established through MARS instead of with the method based on sliding window. In essence, the MARS models in this study are local models. To explain this problem, the model built with the variables in group 1 was chosen as an example. The model is depicted as follows:

Eq. (13)

BF 1 = \max (0, ϵ - 0.986), BF 2 = \max (0, 0.986 - ϵ), BF 3 = \max (0, P_{8} - 0.0941), BF 5 = \max (0, P B I - 12.586), BF 6 = \max (0, 12.586 - PBI), \dots, BF 285 = \max (0, I_{DUR_2} - 153.168),

Eq. (14)

BT = 292.281 - 842.411 \times BF 1 + 838.992 \times BF 2 - 1.527 \times BF 3 - 0.0324524 \times BF 5 + 0.190 \times BF 6 \dots, \dots - 0.0309 \times BF 285,

where a total of 42 BFs were obtained but 36 BFs are omitted here due to space limitations. In the formulas, MARS divides the entire high-dimensional data space into several subspaces, each corresponding to a special form of the model. In other words, for each special value combination of those explanatory variables, only a part of the BFs are active and the rest are constant zero. When each data subspace is mapped to the geographic space, it becomes large or small local patches with similar environmental properties. In each kind of patch, the model has a specific form. Therefore, MARS can reflect the local relationship between target variable and explanatory variables and can adapt to the thermal image sharpening of the complex region.

4.3.

Feasibility of Geographical Statistical Model in Practice

The presented results indicate that the GSM method is capable of generating a fine spatial resolution thermal image from a coarse spatial resolution image ( $\sim 1000 - m$ resolution, e.g., the thermal images of MODIS and AVHRR) accurately. However, the following issue needs further discussion. In this study, the NDVI′, $ϵ$ , and $P_{lct}$ raster data at 57-m resolution used by GSM were generated based on the ETM+ Vis-NIR data. How can these data be obtained in practical applications? For special regions, most satellite-borne moderate-resolution sensors (e.g., ETM+, Operational Land Imager, and Advanced Spaceborne Thermal Emission and Reflection Radiometer) cannot produce Vis-NIR data every day. Thus, when conducting thermal image sharpening on the $J$ ’th day, we may need the Vis-NIR data of moderate resolution sensors on the $(j + n)$ ’th or $(j - n)$ ’th day to generate NDVI′, $ϵ$ , and $P_{lct}$ . Generally, NDVI′ and $ϵ$ are stable within several days, and LCT is stable in a year. Thus, it is feasible to using NDVI′, $ϵ$ , and LCT of the $(j + n)$ ’th or $(j - n)$ ’th day ( $n < 10$ ) to serve as the corresponding variables on the $J$ ’th day. Other variables in the default explanatory variable group of GSM are all stable in the time domain (e.g., slope and elevation), or the data of these variables can be observed or simulated every day (e.g., precipitation, AT, and solar radiation). Therefore, it can be concluded that the GSM sharpening method can eliminate the dependence on synchronous imagery at VIS-NIR bands.

It should be noted that the GSM-based method needs raster layers of many variables involved in land surface energy budget and thermal infrared electromagnetic radiation transmission. The establishment of these layers requires sufficient initial data (DEM, meteorological data, and so on) and much calculation time. In terms of data availability and implementation simplicity, several methods proposed by previous studies (e.g., Refs. 3, 5, 14, and so on) are better than the GSM-based method.

The robustness of this method needs further investigations; specifically, future work should address the variation of the importance of explanatory variables in various regions. In addition, some satellite remote-sensing data products (like MODIS water vapor, MODIS atmosphere profile, and TRMM precipitation products) possess well-defined physical meaning; they are likely suitable for acting as explanatory variables and can enhance the thermal sharpening effect. Once these data are used as explanatory variables, they should be deliberately analyzed in respect of spatial resolution, data quality, and the degree to which they synchronize temporally with the sharpened objects.

5. Conclusions

A sharpening method of thermal image based on GSM was proposed in this study. This method attached more importance to geography law and obtained sharpening results with fewer errors and good visual effects in mountainous areas with complex environment. The GSM-based sharpening method can avoid the blind choice of explanatory variables and reduce the dependence on synchronous imagery at VIS-NIR bands, which makes it more feasible in practice. The method is helpful to the agricultural and forestry studies and practices (e.g., drought monitoring and wildfire danger modeling) in mountainous regions using coarse spatial resolution thermal images (e.g., the images of MODIS and AVHRR). Furthermore, this study found that: (1) the homogeneous sample pixel model is not suitable for sharpening of a thermal image in a mountain area; (2) although macrovariables, such as PBI, AT, and AWVC, can only present spatial variability in macroscale, they are indispensable for thermal image sharpening; and (3) the MARS modeling method is capable of reflecting the locality of the relationship between the target variable and explanatory variables, thereby reducing errors of sharpening results efficiently.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (No. 41201099) and Postdoctoral Research Funds of Henan University (Grant No. BH2012042). We are thankful to the United States Geological Survey (USGS) for providing the ETM+, the National Elevation Dataset, and the National Land Cover Database. We are also thankful to the National Climatic Data Center of National Oceanic and Atmospheric Administration of US for providing meteorological data. We also deeply appreciate the American National Renewable Energy Laboratory for providing the solar radiation data and the atmospheric data related with solar radiation. We also thank the two anonymous reviewers for their insightful comments, which improved this paper.

References

1.

W. F. Zhan et al., “Disaggregation of remotely sensed land surface temperature: literature survey, taxonomy, issues, and caveats,” Remote Sens. Environ., 131 119 –139 (2013). http://dx.doi.org/10.1016/j.rse.2012.12.014 Google Scholar

2.

S. Mukherjee, P. K. Joshi and R. D. Garg, “A comparison of different regression models for downscaling Landsat and MODIS land surface temperature images over heterogeneous landscape,” Adv. Space Res., 54 (4), 655 –669 (2014). http://dx.doi.org/10.1016/j.asr.2014.04.013 ASRSDW 0273-1177 Google Scholar

3.

S. Mukherjee, P. K. Joshi and R. D. Garg, “Regression-Kriging technique to downscale satellite-derived land surface temperature in heterogeneous agricultural landscape,” IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens., 8 (3), 1245 –1250 (2015). http://dx.doi.org/10.1109/JSTARS.2015.2396032 Google Scholar

4.

P. Sismanidis et al., “Evaluating the operational retrieval and downscaling of urban land surface temperatures,” IEEE Geosci. Remote Sens. Lett., 12 (6), 1312 –1316 (2015). http://dx.doi.org/10.1109/LGRS.2015.2397450 Google Scholar

5.

N. Agam et al., “A vegetation index based technique for spatial sharpening of thermal imagery,” Remote Sens. Environ., 107 (4), 545 –558 (2007). http://dx.doi.org/10.1016/j.rse.2006.10.006 Google Scholar

6.

D. Fasbender et al., “Support-based implementation of Bayesian data fusion for spatial enhancement: Applications to ASTER thermal images,” IEEE Geosci. Remote Sens. Lett., 5 (4), 598 –602 (2008). http://dx.doi.org/10.1109/LGRS.2008.2000739 Google Scholar

7.

M. Stathopoulou and C. Cartalis, “Downscaling AVHRR land surface temperatures for improved surface urban heat island intensity estimation,” Remote Sens. Environ., 113 (12), 2592 –2605 (2009). http://dx.doi.org/10.1016/j.rse.2009.07.017 Google Scholar

8.

G. J. Yang et al., “A novel method to estimate subpixel temperature by fusing solar-reflective and thermal-infrared remote-sensing data with an artificial neural network,” IEEE Trans. Geosci. Remote Sens., 48 (4), 2170 –2178 (2010). http://dx.doi.org/10.1109/TGRS.2009.2033180 IGRSD2 0196-2892 Google Scholar

9.

W. F. Zhan et al., “Sharpening thermal imageries: a generalized theoretical framework from an assimilation perspective,” IEEE Trans. Geosci. Remote Sens., 49 (2), 773 –789 (2010). http://dx.doi.org/10.1109/TGRS.2010.2060342 IGRSD2 0196-2892 Google Scholar

10.

C. Jeganathan et al., “Evaluating a thermal image sharpening model over a mixed agricultural landscape in India,” Int. J. Appl. Earth Obs. Geoinf., 13 (2), 178 –191 (2011). http://dx.doi.org/10.1016/j.jag.2010.11.001 Google Scholar

11.

A. Dominguez et al., “High-resolution urban thermal sharpener (HUTS),” Remote Sens. Environ., 115 (7), 1772 –1780 (2011). http://dx.doi.org/10.1016/j.rse.2011.03.008 Google Scholar

12.

M. C. Anderson et al., “Mapping daily evapotranspiration at field to continental scales using geostationary and polar orbiting satellite imagery,” Hydrol. Earth Syst. Sci., 15 223 –239 (2011). http://dx.doi.org/10.5194/hess-15-223-2011 Google Scholar

13.

V. Rodriguez-Galiano et al., “Downscaling Landsat 7 ETM+ thermal imagery using land surface temperature and NDVI image,” Int. J. Appl. Earth Obs. Geoinf., 18 515 –527 (2012). http://dx.doi.org/10.1016/j.jag.2011.10.002 Google Scholar

14.

F. Gao, W. P. Kustas and M. C. Anderson, “A data mining approach for sharpening thermal satellite imagery over land,” Remote Sens., 4 (11), 3287 –3319 (2012). http://dx.doi.org/10.3390/rs4113287 Google Scholar

15.

J. L. Quan et al., “Downscaling remotely sensed land surface temperatures: a comparison of typical methods,” J. Remote Sens., 17 (2), 363 –375 (2013). Google Scholar

16.

J. Nichol, “An emissivity modulation method for spatial enhancement of thermal satellite images in urban heat island analysis,” Photogramm. Eng. Remote Sens., 75 (5), 547 –556 (2009). http://dx.doi.org/10.14358/PERS.75.5.547 Google Scholar

17.

W. Essa et al., “Evaluation of the DisTrad thermal sharpening methodology for urban areas,” Int. J. Appl. Earth Obs., 19 163 –172 (2012). http://dx.doi.org/10.1016/j.jag.2012.05.010 Google Scholar

18.

K. Zakšek and K. Oštir, “Downscaling land surface temperature for urban heat island diurnal cycle analysis,” Remote Sens. Environ., 117 (15), 114 –124 (2012). http://dx.doi.org/10.1016/j.rse.2011.05.027 Google Scholar

19.

F. Wang et al., “An efficient approach for pixel decomposition to increase the spatial resolution of land surface temperature images from MODIS thermal infrared band data,” Sensors, 15 (1), 304 –330 (2015). http://dx.doi.org/10.3390/s150100304 SNSRES 0746-9462 Google Scholar

20.

D. S. Liu and X. L. Zhu, “An enhanced physical method for downscaling thermal infrared radiance,” IEEE Geosci. Remote Sens. Lett., 9 (4), 690 –694 (2012). http://dx.doi.org/10.1109/LGRS.2011.2178814 Google Scholar

21.

X. H. Chen et al., “A combination of TsHARP and thin plate spline interpolation for spatial sharpening of thermal imagery,” Remote Sens., 6 (6), 2845 –2863 (2014). http://dx.doi.org/10.3390/rs6042845 Google Scholar

22.

J. Y. Fang, Z. H. Shen and H. T. Cui, “Ecological characteristics of mountains and research issues of mountain ecology,” Biodivers. Sci., 12 (1), 10 –19 (2004). http://dx.doi.org/10.17520/biods.2004003 Google Scholar

23.

R. G. Barry, Mountain Weather and Climate, 3rd ed.Cambridge University Press, Cambridge (2003). Google Scholar

24.

M. Beniston, “Mountain weather and climate: a general overview and a focus on climatic change in the Alps,” Hydrobiologia, 562 (1), 3 –16 (2006). http://dx.doi.org/10.1007/s10750-005-1802-0 HYDRB8 0018-8158 Google Scholar

25.

M. Iqbal, An Introduction to Solar Radiation, Academic Press, Don Mills, Ontario (1983). Google Scholar

26.

J. Dozier and J. Frew, “Rapid calculation of terrain parameters for radiation modeling from digital elevation data,” IEEE Trans. Geosci. Remote Sens., 28 (5), 963 –969 (1990). http://dx.doi.org/10.1109/36.58986 IGRSD2 0196-2892 Google Scholar

27.

W. P. Kustas et al., “Estimating subpixel surface temperatures and energy fluxes from the vegetation index-radiometric temperature relationship,” Remote Sens. Environ., 85 (4), 429 –440 (2003). http://dx.doi.org/10.1016/S0034-4257(03)00036-1 Google Scholar

28.

N. Agam et al., “Utility of thermal image sharpening for monitoring field-scale evapotranspiration over rainfed and irrigated agricultural regions,” Geophys. Res. Lett., 35 1 –7 (2008). http://dx.doi.org/10.1029/2007GL032195 GPRLAJ 0094-8276 Google Scholar

29.

T. R. Oke, Boundary Layer Climates, 2nd ed.Methuen, New York (1987). Google Scholar

30.

S. M. Bateni and D. Entekhabi, “Relative efficiency of land surface energy balance components,” Water Resour. Res., 48 (4), W04510 (2012). http://dx.doi.org/10.1029/2011WR011357 WRERAQ 0043-1397 Google Scholar

31.

J. H. Friedman, “Multivariate adaptive regression splines,” Ann. Stat., 19 (1), 1 –67 (1991). http://dx.doi.org/10.1214/aos/1176347963 ASTSC7 0090-5364 Google Scholar

32.

Z. H. Qin, A. Karnieli and P. Berliner, “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., 22 (18), 3719 –3746 (2001). http://dx.doi.org/10.1080/01431160010006971 IJSEDK 0143-1161 Google Scholar

33.

J. C. Jiménez-Muñoz and J. A. Sobrino, “A generalized single-channel method for retrieving land surface temperature from remote sensing data,” J. Geophys. Res., 108 (D22), 2015 –2023 (2003). http://dx.doi.org/10.1029/2003JD003480 JGREA2 0148-0227 Google Scholar

34.

J. A. Sobrino, J. El Kharraz and Z. L. Li, “Surface temperature and water vapour retrieval from MODIS data,” Int. J. Remote Sens., 24 (24), 5161 –5182 (2003). http://dx.doi.org/10.1080/0143116031000102502 IJSEDK 0143-1161 Google Scholar

35.

K. B. Mao, Z. H. Qin and J. C. Shi, “A practical split-window algorithm for retrieving land-surface temperature from MODIS data,” Int. J. Remote Sens., 26 (15), 3181 –3204 (2005). http://dx.doi.org/10.1080/01431160500044713 IJSEDK 0143-1161 Google Scholar

36.

J. C. Jiménez-Muñoz et al., “Land surface temperature retrieval methods from Landsat-8 thermal infrared sensor data,” IEEE Geosci. Remote Sens., 11 (10), 1840 –1843 (2014). http://dx.doi.org/10.1109/LGRS.2014.2312032 Google Scholar

37.

J. A. Sobrino, J. C. Jiménez-Muñoz and L. Paolini, “Land surface temperature retrieval from Landsat TM5,” Remote Sens. Environ., 90 (4), 434 –440 (2004). http://dx.doi.org/10.1016/j.rse.2004.02.003 Google Scholar

38.

X. L. Yu, X. L. Guo and Z. C. Wu, “Land surface temperature retrieval from Landsat 8 TIRS—comparison between radiative transfer equation-based method, split window algorithm and single channel method,” Remote Sens., 6 (10), 9829 –9852 (2014). http://dx.doi.org/10.3390/rs6109829 Google Scholar

39.

R. B. Stull, An Introduction to Boundary Layer Meteorology, Kluwer Academic Publishers, Dordrecht, The Netherlands (1988). Google Scholar

40.

W. G. M. Bastiaanssen, “SEBAL-based sensible and latent heat fluxes in the irrigated Gediz Basin, Turkey,” J. Hydrol., 229 (1–2), 87 –100 (2000). JHYDA 70022-1694 Google Scholar

41.

G. Bertoldi et al., “Topographical and ecohydrological controls on land surface temperature in an alpine catchment,” Ecohydrology, 3 (2), 189 –204 (2010). http://dx.doi.org/10.1002/eco.129 Google Scholar

42.

C. Ricotta and G. C. Avena, “Relation between vegetation canopy surface temperature and the Sun-surface geometry in a mountainous region of central Italy,” Int. J. Remote Sens., 18 (14), 3091 –3096 (1997). http://dx.doi.org/10.1080/014311697217251 IJSEDK 0143-1161 Google Scholar

43.

X. Li et al., “Modification of solar radiation model over rugged terrain,” Chin. Sci. Bull., 44 (15), 1345 –1349 (1999). http://dx.doi.org/10.1007/BF02885977 CSBUEF 1001-6538 Google Scholar

44.

K. Yang, G. W. Huang and N. Tamai, “A hybrid model for estimating global solar radiation,” Sol. Energy, 70 (1), 13 –22 (2001). http://dx.doi.org/10.1016/S0038-092X(00)00121-3 Google Scholar

45.

M. P. González-dugo et al., “Canopy temperature variability as an indicator of crop water stress severity,” Irrig. Sci., 24 (4), 233 –240 (2005). http://dx.doi.org/10.1007/s00271-005-0022-8 IRSCD2 1432-1319 Google Scholar

46.

G. S. Buttar et al., “Canopy temperature: a method to estimate plant water stress and scheduling irrigation in cotton and wheat,” J. Agr. Phys., 5 (1), 79 –83 (2005). Google Scholar

47.

ArcGIS resource center, “How to direction works,” (2011) http://help.arcgis.com/en/arcgisdesktop/10.0/help/index.html#//009z00000063000000.htm Google Scholar

48.

P. M. Teillet, B. Guindon and D. G. Goodenough, “On the slope-aspect correction of multispectral scanner data,” Can. J. Remote Sens., 8 (2), 84 –106 (1982). http://dx.doi.org/10.1080/07038992.1982.10855028 CJRSDP 0703-8992 Google Scholar

49.

G. L. Zhu et al., “Evaluation of topographic effects on four commonly used vegetation indices,” J. Remote Sens., 17 (1), 210 –221 (2013). Google Scholar

50.

Z. H. Qin et al., “Estimating of the essential atmospheric parameters of mono-window algorithm for land surface temperature retrieval from Landsat TM6,” Remote Sens. Land Resour., 15 (2), 37 –43 (2003). http://dx.doi.org/10.6046/gtzyyg.2003.02.10 Google Scholar

51.

X. A. Liu et al., “Study on spatialization technology of terrestrial eco-information in China (III): temperature and precipitation,” J. Nat. Resour., 19 (6), 818 –825 (2004). http://dx.doi.org/10.11849/zrzyxb.2004.06.018 JRLEEJ 1059-9053 Google Scholar

52.

X. Li, G. D. Cheng and L. Lu, “Spatial analysis of air temperature in the Qinghai-Tibet plateau,” Arct. Antarct. Alp. Res., 37 (2), 246 –252 (2005). http://dx.doi.org/10.1657/1523-0430(2005)037[0246:SAOATI]2.0.CO;2 Google Scholar

53.

L. W. Zevenbergen and C. R. Thorne, “Quantitative analysis of land surface topography,” Earth Surf. Proc. Land., 12 (1), 47 –56 (1987). http://dx.doi.org/10.1002/(ISSN)1096-9837 ESPLDB 0197-9337 Google Scholar

54.

C. Z. Qin et al., “Mapping soil organic matter in small low-relief catchments using fuzzy slope position information,” Geoderma, 171–172 64 –74 (2012). http://dx.doi.org/10.1016/j.geoderma.2011.06.006 GEDMAB 0016-7061 Google Scholar

55.

Sr. R. A. Pielke, “Land use and climate change,” Science, 310 (5754), 1625 –1626 (2005). http://dx.doi.org/10.1126/science.1120529 SCIEAS 0036-8075 Google Scholar

56.

Sr. R. A. Pielke et al., “An overview of regional land-use and land-cover impacts on rainfall,” Tellus B, 59 (3), 587 –601 (2007). http://dx.doi.org/10.1111/teb.2007.59.issue-3 TELLAL 0040-2826 Google Scholar

57.

J. Zhai et al., “Human-induced land cover changes drive a diminution of land surface albedo in the Loess plateau (China),” Remote Sens., 7 (3), 2926 –2941 (2015). http://dx.doi.org/10.3390/rs70302926 Google Scholar

58.

H. L. Xiao and Q. H. Weng, “The impact of land use and land cover changes on land surface temperature in a karst area of China,” J. Environ. Manage., 85 (1), 245 –257 (2007). http://dx.doi.org/10.1016/j.jenvman.2006.07.016 Google Scholar

59.

J. Pongratz et al., “The impact of land cover change on surface energy and water balance in Mato Grosso, Brazil,” Earth Interact., 10 (1), 1 –17 (2006). http://dx.doi.org/10.1175/EI176.1 Google Scholar

Biography

Pengcheng Qi received his PhD degree in physical geography from Lanzhou University, Lanzhou, China, in 2009. Since 2010, he has been working at Nanyang Normal University, where he is currently an associate professor with the Laboratory of Remote Sensing Monitoring of Natural Disaster. His research interests include remote sensing digital image processing and applications of thermal infrared remote sensing in ecology.

Shixiong Hu received his PhD degree in physical geography and environmental analysis from State University of New York at Buffalo, United States, in 2004. Since 2004, he has been working in the Department of Geography, East Stroudsburg University of Pennsylvania. His research interests include hillslope processes, environmental modeling with GIS, and soil and water conservation.

Haijun Zhang received the BS degree in geographic information systems and the ME degree in cartography and geographic information engineering from Chang’an University, China. He is currently a lecturer in the School of Environmental Science and Tourism, Nanyang Normal University, China. His research interests focus on the integration of remote sensing and GIS in addressing environmental issues.

Guangmeng Guo received his PhD from Institute of Geographic Sciences and Natural Resources Research, Chinese Academy of Sciences, in 2004. He is currently a professor with the Laboratory of Remote Sensing Monitoring of Natural Disaster, Nanyang Normal University. His research interests focus on applications of remote sensing in disaster monitoring.

CC BY: © The Authors. Published by SPIE under a Creative Commons Attribution 4.0 Unported License. Distribution or reproduction of this work in whole or in part requires full attribution of the original publication, including its DOI.

Citation Download Citation

Pengcheng Qi, Shixiong Hu, Haijun Zhang, and Guangmeng Guo "Sharpening method of satellite thermal image based on the geographical statistical model," Journal of Applied Remote Sensing 10(2), 025013 (12 May 2016). https://doi.org/10.1117/1.JRS.10.025013

Published: 12 May 2016

Access the abstract

JOURNAL ARTICLE
17 PAGES

DOWNLOAD PAPER SAVE TO MY LIBRARY

GET CITATION

CITATIONS

Cited by 5 scholarly publications.

Explore citations on Lens.org

KEYWORDS

Thermography

Thermal modeling

Image resolution

Statistical analysis

Satellites

Earth observing sensors

Raster graphics

1.

Introduction

2.

Method

2.1.

Study Area and Data Source

2.2.

Preparation of Images for Sharpening

Eq. (1)

2.3.

Selection of Explanatory Variables and Establishment of Their Raster Layers

2.3.1.

Total solar radiation

Eq. (2)

Eq. (3)

Eq. (4)

2.3.2.

Precipitation before imaging and topographic wetness index

Eq. (5)

2.3.3.

Vegetation index and vegetation coverage

Eq. (6)

Eq. (7)

2.3.4.

Emissivity, water vapor content of air column at imaging, and air temperature at imaging

2.3.5.

Profile curvature, slope degree, and elevation

2.3.6.

Land cover type

2.3.7.

Radiance of visible and near infrared bands

2.4.

Combining Alternative Explanatory Variables

Table 1

2.5.

Selecting the Sample Pixels

2.6.

Multivariate Adaptive Regression Spline

Eq. (8)

Eq. (9)

Eq. (10)

2.7.

Simulation of Thermal Images at 57-m Resolution and Residual Redistributing

2.8.

Validation of the Sharpening Results

Eq. (11)

Eq. (12)

3.

Results

3.1.

Evaluation of the Sharpening Effect

Table 2

Fig. 1

3.2.

Effect of Sampling on Sharpening

3.3.

Effect of Explanatory Variables on Sharpening

4.

Discussion

4.1.

Comparison with Previous Methods

4.2.

Localization of the Model

Eq. (13)

Eq. (14)

4.3.

Feasibility of Geographical Statistical Model in Practice

5.

Conclusions

Acknowledgments

References

Biography

Show All Keywords

Keywords/Phrases

Search In:

Publication Years