A Physically Based Spatial Expansion Algorithm for Surface Air Temperature and Humidity

An algorithm was developed to expand the surface air temperature and air humidity to a larger spatial domain, based on the fact that the variation of surface air temperature and air humidity is controlled jointly by the local turbulence and the horizontal advection. This study proposed an algorithm which considers the advective driving force outside the thermal balance system and the turbulent driving force and radiant driving force inside the thermal balance system.The surface air temperature is determined by a combination of the surface observations and the regional land surface temperature observed from a satellite. The average absolute difference of the algorithm is 0.65 degree and 0.31mb, respectively, for surface air temperature and humidity expansion, which provides a promising approach to downscale the two surface meteorological variables.


Introduction
Air temperature and humidity are the most fundamental elements that human beings interact with in the environment.The heat and steam (moisture) that reach the surface (the surface or active surface of soil, vegetation, rocks, and water, generally called the earth's surface) have complex interaction with the air in the boundary layer, causing new balance and redistribution of heat and steam.The most important and frontier indicator for this new temporal and spatial distribution is the air temperature and humidity at the height of thermometer shelter in weather stations.For a long time, air temperature and humidity are not only the primary items of weather forecast, but also the core input to model the surface sensible heat flux and latent heat flux.The Penman-Monteith equation [1] requires the air temperature and humidity and some other inputs to calculate the evapotranspiration on the basis of energy balance.
The spatial distribution of air temperature and humidity depends on the uniformity of the surface energy balance and the intensity of the horizontal advection [2].So far, the observation means for such important parameters is still limited in very small "point" scale [3].The spatial representativeness of the air temperature and humidity at the height of thermometer shelter in weather stations is about a few hundred square meters [4,5].The density of weather stations varies with each country's capacity on meteorological observations and the corresponding financial budget.Currently in China, there is approximately one weather station for each administrative county, which means one weather station for average area of 5846 square kilometers (close to 6000 MODIS pixel; counties in the east are smaller so the density is higher than the average and counties in the west are larger so the density is lower than the average).Apparently, the air temperature and humidity currently observed by weather stations have often a lack of representativeness for the areas they are located in.As a result, the spatial distribution of the regional or global energy and material flows inversed based on these parameters much deviated from the actual condition.
In the assimilation of land surface model and quantitative remote sensing model, there is an urgent demand for air temperature and humidity data with finer spatial resolution, which could not be achieved by the existing weather stations [6].The high temporal resolution of Fengyun meteorological satellite and the high spatial resolution surface temperature from MODIS, combining with the high time resolution air temperature and air humidity observed from weather stations, were used to obtain air temperature and humidity of high time frequency at 1 km spatial resolution in the study.This has critical effect on the assimilation of land surface model and two-layer remote sensing evapotranspiration model [7] and also benefits the imbalance in two-layer evapotranspiration model brought by the energy imbalance and the horizontal advection.
Many methods, linear and nonlinear, of spatial interpolation and spatial extension have been proposed [8][9][10][11]; however, most of them lack mechanism support for surface interaction, causing the interpolated or extended air temperature and humidity to deviate from the actual situation a lot.In linear or nonlinear distance inverse methods (inverse distance interpolation methods) [8][9][10][11], interpolation is based on two weather stations' observation of temperature and humidity and on the distance between the two stations, where the surface uneven thermal effect between two stations is not considered [12].For example, the vegetation cover of two stations is grass, and the surface between them is dry sand; obviously, the air temperature and humidity over the dry sand would be significantly different with the values of this kind of interpolation.
The spatial expansion methods for surface air temperature and humidity can be summarized into four kinds: (i) linear interpolation methods without physical meaning, (ii) nonlinear interpolation methods without physical meaning, (iii) static feedback methods with physical meaning, and (iv) dynamic feedback methods with physical meaning [13].The spatial extension method for air temperature driven by the thermal radiation of surface temperature [14] expressed the radiant thermal effect of surface thermal driving force on the air temperature; it is a static feedback method with physical meaning.However, this method did not consider the turbulent driving force of air temperature and the advective driving force outside the thermal balance system.This study proposed an algorithm which considers the advective driving force outside the thermal balance system and the turbulent driving force and radiant driving force inside the thermal balance system.

Establishment of the Spatial Extension Algorithm for Air
Temperature.The air temperature from the ground to the height of the thermometer shelter is caused by the advective driving force, turbulent driving force, and radiant driving force.The first two driving forces are supposed to be complementary weighted functions for the hybrid air temperature.Advective driving force depends on the advection intensity, and the advection intensity is determined by the factors of wind speed, temperature difference of inside and outside the system, and so on.The physical mechanism of advective driving force is weighted as hybrid.
Advective driving force can be divided into large-scale advective driving force and local-scale advective driving force.Large-scale advective driving force is caused by large synoptic processes; the driving force is the same within the region of 50 km or larger.Local-scale advective driving force is formed by the local circulation and large-scale advective driving force after being assimilated and denatured by the local underlying surface; the driving force is nearly in proportion to the moisture and heat condition of the underlying surface.Advective driving force mixes with local turbulent driving force and becomes a part of the surface driving mechanism.
For the surface turbulent driving force mechanism, first, the surface performs thermal exchange with the close-surface air temperature layer and changes its temperature and then is weighted as hybrid on the vertical direction in the form of turbulence.These two drive mechanisms are similar.As for the radiant driving force mechanism, the temperature rises through absorbing the long-wave radiation of the surface through the moisture and carbon dioxide in the air; this kind of temperature rising acts on the entire atmosphere layer, and comparing to the entire atmosphere layer, the height from the surface to the thermometer shelter is quite thin; therefore, the radiant driving force is small, and the temperature rising of even thinner layer still acts on the sensor at the height of the thermometer shelter through the turbulence in the atmosphere.Thus, the basic equations which approximately handle this kind of hybrid driving mechanism are In the equations,  is the air temperature of the 50 km * 50 km Fengyun satellite pixel;  ad is the temperature contribution of advective driving force on the basis of the 50 km * 50 km large-scale air average value;  rt is the temperature contribution of radiant and turbulent driving forces on the basis of 50 km * 50 km surface temperature average value.  () is complementary weighted function depending on two variables, the advective and the radiant turbulent driving forces; it is also a weighted function between 0 and 1.  is the parameter which indicates the advection intensity; it is implicated in the temperature during the calculation; the specific value of it is not required.Meanwhile, in ( 2) and ( 3),   and   are two MODIS pixel-scale air temperatures of two weather stations' observed data of temperature. ad and  ad are the temperature contributions of advective driving force in the above two 1 km * 1 km pixels. rt and  rt are the temperature contributions of turbulent driving force of the above two 1 km * 1 km pixels, which can be replaced by the radiant temperature of these two pixels.

Determination of the Complementary Function 𝑓 𝑡 (𝑢).
Assume that, in 50 km pixels scale or larger range (this assumption is supported by the verification later), the temperature contributions of advective driving force are the same; that is,  ad =  ad =  ad , . .., and then subtract ( 2) and (3) to get the complementary functions: or by addition to ( 2) and ( 3): 2  () . (7)

Inversion of MODIS Pixel-Scale Air
Temperature   from the 50 km * 50 km Fengyun Pixel.Assume that, in the range of 50 km * 50 km Fengyun satellite pixel or larger, the complementary function and advection driving force of all the to-be-solved MODIS descendent-scale pixels remain unchanged; there is, where subscript  of   , . . .,  +1 , . . .,  + indicates the to-besolved temperature of weather station MODIS pixel.
(2) Temperatures of the two stations  and  are significantly different with the surface radiant temperatures of their located pixels.
(3) The two stations  and  are inside 50 km pixel and the central location is the best.
(4) If there is only one station in the 50 km pixel or even none, then expand to two or even three 50 km pixel ranges, until there are two weather stations meeting the above conditions.
2.1.5.Further Optimization.The target of further optimization is the consistence comparison and difference correction of thermal driving sources of the area that the temperature of the weather station thermometer shelter is able to represent and the area of the MODIS pixel scale that it is located in.Choose the TM satellite which has higher spatial resolution than MODIS satellite for the consistence comparison and difference correction of thermal driving sources.Approximately, consider the surface temperature of the TM satellite six-wave band pixel scale of the date closest to the date of the to-be-inverted MODIS satellite as the nonadvective thermal driving source of the temperature of the thermometer shelters in the weather stations.Its ratio with the nonadvective thermal driving source of the to-be-inverted MODIS satellite can be expressed as: The ratios in the above equations have different physical meanings: ratio in the left equation is the radiant driving force, and the right one is the air hybrid driving force.According to the above modeling process, for the air layer two meters off the ground, air hybrid drive is predominant.Therefore, equation on the right is adopted.For different weather stations, the ratios are also different.Equations are distinguished with  and  as below: The driving force weighted function equation ( 5) should be revised to Then, we have the corresponding advective driving force polynome which is similar to (7): 2  () . (11b)

Establishment of the Spatial Extension Algorithm for
Humidity.Similar to with the air temperature, humidity is also driven by advection and turbulence.The two driving forces should be complementary weighted functions against the hybrid air humidity; it also depends on the advection intensity, which is determined by the humidity difference inside and outside the system as well as other factors.

Basic Equations.
The basic equations are In the equations, ,   and   , are the average humidity of the 50 km pixel and the two humidity of 1 km pixels with weather station humidity data;  pt ,  pt , and  pt are the turbulent driving forces of 1 km pixel and the turbulent driving forces of the two pixels with weather station humidity data;  ad ,  ad , and  ad are the advective driving forces of the 50 km pixel and the driving forces of the two 1 km pixels with weather station humidity data.Assume that, in the 50 km * 50 km Fengyun satellite pixel, the complementary functions and advective driving forces of all the to-be-solved MODIS descendent-scale pixels remain unchanged.There is no radiant driving force for air humidity.

Obtaining of Turbulent Driving Force for Humidity.
By the surface temperature   of a pixel, its saturated water vapor pressure  *  at this temperature, the apparent thermal inertia   of this pixel, and the apparent thermal inertia  max of the fully wet bear soil pixel, it can be confirmed that the physical meaning of  pt ,  pt , . . .,   / max is the ratio of the two pixels' thermal inertia, which, under the premise that the parent materials of the soil are the same, expresses the water content's ratio of the layer that from the depth of the soil with constant daily temperature changes to the surface, that is, the ratio of the soil layer's water supply ability and which in fact is the ratio of local air humidity driving force.Field experiments show that the bear soil's thermal inertia has linear relationship with the soil's water content; the local driving force for humidity is the vaporized value of soil water content; therefore, the linear relationship between the ratio of the humidity and the ratio of the thermal inertia is on a theoretical and experimental basis.
In addition, it is worth pointing out that in reality the surface is often covered with vegetation of a certain proportion.The thermal inertia calculated by the soil thermal inertia formula when there is vegetation cover is very different with bear soil on the quantity expressing the soil water content; it is called pseudothermal inertia.Practice shows that the ratio between the true and pseudo-thermal inertias expressing the same soil water content reaches 6-8 [6].In other words, the pseudo-thermal inertia is six to eight times more sensitive than the true thermal inertia on expressing the soil water content.
In conclusion, the driving force for local surface humidity can be expressed by In the equation,  tru  ,  fak  ,   , ,  0 ,   ,  1 ,  2 ,  1 , and  2 are, respectively, the soil's true thermal inertia, vegetation pseudo-thermal inertia, vegetation coverage, ratio of true and pseudo-thermal inertias, direct radiation of the sun, surface albedo ratio, thermal inertia calculation start time (usually early morning), thermal inertia calculation end time (usually noon) in seconds, surface temperature at the calculation start time, and the surface temperature at the calculation end time.
Likewise, for the other pixel, Besides that the subscript is changed from pixel  to pixel ; other meanings are all the same with (15).

Obtaining of the Advective Driving Force for Humidity.
Obtain from ( 13) and ( 14), respectively, the advective driving forces for humidity which are 2.2.5.Further Optimization.The target of further optimization is the consistence comparison and difference correction of turbulent humidity driving source of the area that the humidity of the thermometer shelter in the weather station is able to express and the area of the MODIS pixel scale that it is located in.
Choose the TM satellite which has higher spatial resolution than MODIS satellite for the consistence comparison and difference correction of turbulent humidity driving source.Approximately, consider the ratio of the saturated water vapor pressure and thermal inertia under the surface temperature of the TM satellite six-wave band pixel scale of the date closest to the date of the to-be-inverted MODIS satellite as the turbulent humidity driving source of the temperature of the thermometer shelter in the weather station.Its ratio with For different weather stations, the ratios are different.The differences are marked with  and  as below: Driving force weighted function equation ( 5) should be revised to Advective driving force polynome is as follows:  driving force remain unchanged for all the to-be-solved MODIS descendent-scale pixels; then, The   ,  +1 , . . .,  + of all the to-be-solved MODIS descent-scale pixels in the 50 km * 50 km Fengyun satellite pixel can be obtained by In the equations,   ,  +1 , . . .,  + are the descendentscale humidity of the to-be-solved pixels.Surface station observation sites ,  pt , and  pt are together with   and   ; no site selection is required.

Validation of the Algorithm.
Choose MODIS image with a range that contains three weather stations; the image could be larger than the 50 km * 50 km Fengyun satellite pixel.Use the air temperature and humidity of any two stations to extrapolate the air temperature and humidity of the third weather station, and compare then with the actual  observed value of the air temperature and humidity of the third station.Comparative verifications are performed to the three combinations of the three stations.The locations of the meteorological stations in North China Plain are given in Figure 1. Figure 2 shows the interpolated results of air temperature and air humidity on the basis of the method described in Section 2. Values of air temperature and air humidity at the pixel, where weather stations are located, were extracted.Table 1 and Figure 3 show the comparison of the calculated value of air temperature and humidity model of North China Plain and the actual observed value.It can be seen that the average absolute difference of the algorithm is 0.65 k and 0.31 mb, respectively, for surface air temperature and humidity expansion.

Conclusions
Practice shows that using the observed air temperature and humidity data of two adjacent weather stations and corresponding surface temperatures inverted with surrounding MODIS and TM satellite pixels, the air temperature and humidity of surrounding MODIS satellite pixels could be extrapolated.

Figure 1 :
Figure 1: Weather station distribution in North China Plain.

Figure 2 :
Figure 2: Surface air temperature (a) and humidity distribution (b) in North China Plain.

Figure 3 :
Figure 3: Comparison of modeled and observed surface air temperature (a) and water vapor pressure distribution (b) over the meteorological stations.

Table 1 :
Comparison of the calculated value of air temperature and humidity model of North China Plain and the actual observed value.