A method for calculating land surface evapotranspiration based on remote sensing

Abstract

The present invention belongs to the field of geography and ecology, and relates to a method for calculating land surface evapotranspiration based on remote sensing. Land surface evapotranspiration consists of vegetation transpiration, soil evaporation and intercepted evaporation. The present invention mainly optimizes the temperature in vegetation evaporation and the soil moisture in soil evaporation, and mainly includes two parts: using a physical model to calculate the leaf temperature instead of the ambient temperature in the original method to calculate vegetation transpiration; using the soil moisture content of different soil layers to redefine the soil moisture limit in the soil evaporation process. The method of the present invention is simple and practical, and the data involved can be obtained free of charge from the remote sensing platform, China Meteorological Data Network, vortex flux tower and related literature. The present invention belongs to a mechanism model, and no parameter correction is required during extrapolation, and reliable estimation of land surface evaporation can be made quickly and at a lower cost.

Description

Technical Field

[0001] The present invention belongs to the field of geography and ecology, and relates to a method for calculating evapotranspiration of land surface. Specifically, the method uses remote sensing technology to obtain relevant MODIS and FLDAS data, obtains site relative humidity data through the National Meteorological Center, obtains evapotranspiration data through the China Flux Tower and related literature, and optimizes the soil water molecule model and temperature sub-model based on the existing Priestley-Taylor model framework, and uses the optimized model for estimating evapotranspiration of land surface. Background Art

[0002] Evapotranspiration (E) is an important part of the land surface water and heat cycle, which consists of vegetation transpiration, soil evaporation and intercepted evaporation. The maximum evapotranspiration under the most suitable environment is called potential evapotranspiration; the evapotranspiration under the actual conditions of climate, vegetation and land cover is called actual evapotranspiration. More than 50% of solar radiation and 60% of the total annual precipitation are returned to the atmosphere through actual evapotranspiration.

[0003] The relationship between evapotranspiration and available energy can be expressed as:

[0004] LE=R n -GH

[0005] Where, L represents the latent heat of evaporation (MJ / kg), E represents evaporation (mm), LE represents the latent heat flux (W / m 2 ), R n Indicates net radiation (W / m 2 ), G represents soil heat flux (W / m 2 ), H represents the sensible heat flux (W / m 2 ).

[0006] The vortex flux tower is the best way to observe evapotranspiration, but there are only 400 flux towers in the world, with short recording time, high cost, and it is difficult to expand from point scale to surface scale. To obtain spatially continuous evapotranspiration, one way is to increase the density of flux towers and interpolate, but this method is costly, time-consuming and labor-intensive, and can only obtain data after the establishment of the flux tower, which is not very feasible; another method is to interpolate meteorological data through inverse distance weighting, kriging, spline functions and other methods, and then obtain evapotranspiration based on raster data of meteorological data and models such as Penman-Monteith, but the premise of interpolation is that the terrain is flat, the spatial heterogeneity is small, and all elements required for evapotranspiration have enough and continuous station recording data. In fact, in some developing countries and regions, the number of meteorological stations is small and the time span is short, so the interpolation effect is greatly limited, especially in the Qinghai-Tibet Plateau, where the terrain and climate change are extremely large, and there is a lack of meteorological stations in the northwest of the Qinghai-Tibet Plateau, so the interpolation effect will be worse. Therefore, free remote sensing data with high temporal and spatial resolution over a long period of time has become an important data source for estimating high temporal and spatial resolution evapotranspiration. Evapotranspiration estimation based on remote sensing includes statistical models and mechanism models.

[0007] Statistical models can be divided into two categories: (1) Regression model: The regression model establishes a regression relationship between the existing flux tower evapotranspiration data and meteorological data such as temperature and precipitation, and remote sensing index data such as NDVI, and applies the regression relationship to the whole region to estimate evapotranspiration. The regression model has a simple structure and can obtain evapotranspiration through simple raster operations. However, the parameters and structure of the regression model are based on the existing measured data. When applied to new regions or new climate zones, the parameters often need to be recalibrated. (2) Machine learning model: Using methods such as random forests, artificial neural networks, and support vector machines, the model is trained by establishing a nonlinear relationship between meteorological factors, topographic factors and evapotranspiration, and the trained model is extrapolated to other regions. However, the machine learning model ignores the internal physical mechanism and it is difficult to explain the interaction between evapotranspiration and driving factors. It is also limited by the quality, quantity and typicality of the data. When Chen et al. (2014) evaluated the effects of eight models, three machine learning models (artificial neural networks, support vector machines and decision trees) had the worst simulation accuracy.

[0008] Mechanism models include the surface energy balance algorithm for land (SEBAL), the Penman-Monteith (PM) model, and the Priestley-Taylor (PT) model. The surface energy balance model involves some microclimate parameters that can only be obtained through experience during the calculation process, such as dynamic roughness, thermodynamic roughness and other parameters. Therefore, the development of this method is limited. On this basis, Penman combined the surface energy balance model with the principles of aerodynamics to propose a method for calculating evapotranspiration of a wet underlying surface. Monteith introduced the concept of surface resistance and proposed the Penman-Monteith model. Theoretically, the PM model has a clearer mechanism than the PT model, but when remote sensing is used for full-area calculations, it is difficult to quantify parameters such as stomatal conductance and surface conductance in the PM model. Compared with the PM model, the PT model requires fewer parameters and can achieve approximate effects. Therefore, the PT model has been more widely used, especially the Priestley-Taylor Jet Propulsion Lab (PT-JPL) model with a simple structure and clear mechanism (Fisher et al., 2008; Yao et al., 2013; Yao et al., 2015).

[0009] All models, including PT-JPL, have the same problem: all models are calculated based on land surface temperature. However, according to the actual process, the leaf temperature (T) is the real factor that determines the transpiration of vegetation. l ) rather than air temperature (T a )(Leuning et al., 1995; Tuzet et al., 2003), and generally leaf temperature is not equal to air temperature (Blonder and Michaletz, 2018; Michaletz et al., 2016; Michaletz et al., 2015). In addition, the relationship between air temperature and leaf temperature is not fixed. Leaf temperature is sometimes 7-10°C higher than air temperature (Korner et al., 2002; Gates et al., 1964; Lange, 1959), and sometimes even 20°C higher (Doughty and Goulden, 2008). The regulation of leaf temperature with air temperature is actually the result of plant adaptation, which will be beneficial to plant growth, survival and production (Blonder and Michaletz, 2018). Therefore, using air temperature instead of leaf temperature will misjudge the impact of temperature on transpiration, especially in places where there is a large difference between air temperature and leaf temperature.

[0010] In addition to leaf temperature, soil moisture has a very important influence on evapotranspiration (Detto et al., 2006; Krishnan et al., 2006; Wang and Liang, 2008). Many studies have suggested using remote sensing inversion of highly spatially heterogeneous soil moisture to estimate soil evaporation (Gokmen et al., 2012; Jin et al., 2011; Miralles et al., 2011; Yao et al., 2015). However, previous satellite-based soil moisture only considered the surface soil moisture of about 2-5 cm, ignoring the impact of deeper soil moisture on soil evaporation (Anderson et al., 2007; Yao et al., 2015). In recent years, emerging soil moisture data of different soil layers have provided a direct and effective data source for estimating soil evaporation. Vapor pressure deviation can directly describe the dryness of the air. Under the same moisture conditions, the higher the temperature, the greater the VPD. Summary of the invention

[0011] The technical problem to be solved by the present invention is to provide a land surface evapotranspiration calculation method based on remote sensing with higher accuracy.

[0012] The present invention optimizes the temperature sub-model and soil water molecule model in the mechanism model PT-JPL, and applies the model to the calculation of potential evapotranspiration and actual evapotranspiration of the Qinghai-Tibet Plateau, and compares the actual evapotranspiration with the observation data (true value) of the vortex flux tower of the Qinghai-Tibet Plateau for accuracy verification. The accuracy of the model is compared with that of the FLDAS, MOD16A2 and other models, and the effect of the optimization model of the present invention is verified by three indicators: root mean square error, deviation and determination coefficient.

[0013] The technical solution of the present invention is:

[0014] A method for calculating land surface evapotranspiration based on remote sensing, actual evapotranspiration and vegetation surface transpiration E s , soil surface evaporation E s and interception of evaporation E i The following model is used to calculate:

[0015] LE=L(E c +E s +E i )

[0016]

[0017]

[0018]

[0019] Where L represents the latent heat of evaporation (MJ / kg), E represents evaporation (mm), and the product of the latent heat of evaporation and evaporation, LE, represents the latent heat flux (W / m 2 );Ec is the transpiration of dry vegetation surface (mm), E s is the evaporation from dry soil surface (mm), E i Indicates intercepted evaporation (mm); LE c , L.E. s With LE i Respectively represent the latent heat flux of each part; α PT is a dimensionless coefficient, generally taken as 1.26; Δ is the slope of the curve of saturated water vapor pressure versus temperature (kPa / ℃); γ is the psychrometric constant (kPa / ℃); R nc is the net radiation of the vegetation surface (W / m 2 );R ns is the net radiation from the soil surface (W / m 2 ); G is soil heat flux (W / m 2 ).

[0020] From the above formula, the potential evapotranspiration (PE) and the estimated evapotranspiration (E) are:

[0021]

[0022]

[0023] PE in k c , k i With k s Down to E, k c , k i With k s They represent the limiting factors from potential evaporation of vegetation to actual evaporation, from potential interception evaporation of vegetation surface to actual interception evaporation, and from potential evaporation intercepted by soil to actual evaporation. The calculation methods are:

[0024] k c =(1-f wet )f g f T f M

[0025] k s =f wet +(1-f wet )f sm

[0026] k i =f wet

[0027] Among them, f wet represents the interception surface ratio; f g 、fT 、f M They represent the stress coefficients of greenness, temperature and water on the plant surface respectively; f sm is the soil surface water stress coefficient. The value range of the above five indicators is 0-1. When the value is 1, it means there is no environmental stress; the smaller the value, the greater the environmental stress.

[0028] The calculation method of each factor is as follows:

[0029] (1) Temperature stress coefficient (f T )

[0030] After improvement T The calculation method is as follows:

[0031]

[0032] Among them, T l represents the blade temperature; T opt represents the optimum temperature; f T The value range is 0-1.

[0033] (2) Interception surface ratio (f wet )

[0034]

[0035] Where R h Relative humidity is calculated as follows:

[0036]

[0037] Where p is the atmospheric pressure, q is the absolute humidity, T a is the air temperature.

[0038] (3) Leaf temperature (T l )

[0039]

[0040] Among them, T a is the air temperature; R nc is the net radiation falling on the vegetation canopy; c P is the specific heat of air at constant pressure; δ b is the blade boundary conductance, Among them, μ represents the wind speed and d is related to the blade width.

[0041] (4) Vegetation greenness (f g )

[0042] Vegetation greenness refers to the part of vegetation transpiration that is absorbed by green vegetation. The calculation method of vegetation greenness is:

[0043]

[0044] Among them, f g Indicates vegetation greenness; f APAR represents the proportion of land surface solar radiation absorbed by living vegetation; f IPAR Represents the proportion of all solar radiation intercepted by vegetation. The specific expressions of each factor are:

[0045] f APAR =1.36SAVI-0.048

[0046] f IPAR =NDVI-0.05

[0047] SAVI=1.5×(ρ NIR -ρ Red ) / (ρ NIR +ρ Red +0.5)

[0048] DNVI=(ρ NIR -ρ Red ) / (ρ NIR +ρ Red )

[0049] Where SAVI is the soil corrected vegetation index, ρ NIR and ρ Red are the reflectances of the near-infrared band and the red band, respectively.

[0050] (5) Vegetation water stress coefficient (f M )

[0051] The expression for vegetation water limitation is:

[0052]

[0053] in, Indicates the maximum value of the percentage of solar radiation absorbed by the green vegetation in the pixel.

[0054] (6) Soil moisture stress coefficient (f sm )

[0055] The soil water stress optimized by this method is:

[0056]

[0057] Where SM represents soil moisture content, SM max With SM min They represent the maximum and minimum values ​​of soil moisture content respectively.

[0058] (7) Psychrometric constant (γ)

[0059]

[0060] Among them, C p is the specific heat of air at constant pressure, which means the energy required to increase the temperature of a unit volume of air by 1°C under a certain air pressure; P is the atmospheric pressure, which is related to the altitude; λ is the latent heat of evaporation; ε is the ratio of the molecular weight of water vapor to the molecular weight of dry air.

[0061] (8) Vegetation surface (R nc ) and soil surface (R ns ) can absorb the following solar radiation:

[0062] R nc =FVC×R n

[0063] R ns =(1-FVC)×R n -G

[0064] Among them, FVC represents vegetation coverage, R n is the vegetation net radiation, and G is the soil heat flux. BRIEF DESCRIPTION OF THE DRAWINGS

[0065] Figure 1 , comparative verification of remote sensing simulation and meteorological station measured relative humidity.

[0066] Figure 2 , validation of latent heat flux, comparison of true and simulated values ​​between different models: a) Fisher; b) MOD16A2; c) TerraClimate; d) FLDAS; e) PT-Hybrid; f) this paper. DETAILED DESCRIPTION

[0067] The present invention takes the Qinghai-Tibet Plateau as a research object.

[0068] 1 Overview of the study area

[0069] The Qinghai-Tibet Plateau is known as the "Roof of the World". It is located in the southwest of China, with longitude 26°03'N~40°00'N and latitude 73°06'E~105°08'E. The average altitude is above 4,000 meters and covers an area of ​​2572.4×10 3 km 2 , is the source of the Yangtze River, the Yellow River, the Lancang River and many other rivers, and 37% of China's surface runoff comes from the Qinghai-Tibet Plateau. The important water conservation function of this region makes the assessment of evapotranspiration very important for accurately understanding the water and heat cycle of the entire ecosystem.

[0070] 2 Data and Data Preprocessing

[0071] 2.1 Meteorological data

[0072] Meteorological data comes from the China Meteorological Data Network's China Surface Climate Data Daily Dataset (V3.0), which collects eight types of information, including evaporation, temperature, and relative humidity, from nearly 700 meteorological stations across the country. Due to research needs, this paper uses MATLAB to calculate the arithmetic mean of the daily relative humidity data of 65 meteorological stations in the Qinghai-Tibet Plateau region from 2000 to 2019 to obtain monthly data.

[0073] 2.2 Remote sensing data

[0074] Remote sensing data include MODIS raw data and FLDAS and TerraClimate reanalysis datasets. MODIS is an optical sensor carried on the Terra satellite (10:30) and the Aqua satellite (1:30). The sensor has 36 bands with wavelengths ranging from 0.415 μm to 14.235 μm and spatial resolutions ranging from 250 m to 1000 m (Running et al., 1994; Wang et al., 2011). MODIS provides a variety of data products such as NDVI (MOD13Q1), land surface temperature (MOD11A2), and leaf area index (MOD15A2), which have become an important data source for estimating evapotranspiration (Chen et al., 2014; Hao et al., 2019; Yao et al., 2015). The FLDAS dataset (NcNally et al., 2017) is designed to help developing countries with limited data to conduct food security assessments. FLDAS uses the Noah 3.6.1 ground model. Its CHIRPS-6 hourly precipitation has been reduced using the NASA Land Surface Data Toolkit. It is part of the NASA Land Information System (LIS) and provides a variety of monthly data products such as evapotranspiration, downward solar radiation, soil heat flux, and water content in different soil layers from 1982 to the present. The TerraClimate dataset (Abatzoglou et al., 2018) provides the maximum air temperature, minimum air temperature, and wind speed at 10 m over the global land surface.

[0075] The data products and related descriptions are shown in Table 1. The near-infrared band NIR and the red band Red of MOD13Q1 can be used to calculate the SAVI and NDVI required for vegetation greenness, and the enhanced vegetation index EVI can be used to solve the vegetation coverage. The photosynthetically active radiation component FPAR of MOD15A2 and the net radiation of FLDAS can be used to calculate the total amount of solar radiation absorbed by vegetation, and then the water limitation factor of vegetation can be calculated. SRTMDEM elevation is used to solve the atmospheric pressure, and then solve the psychrometric constant, and it is also used with FLDAS moisture content and TerraClimate air temperature to solve relative humidity. TerraClimate wind speed is used to solve leaf temperature. FLDAS soil moisture data is used to solve the moisture limitation during soil evaporation. Evaporation from MOD16A2, latent heat of evaporation from TerraClimate and evaporation from FLDAS are used to compare and verify the accuracy of the model.

[0076] In order to compare the accuracy of FLDAS soil moisture content, this paper also collected a downscaled dataset (SoilMoisture in China dataset, hereinafter referred to as CN_SM) and a global data product of the European Space Agency (ESACCI SM v06.1, hereinafter referred to as ESA_SM). The data sources and data units are shown in the table.

[0077] Table 1 Data sources

[0078]

[0079]

[0080] 3 Flux data

[0081] The accuracy of the evapotranspiration algorithm was verified using ground-measurements of eddy covariance (ECOR) data from the China Flux Tower. The evapotranspiration data came from the China Flux Tower and related literature. The data from the China Flux Tower were downloaded from relevant websites. The data included three observation stations: Haibei Shrub, Haibei Wetland, and Dangxiong. Each observation station provided half-hourly observation records, daily observation records, monthly observation records, and annual observation records. The observation data for Tanggula, Maqu, and Nagqu came from relevant literature. The specific information and related literature are shown in Table 2.

[0082] Table 2 ECOR ground observations from flux towers and related studies in China

[0083]

[0084] Among them, the LE data provided by ECOR has the problem of energy imbalance, that is, the difference between solar radiation and soil heat flux is not equal to the sum of latent heat flux and sensible heat flux. That is:

[0085] R n -G≠LE+H

[0086] Therefore, before using LE, a Bowen ratio correction is required to convert the observed LE into the actual LE (Twine et al., 2000). The correction method is as follows:

[0087]

[0088]

[0089] Among them, LE r R represents the latent heat flux after Bowen ratio correction, LE and H represent the observed latent heat flux and sensible heat flux, respectively. n represents the net radiation of the land surface, and G represents the soil heat flux. In order to meet the needs of long-term series research, the flux tower observation data and the data obtained from related literature were finally converted to the monthly scale.

[0090] 4 Model structure

[0091] After Priestley and Taylor conducted experimental studies in open water and arid areas, they determined that the basic structure of the PT model is:

[0092]

[0093] L is the latent heat of evaporation, which is 2.45 MJkg -1 ; E is land surface evapotranspiration; α PT is a dimensionless coefficient, ranging from 1.25 to 1.31, and generally 1.26; Δ is the slope of the curve of saturated water vapor pressure versus temperature; γ is the psychrometric constant; R n is the net radiation on the vegetation surface.

[0094] Since the PT model calculates potential evapotranspiration, when converting potential evapotranspiration into actual evapotranspiration, the effects of environmental factors such as temperature, moisture and vegetation greenness on LE are not taken into account. Therefore, Fisher et al. (2008) proposed the PT-JPL model, which assumes that actual evapotranspiration E is related to canopy surface transpiration (E c ), soil surface evaporation (E s ) and interception loss (E i ), the specific calculation method is as follows:

[0095] LE=L(E c +E s +E i )

[0096]

[0097]

[0098]

[0099] Where, L is the latent heat of evaporation (MJ / kg); E is evaporation (mm); LE is the product of latent heat of evaporation and evaporation, which represents the latent heat flux (W / m 2 );E c is the canopy surface transpiration (mm); E s is soil surface evaporation (mm); E i Indicates intercepted evaporation (mm); LE c , L.E. s With LE i Respectively represent the latent heat flux of each part; α PT is a dimensionless coefficient, generally taken as 1.26; Δ l , Δ a With Δ s are the slopes of the curves of saturated water vapor pressure versus leaf temperature, air temperature, and soil temperature (kPa / ℃); γ is the psychrometric constant (kPa / ℃); R nc is the net radiation of the canopy surface (W / m 2 );R nc is the net radiation from the soil surface (W / m 2 ); G is soil heat flux (W / m 2 ).

[0100] From the above formula, the potential evapotranspiration (PE) and the estimated evapotranspiration (E) are:

[0101]

[0102]

[0103] PE in k c , k i With k s Down to E, k c , k i With k s They represent the limiting factors from potential evaporation of vegetation to actual evaporation, from potential interception evaporation of vegetation surface to actual interception evaporation, and from potential evaporation intercepted by soil to actual evaporation. The calculation methods are:

[0104] k c =(1-f wet )f g f T f M

[0105] k s =f wet +(1-f wet )f sm

[0106] k i =f wet

[0107] Among them, f wet represents the interception surface ratio; f g f T f M They represent the restriction coefficients of greenness, temperature and water on the plant surface respectively; f sm is the moisture limitation coefficient of the soil surface. The value range of the above five indicators is 0-1. When the value is 1, it means there is no environmental stress; the smaller the value, the greater the environmental stress.

[0108] The specific models involved in the above parameter settings are calculated below:

[0109] (1) Temperature stress coefficient (f T )

[0110] The model itself T The calculation method is:

[0111]

[0112] Among them, T a is the air temperature; T opt Indicates the optimum temperature. Since it is the temperature of the vegetation surface rather than the air temperature that actually affects the transpiration of vegetation, the present invention uses the temperature of the vegetation surface to replace the air temperature to calculate the temperature stress coefficient. T The calculation method is as follows:

[0113]

[0114] Among them, T l represents the blade temperature; T opt represents the optimum temperature; f T The value range is 0-1.

[0115] (2) Interception surface ratio (f wet )

[0116] The interception surface ratio is related to the relative humidity of the air:

[0117]

[0118] Among them, R hRelative humidity. Relative humidity is an important data source for calculating water vapor pressure difference, but many remote sensing products currently only provide absolute humidity. Therefore, it is necessary to convert absolute humidity to relative humidity. The derivation equation used is:

[0119]

[0120] Where w is the actual water vapor pressure, w s is the saturated water vapor pressure, p is the atmospheric pressure, q is the absolute humidity, T a is the air temperature.

[0121] (3) Leaf temperature (T l )

[0122]

[0123] Among them, T l represents the blade temperature, T opt Indicates the optimum temperature. T The value range is 0-1.

[0124] Leaf temperature affects vegetation transpiration and vegetation photosynthesis, and vegetation transpiration in turn regulates leaf temperature. The change of leaf temperature is nonlinearly related to factors such as air temperature and land surface solar radiation. This paper adopts the method of Dong et al. (2017) to convert leaf temperature to air temperature:

[0125]

[0126] Among them, T a is the air temperature; R nc is the net radiation falling on the vegetation canopy; c P is the specific heat of air at constant pressure; δ b is the blade boundary conductance, Among them, μ represents the wind speed and d is related to the blade width.

[0127] (4) Vegetation greenness (f g )

[0128] Vegetation greenness refers to the part of vegetation transpiration that is absorbed by green vegetation. The calculation method of vegetation greenness is:

[0129]

[0130] Among them, f g represents vegetation greenness, f APAR represents the proportion of land surface solar radiation absorbed by green vegetation, f IPAR Represents the proportion of all solar radiation intercepted by vegetation. The specific expressions of each factor are:

[0131] f APAR=1.36SAVI-0.048

[0132] f IPAR =NDVI-0.05

[0133] SAVI=1.5×(ρ NIR -ρ Red ) / (ρ NIR +ρ Red +0.5)

[0134] NDVI=(ρ NIR -ρ Red ) / (ρ NIR +ρ Red )

[0135] Where SAVI is the soil corrected vegetation index, ρ NIR and ρ Red are the reflectances of the near-infrared band and the red band, respectively.

[0136] (5) Vegetation water stress coefficient (f M )

[0137] The expression for vegetation water limitation is:

[0138]

[0139] in, Indicates the maximum value of the percentage of solar radiation absorbed by the green vegetation in the pixel.

[0140] (6) Soil moisture stress coefficient (f sm )

[0141] In the original method, f sm The soil water stress coefficient is defined as:

[0142] f sm =R h VPD

[0143] Among them, R h represents relative humidity and VPD is the vapor pressure deviation.

[0144] R hBoth VPD and VPD can reflect the dryness of the air, and can be used as an approximate substitute but cannot directly reflect the soil water stress condition. In order to improve the reliability of the evaluation results, we directly use the remote sensing soil water data developed in recent years to quantify soil water stress. There is the following relationship between soil water, VPD and soil water stress (water stress coefficient): the higher the soil water content, the smaller the water stress (the larger the water stress coefficient); the higher the VPD, the greater the demand for soil water, and the smaller the soil water stress (the larger the water stress coefficient). Therefore, soil water stress and VPD are both proportional to the soil water stress coefficient. Under the synergistic effect of the two, this method defines the soil water stress equation as:

[0145]

[0146] Where SM represents soil moisture content, SM max With SM min They represent the maximum and minimum values ​​of soil moisture content respectively.

[0147] Since the research object is grassland with shallow root depth, the soil water content is obtained from 0-10cm and 10-40cm. The water content of the 0-40cm soil layer is expressed as a weighted average:

[0148] SM=0.25×SM1+0.75×SM2

[0149] Among them, SM1 is the surface soil water of 0-10 cm, and SM2 is the soil water of 10-40 cm.

[0150] SM max With SM min Represent the maximum and minimum values ​​of soil moisture content respectively. This paper calculates the maximum value of soil moisture in 20 years as SM max , the minimum soil moisture is taken as SM min .

[0151] (7) Psychrometric constant (γ)

[0152] The expression of the psychrometric constant is:

[0153]

[0154] Among them, C p is the specific heat of air at constant pressure, which means the energy required to increase the temperature of a unit volume of air by 1°C under a certain air pressure; P is the atmospheric pressure, which is related to the altitude; λ is the latent heat of evaporation; ε is the ratio of the molecular weight of water vapor to the molecular weight of dry air.

[0155] (8) Vegetation surface (R nc ) and soil surface (R ns ) can absorb the solar radiation respectively

[0156] R nc =FVC×R n

[0157] R ns =(1-FVC)×R n -G

[0158] Among them, FVC represents vegetation coverage, R n is the vegetation net radiation, and G is the soil heat flux.

[0159] The expression for calculating vegetation coverage based on EVI is:

[0160]

[0161]

[0162] G = a g (1-FVC)R n

[0163] Among them, EVI stands for Enhanced Vegetation Index, EVI max With EVI min Respectively represent the maximum and minimum values ​​of the enhanced vegetation index of the pixel; K, L, C1, and C2 are 2.5, 1, 6.0, and 7.5, respectively; a g is a constant of 0.18, FVC is the vegetation coverage, R n is the net radiation.

[0164] 5 Data analysis and evaluation methods

[0165] The present invention uses root mean square error, deviation and R 2 Three statistical analysis methods were used to evaluate the simulation effect of evapotranspiration.

[0166] The root mean square error (RMSE), also known as the standard error, explains how concentrated the simulated values ​​are around the observed values. It is the square root of the ratio of the sum of the squares of the deviations between the observed values ​​and the simulated values ​​to the number of observed values.

[0167]

[0168] Among them, RMSE represents the root mean square error, y i is the ith observation value, is the ith simulated value, and N is the number of observed values.

[0169] Bias describes the degree of deviation between the observed value and the simulated value, reflecting the fitting ability of the model itself.

[0170]

[0171] Among them, Bias represents the deviation, Represents the simulated value, y i represents the observed value, and N is the number of observed values. The unit of the deviation is the same as the observed value.

[0172] R 2 Also called goodness of fit, it represents the percentage of observed values ​​that can be explained by simulated values.

[0173]

[0174] Among them, R 2 The range is 0-1, and the larger the value, the higher the explained percentage. Represents the simulated value, y i represents the observed value, N is the number of observed values, Represents the mean of all observations.

[0175] 6. Result verification and analysis

[0176] 6.1 Relative humidity verification

[0177] This section mainly verifies the accuracy of the relative humidity calculation equation. The verification data uses the relative humidity station data of the National Meteorological Data Center. The data includes 14,400 data records from 2000 to 2019. The comparative analysis results are as follows Figure 1 As shown, the root mean square error RMSE is 8.55% and the determination coefficient R 2 The value is 0.401, and the deviation is 2.80%. This proves that the conversion from absolute humidity to relative humidity is reliable.

[0178] 6.2 Soil moisture verification

[0179] The results are shown in Table 3. We calculated the correlation between soil moisture and observation data in the FLDAS dataset (FLDAS_SM), SMC dataset (SMC_SM) and ESACCI SM v06.1 (ESA_SM), and obtained the average correlation values ​​of different models by weighted average method. Before comparison, we converted the daily ESA_SM values ​​into monthly values, but there were still many data gaps. Therefore, for the objectivity of the results, we only used their intersection for comparison and verification. The vacancy values ​​of each station are different, and the number of samples (n) of each station is not exactly the same as that of other stations. It should also be noted that at the Qingshuihe station, it has only one valid record, so accuracy verification is not performed here. Among the comparison results of the three, the FLDAS_SM dataset has the highest accuracy at Guide Station, and the ESA_SM has the worst accuracy at Yushu. Overall, the average coefficients corresponding to FLDAS_SM, SMC_SM, and ESA_SM are 0.46, 0.34, and 0.40, respectively.

[0180] Table 3 Comparison of correlation between three soil moisture products and observed data

[0181]

[0182] 6.3 Comparison of different evapotranspiration models

[0183] In order to verify the effect of the model of the present invention, MOD16 (Mu et al., 2007; Mu et al., 2011), PT-JPL (Fisher et al., 2008), PT-Hybrid (Yao et al., 2015) and FLDAS (McNally et al., 2017) are used for comparison. 2 Statistical analysis was performed.

[0184] The statistical results of latent heat flux observed by various remote sensing models and flux towers are shown in Table 4 and Figure 2 Table 4 compares the deviations and root mean square errors of the six stations: Tanggula, Maqu, Naqu, Haibei shrub, Haibei wetland and Dangxiong. Figure 2 The R of the six sites was compared. 2 .

[0185] From Table 4, Dangxiong has the best simulation effect, with a deviation of 2.69W / m 2 , the root mean square error is 10.15W / m 2; followed by Haibei Wetland and Haibei Shrub; the simulation results of Maqu, Naqu and Tanggula are poor, because the data of these three stations are from relevant literature, and there is no data such as sensible heat flux and soil heat flux based on the station, so no Bowen ratio correction is performed. From the comparison between the models, the simulation results of this paper are the best, FLDAS has better results, and MOD16 and PT-Hybrid have the worst results. Since FLDAS itself is based on the development of African agriculture (McNally et al., 2017), it is suitable for arid and semi-arid areas. MOD16 is a global general model based on PM, involving parameters such as stomatal conductance that are difficult to quantify, which ultimately manifests in poor local simulation results. When PT-Hybrid changes meteorological data, the parameters of the model often need to be recalibrated, so when simulating the data in this paper, its effect is poor.

[0186] from Figure 2 It can be seen that the overall R of the PT-JPL model is 2 The overall R of MOD16 is 0.727. 2 The overall R of TerraClimate is 0.547. 2 The overall R of PT-Hybrid is 0.878. 2 The overall R 2 is 0.903. From the overall distribution trend, the PT-JPL model, MOD16 and PT-Hybrid all have a relatively obvious underestimation trend. The simulated values ​​of FLDAS are consistent with the observed values, and the simulation effect is closest to the simulation results of this paper. Only a few individual values ​​of FLDAS are better than this paper. However, since FLDAS is more suitable for arid and semi-arid areas, the model will have a certain underestimation at the Haibei Wetland Station. This paper directly quantified the soil water content, and still has a good simulation effect in the Haibei Wetland.

[0187] Table 4 PT-JPL, MOD(MOD16A2), TC(TerraClimate), FLDAS, PT H Comparison of RMS error and determination coefficient between (PT-Hybrid) and the present invention

[0188]

[0189] The present invention optimizes a calculation method for evapotranspiration, replaces the air temperature in vegetation transpiration with leaf temperature, and replaces the water limiting factor in soil evaporation with normalized soil moisture. The present invention has the following advantages:

[0190] 1. Data is easy to obtain. The calculation of evapotranspiration depends entirely on free remote sensing data, meteorological data and flux tower data, which can be downloaded free of charge from relevant websites.

[0191] 2. The calculation results are more reliable. The present invention uses the soil moisture content of different soil layers to quantify the effect of soil moisture on soil evaporation; uses the mechanism model to calculate the nonlinear relationship between leaf temperature and air temperature, and uses leaf temperature to calculate the temperature limit in vegetation transpiration. The theoretical basis and result accuracy of this calculation method are more reliable.

[0192] 3. The present invention adopts remote sensing method to calculate evapotranspiration, which can reflect the potential evapotranspiration and actual evapotranspiration status at different time and space scales in real time, and can comprehensively reflect the moisture status of a region, which is conducive to accurately quantifying the water cycle process.

Claims

1. A method for calculating land surface evapotranspiration based on remote sensing, characterized in that: Actual evapotranspiration (E) and canopy surface transpiration (E c ), soil surface evaporation (E s ) and interception evaporation (E i ), the specific calculation method is as follows: THE=L(E c +E s +E i ) Where, L is the latent heat of evaporation (MJ / kg); E is evaporation (mm); LE is the product of latent heat of evaporation and evaporation, which represents the latent heat flux (W / m 2 );E c is the canopy surface transpiration (mm); E s is soil surface evaporation (mm); E i Indicates intercepted evaporation (mm); LE c , L.E. s With LE i Respectively represent the latent heat flux of each part; α PT is a dimensionless coefficient, generally taken as 1.26; Δ l , Δ a With Δ s are the slopes of the curves of saturated water vapor pressure versus leaf temperature, air temperature, and soil temperature (kPa / ℃); γ is the psychrometric constant (kPa / ℃); R nc is the net radiation of the canopy surface (W / m 2 );R ns is the net radiation from the soil surface (W / m 2 ); G is soil heat flux (W / m 2 ); From the above formula, the potential evapotranspiration (PE) and actual evapotranspiration (E) are: PE in k c , k i With k s Down to E, k c , k i With k s They represent the limiting factors from potential evaporation of vegetation to actual evaporation, from potential interception evaporation of vegetation surface to actual interception evaporation, and from potential evaporation intercepted by soil to actual evaporation. The calculation methods are: k c =(1-f wet )f g f T f M k s =f wet +(1-f wet )f sm k i =f wet Among them, f wet represents the interception surface ratio; f g f T f M They represent the stress coefficients of greenness, temperature and water on the plant surface respectively; f sm is the moisture stress coefficient of the soil surface; the value range of the above five indicators is 0-1. When the value is 1, it means there is no environmental stress; the smaller the value, the greater the environmental stress; The calculation and improvement methods of each factor are as follows: (1) Temperature stress coefficient (f T ) Among them, T l represents the blade temperature; T opt represents the optimum temperature; f T The value range is 0-1; (2) Interception surface ratio (f wet ) Among them, R h Indicates relative humidity; (3) Leaf temperature (T l ) Leaf temperature (T l ) is expressed as: Among them, T a is the air temperature; R nc is the net radiation falling on the vegetation canopy; c P is the specific heat of air at constant pressure; δ b is the blade boundary conductance, Among them, μ represents the wind speed, and d is related to the blade width; (4) Vegetation greenness (f g ) Vegetation greenness refers to the part of vegetation transpiration that is absorbed by green vegetation. The calculation method of vegetation greenness is: Among them, f g Indicates vegetation greenness; f APAR represents the proportion of land surface solar radiation absorbed by living vegetation; f IPAR represents the proportion of all solar radiation intercepted by vegetation; (5) Vegetation water stress coefficient (f M ) The expression for vegetation water limitation is: in, Indicates the maximum value of the solar radiation absorbed by the green vegetation in the pixel; (6) Soil moisture stress coefficient (f sm ) Where SM represents soil moisture content, SM max With SM min They represent the maximum and minimum values ​​of soil moisture content respectively; VPD represents vapor pressure deviation; (7) Psychrometric constant (γ) Among them, C p is the specific heat of air at constant pressure, which means the energy required to increase the temperature of a unit volume of air by 1°C under a certain pressure; P is the atmospheric pressure, which is related to the altitude; λ is the latent heat of evaporation; ε is the ratio of the molecular weight of water vapor to the molecular weight of dry air; (8) Solar radiation: vegetation surface (R nc ) and soil surface (R ns ) R nc =FVC×R n R ns =(1-FVC)×R n -G Among them, FVC represents vegetation coverage, R n is the vegetation net radiation, and G is the soil heat flux.

Images