CN102354348A - Watershed scale soil moisture remote sensing data assimilation method - Google Patents

Abstract

The invention discloses a watershed scale soil moisture remote sensing data assimilation method belonging to the field of remote sensing data assimilation methods. The method comprises the following steps of: (A) preparing for data; (B) constructing a watershed soil moisture assimilation observation operator; (C) constructing a distributed hydrologic model assimilation platform; and (D) constructing a watershed soil moisture remote sensing data assimilation scheme based on a distributed hydrologic model and particle filter assimilation algorithm. In the invention, a novel distributed watershed hydrologic model which is capable of effectively fusing microwave remote sensing information and has a certain physical basis is constructed by utilizing soil water hydrodynamic method and combining a saturation excess runoff principle, a hydrological simulation detection result of the Yi River watershed in a typical semi-arid and semi-humid area shows that the novel distributed watershed hydrologic model has a better daily runoff simulation effect and stable surface soil moisture simulation precision and can be used as a watershed soil moisture remote sensing data assimilation model operator. By using the watershed scale soil moisture remote sensing data assimilation method, a watershed scale soil moisture assimilation data set in temporal and spatial distribution can be acquired effectively.

Description

Basin yardstick soil moisture remotely-sensed data assimilation method

Technical field

The present invention relates to a kind of remotely-sensed data assimilation method, the basin topsoil humidity retrieval information and the specific data assimilation method of saying so more specifically and utilizing remote sensing to obtain.

Background technology

The research of basin yardstick and field yardstick soil water regime helps deepening people and offers help to the understanding of hydrologic process and for the final complex interaction effect of grasping surface water and groundwater resource.In the hydrographic water resource field, soil moisture space distribution accurately can be used as the initial value of the required soil moisture variable of hydrological distribution model flood simulation, improves the accuracy of flood forecasting; In addition, will be through the soil moisture data substitution hydrological model accurately continuously of overcorrect, for the precision that improves long-term hydrologic process simulation and Runoff Forecast and estimate that accurately surface water resources amount has potential value.Therefore invention with the soil moisture variable of basin yardstick as main research object.

The research of basin soil moisture be unable to do without two kinds of means of model and observation.At first, hydrological model is the effective tool that calculates and analyze soil moisture content, and is significant in the relevant time-space process of simulated soil moisture.With the Xinanjiang model is the forecasting model that soil water regime has been set up on the basis, is based on the lump type model of notion, can only provide the average soil humidity information in basin roughly, can't depict the inner detailed space hydrologic process in basin.In recent years, along with the development of computer science and infotech, the distributed basin hydrological model under RS and the GIS technical support becomes the focus of research.The DEM treatment technology of RS and GIS provides hydrometeorological parameter and underlying surface vegetation, soil information and the basin terrain parameter of the required different spatial and temporal resolutions of hydrological distribution model, like parameters such as the gradient, aspect, flow path and borders, basin; Distributed basin hydrological model in conjunction with the foundation of Spatial Information Technology such as RS, GIS, DEM; Can be on littler hydrological simulation unit; The accurately hydrology cyclic process under Simulation of Complex weather conditions and the surface condition, thereby hydrological distribution model can provide the space output of multiple hydrology variable in the basin.This has remarkable difference with the lump type hydrological model of mainly simulating basin water delivering orifice flow through notion parameter and soil moisture variable.Just because of above characteristics; Hydrological distribution model not only can be moved for relevant hydrology ecological Studies such as basin product, sediment transport, nutrients are defeated, Pollutants Diffusion and mankind's activity provide advanced calculating and analog platform to the influence of water cycle etc.; And, important purposes is arranged like pre-warning time aspect with weather forecast Mode Coupling prolongation flood damage at flood forecasting.But; Also there is certain deficiency in hydrological distribution model; For example; Because using, the dynamics scheme that receives laboratory model is still waiting improvement on the basin; How to set up the stronger hydrological distribution model of physical basis, obtaining more, the spatial and temporal pattern of the hydrologic parameters such as soil water content of objective reality is the important directions of distributed model research; In addition; Owing to receive the influence of model error of input data and model self structure error; How certain uncertainty that the analog result of hydrological distribution model also exists quantizes and reduces the uncertainty of hydrological model prediction, also is the forward position focus of hydrological science current research.Utilizing remotely-sensed data to pass through the simulation error accumulation that data assimilation method reduces model, is the very promising method of obtaining high precision and high resolving power soil moisture data.

Secondly, soil moisture can be obtained through observing directly, and the direct monitoring method of soil moisture is broadly divided into two kinds of ground investigations and airborne/satellite remote sensing.Monitoring method based on ground investigation is different with means according to obtaining data mode; Can be divided into earth boring auger fetch earth weight method, oven drying method, neutron appearance method, electric-resistivity method, TDR method again; These methods adopt ocean weather station observation on the space; Sample rate is slow, manpower, material resources and financial resources consumption is big; Be difficult to satisfy the needs of the continuous dynamic monitoring of large tracts of land; Receive the restriction of its spatial sampling density and the influence of atmospheric water energy process random perturbation, the real space general layout of soil moisture also can't accurately be observed.Compare with ground investigation, remote sensing technology has the characteristics of large tracts of land simultaneous observation, ageing, economy, is the effective way that solves soil moisture observation.Along with remote sensing technology develops to high spatial resolution, high spectral resolution direction, the effect that remote sensing technology becomes more and more important performance in the dynamic monitoring of soil moisture large tracts of land.

The research of remote sensing monitoring soil moisture starts from twentieth century end of the sixties; Over nearly 40 years; The method of remote sensing monitoring soil moisture is also constantly being improved and is being upgraded; Computing method based on different remote sensing principles have appearred, as: thermal inertia method, normalized differential vegetation index method, state of temperature index method, temperature vegetation drought index (TVDI) method, microwave remote sensing method etc.Can be used for monitoring soil moisture situation wave band and comprise visible light and near infrared, thermal infrared wave band and microwave region, wherein visible light and near-infrared band data owner will calculate inverting through above-mentioned normalized differential vegetation index method, state of temperature index method, TVDI method and obtain soil moisture; The thermal infrared wave band data can adopt the thermal inertia method to obtain humidity information; Microwave data then mainly concerns, obtains soil moisture based on the method for the microwave radiation mode inverting of physics through the empirical statistics of brightness temperature and soil moisture.Compare with infrared band with visible light, utilize microwave data inverting soil moisture to have solid physical basis, its application is extensive far beyond the data of visible light and infrared band, and inverting reliability and precision are also higher.With the advantage of microwave data inverting soil moisture be mainly reflected in following some; 1) physical characteristics of microwave.The microwave remote sensing soil moisture inverting has solid physical basis: atural object microwave emissivity depends mainly on the specific inductive capacity of object, and soil dielectric constant depends mainly on the moisture of soil.In microwave region; The specific inductive capacity of water is approximately 80; Dry ground only is 3; Variation along with moisture; The emissivity of soil from 0.6 (30% volume soil moisture) of wet soil to changing 0.9 (the 9% volume soil moisture) of dry ground; Have bigger contrast between them, so microwave imagery is very responsive to moisture.2) regional adaptability.In China south, because weather is moistening, cloudy rainy, limited in time obtaining of remote optical sensing data, influenced the completion on time of survey tasks.The high-resolution radar remotely-sensed data source that the soil moisture investigation urgent need of cloudy area of heavy rainfull, China south has round-the-clock earth observation ability.3) because the soil moisture variation mostly occurs at rainy weather, at the time window of this important humidity changing features of rainy weather, microwave data is unique feasible soil moisture investigation method.Microwave remote sensing observation soil moisture also has certain limitation, owing to receive the influence of vegetation and the restriction of microwave wavelength, satellite-borne microwave remote sensing at present can be surveyed the soil moisture that obtains and only be the moisture on top layer, ground (only several centimetres).The advantage that therefore will combine remotely-sensed data and hydrological model, must through the data assimilation system with the moonscope data integration in hydrological distribution model, improve simulation, and the soil profile moisture of physics unanimity be provided the basin water cycle process.

Summary of the invention

1, invents the technical matters that will solve

The purpose of this invention is to provide basin yardstick soil moisture data assimilation method; That be suitable for soil moisture simulation through structure and have a hydrological distribution model that holard process kinetics is described; On the basis of River Basin Hydrology process simulation and checking; It is developed into the hydrological distribution model assimilation platform that can effectively assimilate the remote sensing soil humidity information, can be to obtain the basin soil moisture space-time assimilation data set of higher time precision and spatial resolution.

2, technical scheme

The object of the invention is mainly realized through following steps:

Basin yardstick soil moisture data assimilation method the steps include:

A) data are prepared;

B) make up basin soil moisture assimilation observation operator

C) make up hydrological distribution model assimilation platform;

D) made up basin soil moisture remotely-sensed data assimilation scheme based on hydrological distribution model and particle filter assimilation algorithm.

Steps A) data in are prepared, and comprising:

1) hydrology, weather data;

2) DEM (digital elevation model) data;

3) vegetation, soil parameters storehouse; These underlying surface supplemental characteristics of soil, basin utilization/cover and soil types; And, set up the required soil parameters spatial database of model according to the typical soil cross-sectional data that " the soil will " of location, basin has been collected each subtype of soil of study area;

4) satellite image data comprise initiatively microwave ENVISAT-ASAR remotely-sensed data and MODIS image data;

5) runoff and soil moisture measured data day by day.

Step B) concrete steps of Gou Jianing are following: adopt new ASAR initiatively microwave data and MODIS image data respectively; Inverting basin topsoil humidity; And soil moisture actual observation data in the basin are used to estimate precision test as a result; Provided the assessment result that quantizes; Set up the basin topsoil humidity time series of remote sensing appraising, and made up basin soil moisture assimilation observation operator based on this remote sensing appraising result.

1.MODIS visible spectral remote sensing data inversion basin topsoil humidity

The present invention adopts the surface soil humidity of TVDI method (temperature vegetation drought index) inverting estimation MODIS data 1km spatial resolution.The TVDI index is considered to and surface soil humidity linear dependence by obtaining in the image LST-NDVI feature space, and TVDI is calculated by formula 1 and tries to achieve,

$\mathrm{TVDI}=\frac{\mathrm{LST}-{\mathrm{LST}}_{\min }}{{\mathrm{LST}}_{\max }-{\mathrm{LST}}_{\min }}---\left(1\right)$

LST is a surface temperature; LST _MinBe the minimum surface temperature of identical NDVI value, corresponding is the wet limit of LST-NDVI feature space; LST _MaxBe the maximum surface temperature of identical NDVI value, corresponding is the dried limit of LST-NDVI feature space.

The scope of TVDI value is from 0 to 1, and the corresponding TVDI value of doing on the limit of point is 1, and corresponding TVDI value is 0 on the wet limit.Can't comprise the picture dot space fully owing to do wet limit in the actual computation, the TVDI exponential quantity can exceed this scope.The TVDI value is low more, and expression soil moisture is big more, the closer to wet limit.Wet limit is that corresponding soil moisture is the isoline of field capacity, does the corresponding soil moisture theoretical boundary zero in limit.

2.ASAR active microwave data inverting basin topsoil humidity

Based on improved ERS method for normalizing, the topsoil humidity time series of study area of having utilized the ENVISAT ASAR image inverting obtain, and utilize existing soil moisture observational data that inversion result is verified, provided the quantitative evaluation result.

The ASAR remotely-sensed data is through radiant correction, and geometry correction is after the pre-service such as hot spot filtering; For consistent with the spatial resolution of MODIS data; To be used for the assimilation test, the ASAR image is resampled to the 1000m raster resolution, as the input of soil moisture appraising model.

Among the present invention,, the correlation parameter of appraising model has been done certain correction, obtained being suitable for the model parameter of the ASAR soil moisture estimation in basin, semiarid study area Yihe according to study area characteristic and actual measurement soil humidity data.Detailed ASAR data humidity estimation steps is following:

Scattering coefficient normalization

Different face of land back scattering that cover and incident angle exist the linear dependence relation, and IEM radiation delivery model simulation experiments has also been verified the existence of this relation, and face of land back scattering and incident angle relation can use following formula 2 to express,

σ ⁰＝pθ+q (2)

Here, θ representes incident angle, and p, q cover relevant empirical parameter with the soil, and concrete numerical value can be by finding among the Loew et al..Although the incident angle of each scape ASAR WS pattern image has nothing in common with each other, the zones of different incident angle variation range in the same scape image is also bigger, can through following formula, the backscattering coefficient of all images is normalized to the incident angle θ of reference ₀, θ here ₀Be incident angle value 23 degree of ERS, concrete method for normalizing is following,

${\mathrm{\σ}}_{\mathrm{nor}}^0=\frac{p{\mathrm{\θ}}_0+q}{\mathrm{p\θ}+q}{\mathrm{\σ}}^0---\left(3\right)$

Wherein, σ ⁰Be the original backscattering coefficient of image,

Be the image backscattering coefficient after the normalization.Inverted parameters

Numerous researchs show; The correlationship of backscattering coefficient and face of land specific inductive capacity has following positive correlation; This empirical relationship is mainly drawn by the experience database at 20 observation experiment stations that are located at Germany; And being used to the soil moisture inverting of the radar data of ERS and mesoscale, the precision (RMSE) of inverting of utilizing the topsoil humidity that this method obtains is at 0.04m ³/ m ³To 0.07m ³/ m ³Between.Utilize last one to go on foot the ASAR back scattering that the normalization aftertreatment obtains among the present invention, convolution 4 is inquired into the specific inductive capacity that obtains soil, and obtains the soil moisture value on top layer with this final inverting,

∈ _r＝a+bσ ⁰+c(σ ⁰) ² (4)

∈ wherein _rBe specific inductive capacity; A, b, c are empirical parameter; Relevant with the soil cover type; Owing to receive study area difference, according to different soil cover types, the present invention does certain correction to the parameter of this formula; Its occurrence such as table 1; Farmland accounting example 69% in the study area land use pattern, the meadow accounts for 16%, and all the other ratios are 15%.All the other unclassified types are calculated according to the parameter in farmland.

Table 1: the ASAR remote sensing soil moisture estimation parameter of using among the present invention

Dielectric model

Because specific inductive capacity and soil moisture have correlativity preferably, therefore the specific inductive capacity result who acquires through said method can utilize the soil humidity information that dielectric model is translated into to be needed.

Developed at present the relation that a lot of theories and semiempirical model are used to describe specific inductive capacity and soil moisture under the natural terrain condition.What wherein use always the most in the remote sensing field is the specific inductive capacity model that Hallikainen et al. proposes; Soil dielectric constant and the funtcional relationship between the soil moisture that this model provides are following; And provided specific inductive capacity and soil moisture correlationship (the sand grains proportion 51.5% that obtains by this modeling; Clay proportion is 13.4%)

∈ _r＝(a ₀+a ₁S+a ₂C)+(b ₀+b ₁S+b ₂C)m _v+(c ₀+c ₁S+c ₂C)m _v ² (5)

Wherein, m _vBe calculative volumetric(al) moisture content, S and C are the sand grains of soil and the percentage by weight of clay, and promptly soil texture parameter can obtain a by the basin soil physics parameter library that we set up _i, b _iAnd c _iBe the empirical constant of model, provide by document.

Based on above-mentioned model, utilize basin soil dielectric constant image and the existing soil texture parameter obtained, inquire into the topsoil humidity estimation result who obtains final ASAR remotely-sensed data by grid through iterative calculation method (error threshold is 0.01).

Step C) made up the hydrological model of simplifying that can be used for basin runoff and the simulation of topsoil humidity, and with its platform as the assimilation of basin soil moisture remotely-sensed data.

The modeling thinking:

Model is regarded the basin by the soil cylinder of different-thickness as is formed, and the soil horizon below is saturated underground water.The thickness of soil cylinder is given by the actual soil investigation data in basin, mainly confirms according to the representative section investigation data on the spot of different soils type, so the cylinder thickness of different soils type and inequality.The special heterogeneity of this characteristics decision basin soil water storage capacity; Runoff yield under saturated storage and watershed concentration through lattice point soil cylinder are calculated; Simulation obtains basin general export flow; And, describe soil water content and dynamically migration on the space scale of basin in detail based on vertical unsaturated soil hydrodynamics method and soil physics parameter library.

The soil horizon in root district mainly takes place and is decided by in the production process of soil moisture transmission and rainwash, so the description of this soil mass is most important.In the model, be to be no more than 10 layers soil horizon with the soil cylinder in basin according to identical layering criteria, maximum depth of soil is 3.43 meters, thereby the soil cylinder on the lattice point is because the difference of thickness of soil has the different soil numbers of plies.The soil texture parameter of each layer of soil comes from the soil parameters storehouse that the soil investigation data is set up, and equally also is that the investigation data on the spot according to the representative section of different soils type obtains.On the basis that obtains the soil layering parameter; The quantity of precipitation on utilization to ground and runoff yield difference obtain the average milliosmolarity down of soil surface; Based on unsaturated soil hydrodynamic force equation; Simulate vertical moisture transmission of moisture state and the interlayer of each layer soil, model has also been considered the influence to soil moisture of vegetation root system water sorption and evaporation of soil moisture simultaneously.

The concrete construction framework of model as shown in Figure 1; Wherein k1...k3 representes the soil hydraulic conductivity; θ 1... θ 3 expression soil water contents; Φ 1... Φ 3 expression soil water potentials; Q1...Q3 representes soil water flow; Dz1...dz3 representes soil horizon thickness, and the hydrologic process that provides among Fig. 1 all is contained among the new model.

Soil moisture is the key variables in the River Basin Hydrology process, and studying in great detail and understanding of its Spatial Distribution Pattern is significant for applied research fields such as land gas exchange, environment and the ecological hydrology and subject.Hydrological model is in the past mainly paid close attention to the calculation of the calculated example of rainfall-runoff process and water resources quantity population equilibrium like the outlet run-off; Shortage is carried out detailed description to the moisture movement of soil profile; Therefore hydrological model often only calculates the soil moisture degree of saturation of generalization in the basin; And the average moisture content information of one deck of generalization or multilayer soil; And can't provide really can with the quantification soil humidity information of measured data contrast, like the soil water content of top layer 10cm.

As the model operator of basin soil moisture remote sensing assimilation system, must calculate and the corresponding top layer of remote sensing soil moisture estimation result 5cm-10cm soil moisture simultaneously, could further realize soil moisture assimilation test.Comprehensive above 2 points, topsoil humidity is calculated and the hydrological model of assimilation therefore to be necessary new being suitable for of research and establishment.Table 2 has been listed the physical process description of new model that the present invention proposes and above common and up-to-date hydrological model.

Table 2: the hydrologic process guide look of different hydrology modelings

This model simplification River Basin Hydrology runoff yield simulation, on this basis vertical hydrologic process is carried out comparatively detailed modeling, mainly described soil moisture transmission, the holard with phreatic alternately and the moisture extraction effect of root system.Model comprises the soil type map and the typical soil sectional parameter data in basin through local soil investigation historical summary storehouse, sets up a soil space parameter library, and has got in touch the dynamic calculation of soil moisture and face of land runoff yield based on this database.Below each step provide hydrologic process and the description thereof relevant of simulating in the model with the soil moisture process.Table 3 has provided all hydrology variablees that model can simulate and required meteorology drives variable.

Table 3: in this model the hydrologic process argument table

1) canopy is held back

The present invention adopts improved concept nature Aston exponential model to describe canopy and holds back.Crown canopy interception computing formula is following:

$I_{\mathrm{int}}=C_m\·(1-e^{-{\mathrm{kPC}}_m})+e_w\·T---\left(6\right)$

C wherein _mBe the canopy interception capacity, k is the rainfall interception coefficient, and P is rainfall amount (mm), and T is rainfall duration (hr), e _wBe the crown canopy rate of evaporation, so e _wT representes tree crown evaporation and trunk evaporation sum.

Canopy interception capacity C _mThe method that adopts Von Hoyningen-Huene to be proposed is calculated:

C _m＝C _p·(0.935+0.498·LAI-0.00575·LAI ²) (7)

C wherein _pBe vegetation cover degree, LAI is the leaf area index day by day of lattice point cell-average

The rainfall interception coefficient k is the function of vegetation leaf area index, is provided by following formula 8:

k＝0.046·LAI (8)

Crown canopy rate of evaporation e _wRelevant with factors such as the degree of drying at rainfall initial stage and wind speed; When no measured data is verified; Usually get following empirical parameter: the subtropical zone gets 0.031～0.047, and the temperate zone is moistening, half humid region gets 0.063～0.093, and temperate zone arid, semiarid zone get 0.170～0.200.

2) face of land runoff yield

Study area of the present invention is humid region and semi-arid-semi-humid climate district, so model employing runoff yield under saturated storage mechanism calculating rainwash, in the hydrological simulation process, has only when the quantity of precipitation that drops down onto soil surface is dark greater than the soil lack of water, just can produce rainwash.

In model, for each lattice point unit independently, water balance is calculated by 9 formulas down and is provided:

P-Es-Et-G＝Δw+Δc+R (9)

Wherein, P is precipitation (comprising rainfall and snowfall), and Es is an evaporation of soil moisture, and Et is a plant transpiration, and G is a groundwater discharge, and Δ w is the increment of the holard, and Δ c is that canopy is held back, and R is a rainwash.

3) soil moisture transmission

In this model, Richard ' the s equation of having considered gravity and absorption affinity is used to describe vertical unsaturated soil moisture movement of soil profile, numerical solution adopts method of finite difference.Because the gradient in study area basin is less, the cross flow of (1km) holard is less relatively under the space scale of model running, so ignores the influence that lateral current is simulated humidity in the soil in the model.Model is divided into 6～10 layers through unified branch layered scheme (as shown in table 4) with the soil profile of unit lattice point, and the thickness of soil that the surface is 3 layers (near 10cm) is complementary with remote sensing appraising result's significant depth.Soil is counted as homogeneous in the graticule mesh of modeling, also regards consistent as in the inner soil property of sub-soil horizon, and the model hypothesis rock face degree of depth becomes positive correlation with the root region soil degree of depth on the lattice point in addition.

Table 4: model soil layering scheme

The soil layering scheme that table 4 provides is with reference to the soil profile lift height computing formula of the land face procedure schema CLM that is suitable for the simulation of topsoil humidity.In each independent soil horizon, the distribution of property parameters such as this model hypothesis soil texture is a homogeneous.Its concrete soil parameters obtains from the basin soil space parameter library of soil profile survey data and foundation thus on the spot.

Relate to Richard ' s equation in the soil moisture transmission, be based on following two equation inferences, unsaturated soil water flux law Buckingham-Darcy equation (formula 10) and continuity equation (formula 11):

The Buckingham-Darcy equation:

$q=-K\left(h\right)(\frac{\∂h}{\∂z}+1)---\left(10\right)$

Wherein q is unsaturated water flux density, and h is a holard matric potential, and K (h) is a hydraulic conductivity, and z is the degree of depth of vertical direction.

Continuity equation:

$\frac{\∂\mathrm{\θ}}{\∂t}+\frac{\∂q}{\∂z}=0---\left(11\right)$

Based on the hypothesis of vertical flow, the q of formula 11 is replaced by formula 10, tries to achieve to separate to be one dimension Richard ' s equation expression formula:

$\frac{\∂\mathrm{\θ}}{\∂t}=\frac{\∂}{\∂z}\left(K\left(h\right)\frac{\∂h}{\∂z}\right)+\frac{\∂K\left(h\right)}{\∂z}---\left(12\right)$

Wherein t is the time, and θ is the soil volumetric(al) moisture content.

Root soil moisture extracts, and needs in the holard cubage to consider the root water uptake effect, and this model adds it in Richard ' s equation to as the S item.The present invention uses following method to consider the root tissue water status absorption of different soils layer:

$\frac{\∂\mathrm{\θ}}{\∂t}=\frac{\∂}{\∂z}\left(K\left(h\right)\frac{\∂h}{\∂z}\right)+\frac{\∂K\left(h\right)}{\∂z}-S---\left(13\right)$

Wherein S representes soil root water uptake or evaporation of soil moisture, is described by following formula 14

$S_i=E\frac{r_i{w}_i}{\mathrm{\Σ}{r}_i{w}_i}---\left(14\right)$

Wherein i is the sequence number (i>=2) of simulated soil layer, for ground floor soil (i=1), S _iBe topsoil surface evaporation, E _VegBe the transpiration of plant, r _iBe Root Distribution parameter in the soil horizon, relevant with vegetation pattern, w _iBe holard saturation coefficient, when humidity is in wilting point, equal 0, complete is 1 when saturated, is provided by following formula 15:

Wherein

Be the flow of water of wilting point, equal-1.5 * 10 ⁵,

Be the saturated soil flow of water,

It is holard matric potential.

4) potential evapotranspiration is sent out

The present invention adopts two kinds of methods to calculate potential evapotranspiration and sends out, and when there was the evapotranspiration research station in the basin, potential evapotranspiration was sent out calculating and obtained by the observation data interpolation; When the basin does not have the observation of evapotranspiration, utilize meteorological measuring to calculate according to the Penman-Monteith formula, Penman-Monteith formula 16 is provided by following formula:

$E=\frac{1}{\mathrm{\λ}}\left[\frac{\mathrm{\ΔA}+{\mathrm{\ρ}}_a{c}_{p}D/r_a}{\mathrm{\Δ}+\mathrm{\γ}(1+r_s/r_a)}\right]---\left(16\right)$

Wherein E is an evaporation capacity, and λ is a steam latent heat, and Δ is the slope of saturation vapour pressure temperature curve, and A is a net radiance, ρ _aBe the density of soft air, c _pBe the air latent heat of vaporization under the normal pressure, D is the difference of vapour pressure and saturation vapour pressure, r _aBe the aerodynamic drag coefficient, γ is the psychrometer constant, r _sIt is the resistance coefficient on crop blade face.

5) plant transpiration

In the practical application of plant transpiration, the use of empirical method and accurate physical method is more, because the simplicity and the validity of Richtie model have been used improved Ritchie model among the present invention.Transpiration rate calculates as follows:

$E_{\mathrm{veg}}=E\frac{\mathrm{LAI}}{3}\left(\frac{1+p}{2}\right),0<\mathrm{LAI}\&le;3.0---\left(17\right)$

$E_{\mathrm{veg}}=E\left(\frac{1+p}{2}\right),\mathrm{LAI}>3.0---\left(18\right)$

E wherein _VegBe the transpiration rate of plant, E is the potential evapotranspiration amount of sending out, and p is the vegetation coverage on the lattice point, and LAI is the leaf area index day by day on the lattice point.

Vegetation coverage p calculates through leaf area index LAI day by day, and formula 19 is following:

p＝1-e ^-α·LAI (19)

In the formula: α is the extinction coefficient of face of land solar radiation, and default value is 0.45.

6) evaporation of soil moisture

In model, adopted a semi-empirical approach based on soil physical properties, considered two Main Stage of soil water evaporation respectively, promptly supplied with abundance when soil moisture, evaporation capacity is higher and relatively stable; When the soil moisture reduction, the evaporation of soil table is along with soil water content is non-linear decline:

$E_{\mathrm{soi}}=E\frac{\mathrm{dz}}{\mathrm{dze}+\exp (2.374-0.00713\&CenterDot;\mathrm{dz})}\exp \left(\frac{2.5(\mathrm{\&theta;}-f_c)}{f_c-w_p}\right),\mathrm{\&theta;}<f_c---\left(20\right)$

$E_{\mathrm{soi}}=E\frac{\mathrm{da}}{\mathrm{dz}+\exp (2.374-0.00713\&CenterDot;\mathrm{dz})},\mathrm{\&theta;}\&GreaterEqual;f_c---\left(21\right)$

E wherein _SoiBe the evaporation of soil moisture amount, E is that potential evapotranspiration is sent out, and dz is the thickness (mm) of topsoil, and θ is the moisture (volumetric water content) of topsoil, f _cBe the field capacity (volumetric water content) of topsoil, w _pIt is the wilting coefficient (volumetric water content) of topsoil.

7) underground water and unsaturated soil water are mutual

In order to consider the variation of soil profile water cut more accurately; Invention combines to simulate the soil water content and the flow of water of the root district deep layer that obtains; Underground water form general model SIMGM based on new has done simple description to this process, the variation of and underground water table mutual with simulated soil water and underground water.For in the simulation of basin, considering reciprocation; Model has mainly been made following three hypothesis; The one, suppose to exist in the water-bearing zone under the deep soil saturated groundwater reservoir with certain depth of water; The 2nd, the hydraulic conductivity between soil horizon and the groundwater reservoir is along with degree of depth increase is index decreased, and the 3rd, the groundwater reservoir on the space, basin is identical deeply.Phreatic supply is calculated and is given by the following formula:

Wherein Q is the increment of groundwater (mm/s), when flowing to underground reservoir under the current direction on the occasion of, otherwise be negative value,

Be the degree of depth of underground water table,

Be subsoil water-based gesture, z _BotBe the node degree of depth of subsoil, K _aBe the hydraulic conductivity in water-bearing zone, provide by following formula:

$K_a=\frac{{\&Integral;}_{z_{\mathrm{bot}}}^{z_{\&dtri;}}{k}_{\mathrm{bot}}{e}^{-f(z-z_{\mathrm{bot}})}\mathrm{dz}}{z_{\&dtri;}-z_{\mathrm{bot}}}=\frac{k_{\mathrm{bot}}(1-e^{-f(z_{\&dtri;}-z_{\mathrm{bot}})})}{f(z_{\&dtri;}-z_{\mathrm{bot}})}---\left(23\right)$

Wherein f is the time-delay factor, is used for the calculating of hydraulic conductivity.

8) conflux in the face of land

Confluxing in the face of land in this model, what use is the period unit line method of confluxing, and the unit line parameter obtains through the parameter calibration of model, unit line be input as the mean value that basin grid runoff yield calculates.Unit line is concrete to be defined as follows, on given basin, and the unit effective precipitation that the uniform rainfall of spatial and temporal distributions produces in the unit period, the direct runoff hydrograph on the outlet formed ground of section, basin, note is made UH.Tp when the nemaline index of control unit has flood peak to stagnate, unit line crest discharge qp and unit line last T, are called the unit line three elements, like Fig. 2.The relation of outlet run-off Q and UH ordinate value is following,

$Q_t=\underset{j=k_1}{\overset{k_2}{\mathrm{\&Sigma;}}}{r}_j{q}_{t-j+1}---\left(24\right)$

In the formula, r _jBe the period effective precipitation, use corresponding face of land runoff yield dark (mm) in the Model Calculation, Q _tBe the basin outlet section rainwash flow (m of period Mo ³/ s), q is unit line period later and decadent stage of a school of thought amount, the i.e. ordinate value of UH; Be the parameter of confluxing that needs calibration in the model, t is the flow path surface sequential, t=1; 2 ..., m+n-1; Hop count when wherein m is net rainfall, i.e. hop count during the simulation of model, hop count when n is unit line; K is the bound of accumulation period, and its value is by m, and the n value determines; Computing method are following

${\mathrm{<figure-callout\; id="1"\; label="k"\; filenames="HDA0000101114300000041.png,HDA0000101114300000051.png"\; state="\{\{state\}\}">k</figure-callout>}}_1=\left\{\begin{array}{cc}1& t<n\\ t-n+1& t\&GreaterEqual;n\end{array}\right.---\left(25\right)$

$k_2=\left\{\begin{array}{cc}t& t<m\\ m& t\&GreaterEqual;m\end{array}\right.---\left(26\right)$

The segment length can appoint and gets during unit in the unit line, like 1h, 2h, 6h etc.Because this model is mainly used in day by day the footpath flow field simulation of yardstick, therefore the time segment length of the unit line that uses is 24h.Between 4～8, the n value is set at 7 to hop count n value in the present invention usually during unit line.

9) underground confluxing

Ground water movement is the hydraulic problem of seepage flow, and Boussinesq has set up the differential equation of the one dimension tangential movement in the saturated zone according to continuity equation and Darcy's law, and linearization has obtained diffusion equation to equation, its separate try to achieve into:

Q _t＝Q ₀e ^-αt (27)

In the formula, Q ₀Be the run in depth initial value, α is phreatic coefficient of extinction, is the empirical value of fixing, and is made as 0.984 in the model, Q _tFor t run in depth constantly effluents, try to achieve through the flow of adjacent period.

Therefore model adopts following formula to simulate effluenting of run in depth, and combines to go up that groundwater reservoir effluents on the lattice point that joint calculates, with the input parameter of its basin mean value as underground runoff concentration calculation.

Step D) structure mainly comprises model operator, observation operator, drives weather data and assimilates four parts of algorithm based on the basin soil moisture remotely-sensed data assimilation scheme of hydrological distribution model and particle filter assimilation algorithm.

The model operator of this assimilation scheme is the new hydrological distribution model that makes up among the step C; The observation operator is a soil moisture of using ASAR and MODIS data estimation to obtain among the step B; The assimilation algorithm is particle filter algorithm (PF), and meteorological driving data comprises precipitation, evaporation etc.The observation operator is the soil moisture that ASAR and MODIS data estimation obtain, and the assimilation algorithm is a particle filter algorithm, and meteorological driving data comprises precipitation, evaporation etc., sees Fig. 3 for details.

Particle filter assimilation algorithm

Particle filter is the algorithm that is used for non-linear non-gauss' condition spatial model optimal estimation problem, and the mathematical description of its concretism is following: for stochastic process stably, the posterior probability density of etching system is p (x when supposing k-1 _K-1| z _1:k-1), according to certain principles of selected n random sample point, after k obtained observation information constantly, through state and time renewal process, the posterior probability density of n particle can be approximately p (x _k| z _1:k), reach the effect of optimum Bayesian Estimation.The core of PF algorithm is based on the Monte Carlo, and (Monte Carlo methods, MC) method are similar to the state probability density function through one group of population of in state space, propagating, obtain the minimum variance of state, so also title order MC method of PF.Although just true distribute a kind of approximate of the probability distribution in the algorithm, because imparametrization, the posterior probability that PF can more accurately express based on observed quantity and controlled quentity controlled variable distributes.For linear gauss' condition spatial model, the optimal filtering method is Kalman filtering.

Because River Basin Hydrology system and soil water movement have high complexity and nonlinear characteristic; Therefore must use the PF algorithm based on the soil moisture assimilation research of hydrological distribution model; To guarantee the rationality of assimilation test, this also is that the present invention adopts the reason of particle filter as the assimilation operator.

The state space equation of dynamic system can be expressed as

x _k＝f _k(x _k-1，v _k-1) (28)

z _k＝h _k(x _k，u _k) (29)

x _kThe expression system state, the topsoil humidity of simulating for hydrological model in this research, z _kExpression observation is the topsoil humidity that is obtained by MODIS and the estimation of ASAR remotely-sensed data, u _k, v _kBe independently observation noise and system noise.Suppose x _kObey the single order Markov process, k-1 posteriority distribution function p (x constantly _K-1| z _1:k-1) known, the model state PDF of k can try to achieve through following formula constantly,

p(x _k|z _1:k-1)＝∫p(x _k|x _k-1)p(x _k-1|z _1:k-1)dx _k-1 (30)

P (x _k| x _K-1) by system equation 30 and known system noise v _kDefine, at k observed quantity z constantly _kBe after the remote sensing appraising result obtains, upgrade through the state of Bayes's equation to system,

$p\left(x_k\right|z_{1:k})=\frac{p\left(z_k\right|x_k)p\left(x_k\right|z_{1:k-1})}{p\left(z_k\right|z_{1:k-1})}---\left(31\right)$

Normaliztion constant p (z wherein _k| z _1:k-1) can provide by following formula,

p(z _k|z _1:k-1)＝∫p(z _k|x _k)p(x _k|z _1:k-1)dx _k (32)

P (z in the formula _k| x _k) by observation equation 29 and known observation noise u _kDefinition.Equation 30; 31 constitute optimum Bayes separates; But its analytic solution are only set up limited model; The estimated accuracy of approach methods such as EKF, GMF and IMM is limited; Particle filter is based on Monte Carlo thought; (Importance Sampling is IS) from significance distribution function q (x for the method for employing importance sampling _0:k| z _1:k-1) middle independent draws sample

The probability distribution function of model state (Posterior Distribution Function, PDF) approach into

$p\left(x_{0:k}\right|z_{1:k})=\underset{i=1}{\overset{N_s}{\mathrm{\&Sigma;}}}{w}_k^i\mathrm{\&delta;}(x_{0:k}-x_{0:k}^i)---\left(33\right)$

Where the weighting to be normalized, that satisfy

is recursive estimation, the importance of the function selected

q(x _0:k|z _1:k)＝q(x _0:k-1|z _1:k-1)q(x _k|x _0:k-1，z _1:k) (34)

Therefore from q (x _k| x _0:k-1, z _1:k) middle sample drawn

Its weight

Computing method are following, and among the present invention, right of formula 35 promptly is modeling and remote sensing appraising soil moisture result's error ratio,

$w_k^i=w_{k-1}^i\frac{p\left(z_k\right|x_k^i)p\left(x_k^i\right|x_{k-1}^i)}{q\left(x_k^i\right|x_{k-1}^i,z_k)}---\left(35\right)$

Above formula has been formed order IS (Sequential IS; SIS) method; To the power degradation phenomena that occurs after the iteration; Gordon proposes the sample resampling; The high particle of breeding weights of importance; Suppress degradation phenomena thereby eliminate the low particle of weight, the most frequently used method for resampling is SIR (Sampling Importance Resampling).General particle filter just is made up of SIS and resampling method.

The calculation process of particle filter algorithm is following, and 1. the soil moisture that obtains of the simulation of pair model is carried out disturbance, produces population, promptly from q (x _k| x _0:k-1, z _1:k) in randomly draw N limited sample, N is a particle number; 2. the obtaining constantly of remotely-sensed data, the soil moisture result who uses remote sensing appraising to obtain, the weight of calculating corresponding particle based on formula 37; 3. the weight of particle being done normalization handles; 4. resample, utilize importance sampling SIR method that the particle weight is carried out resampling; 5. after the weight sampling, 33 pairs of model states of convolution carry out last update calculation, the soil moisture result after obtaining assimilating.

3, beneficial effect

The present invention compared with prior art, its beneficial effect mainly embodies as follows:

(1) to this weak link of holard simulation in the existing hydrological model, solved hydrological model and can't provide topsoil humidity quantitative information, and the problem that can't directly assimilate the model operator of research as remote sensing soil moisture.Utilize holard dynamic method to combine the runoff yield under saturated storage principle; Make up the new microwave remote sensing information that can effectively merge, had the distributed basin hydrological model of certain physical basis; Hydrological simulation check through basin, typical semiarid subhumid Yihe; The result shows; Day by day the runoff simulate effect better, topsoil humidity simulation precision has stability, can be used as the model operator of basin soil moisture remotely-sensed data assimilation.

(2) set up basin yardstick soil moisture assimilation scheme, utilized the soil moisture of earth observation satellite TERRA-MODIS and ENVISAT-ASAR to estimate that the result has carried out four-dimensional data assimilation to the topsoil humidity variables in the hydrological model based on hydrological distribution model.Research through in continuous hydrological simulation in basin, semiarid subhumid Yihe and particle filter assimilation test shows; The assimilation result has effectively merged the spatial framework information of remote sensing humidity; And having improved the soil moisture simulation precision of hydrological model, this scheme can effectively be used for the obtaining of basin yardstick soil moisture assimilation data set of spatial and temporal distributions.

The present invention has made up a kind of soil moisture simulation of being suitable for; Has the hydrological distribution model that soil water process kinetics is described; On the basis of River Basin Hydrology process simulation and checking; It is developed into the hydrological distribution model assimilation platform that can effectively assimilate the remote sensing soil humidity information; Estimate that with ASAR and MODIS soil moisture the result is as the observation operator; In carrying out soil moisture numerical simulation process; Adopt particle filter assimilation ASAR or MODIS estimation to obtain the soil moisture result; And the result after will assimilating returns distributed model continuation calculating generation; Correct the model running track; Iterative cycles obtains basin soil moisture space-time assimilation data, finally improves the soil moisture simulation precision.

Test findings shows in the checking district; The new hydrological distribution model that the present invention makes up is satisfactory to the simulate effect of runoff and soil moisture day by day; And has an advantages of higher stability; And; Soil moisture result after the assimilation of employing particle filter algorithm can better correct the trend of underestimating of modeling value, soil moisture is changed more be tending towards reasonable.Thereby accurately and at large obtain the basin soil water space-time data set of the harmonious property of physical process and power mechanism, for statistics and assessment, valley environment analysis and the flood forecasting of the saturated situation of basin yardstick soil moisture provides reliable basis.

Figure of description

Fig. 1 is the structure thinking figure of hydrology assimilation model;

Fig. 2 is a unit line key element synoptic diagram;

Fig. 3 is the basin topsoil moisture assimilation scheme based on hydrological distribution model;

Fig. 4 makes up flow process for the observation operator of remote sensing soil moisture assimilation;

Fig. 5 is inverting of MODIS soil moisture and actual measurement comparison diagram;

The actual measurement checking of Fig. 6 ASAR image soil moisture estimation (2006, DOY148)

Fig. 7 is that rate is regularly and the checking phase runoff analogue value and measured value comparison diagram day by day in the basin, Yihe for model;

Fig. 8 is the contrast (fixed observer point) of soil moisture (10cm) analog result and measured data;

Fig. 9 is the contrast (on-fixed observation station) of soil moisture (10cm) analog result and measured data;

Figure 10 is the contrast (fixed observer point) of soil moisture (20cm) analog result and measured data;

Figure 11 is the contrast (on-fixed observation station) of soil moisture (20cm) analog result and measured data;

Figure 12 is a basin yardstick soil moisture remotely-sensed data assimilation process flow diagram;

Figure 13 is the ideal test assimilation result of different lattice points in the study area (is example with actual measurement soil moisture observation station);

Figure 14 is simulation of topsoil humidity and assimilation result's actual measurement checking (DOY:2006148);

Figure 15 is continuously assimilation result and the comparison diagram of surveying humidity.

Embodiment

Through following examples the present invention is further described:

Embodiment

1. study area overview

The basin, Yihe is positioned at south, Shandong Province, North China, and its catchment area is 2413km2, basin mean annual precipitation 830mm, and Rainfall in Flood Seasons 616mm accounts for 75% of annual precipitation, belongs to half moistening semiarid region.

2. study area basic data

A. the hydrology, weather data

Collected basin, Yihe 2001～2003,2006,2007 years observation of quantity of precipitation day by day and the evaporation from water surface observational data of totally 10 weather stations, basin flow verification website is hydrometric station, Linyi.

The B.DEM data

Compiled the digital elevation that simulation and assimilation research institute need: the DEM of the 90m resolution that provides by the SRTM plan, sample respectively to 1000m and 240m spatial resolution, adopt utm projection, the WGS84 coordinate system was with on the 50th minute.

C. soil, vegetation data

The pedological map in basin derives from the basin digitizing pedological map of 1: 400 ten thousand engineer's scale and the basin digitizing vegetation map of 1: 100 ten thousand engineer's scale; Utilize ARCGIS software that basin vegetation and the soil types digital vector figure that collects carried out the rasterizing processing, generate the soil types grid map and the vegetation pattern grid map of the 1 kilometer spatial resolution in basin respectively.In order to be complementary the projection that The data is identical with DEM, coordinate system and spatial resolution with dem data.

D. satellite image data

Basin, Yihe in 2006 all retrievable MODIS L1B image and ENVISAT-ASAR images have been collected.In order to set up the soil moisture assimilation system of correct lattice point coupling, remotely-sensed data has all adopted the spatial resolution identical with DEM (wherein the ASAR data space being resampled to 1000m resolution), projection and coordinate system.

E. runoff and soil moisture measured data

Collected the runoff observation data at hydrometric station, Linyi on the basin, Yihe; Its latitude and longitude coordinates is 35 ° of 01 ' N; 118 ° of 24 ' E; The runoff measured data of obtaining comprised hydrometric station, Linyi 2001～2003 years; The actual measurement in 2006～2007 years is the runoff data day by day; Wherein 2001～2003 annual datas are used for the parameter calibration of hydrological distribution model, and 2006～2007 annual datas are used for modelling verification;

In the study area basin; Having 8 soil moisture actual measurement websites can provide data to be used for checking; Mainly chosen Pingyi wheatland wherein in the instance as main fixedly location observation station; The normal observation time of this website is distributed in 8,18 and 28 days of upper, middle and lower ten days every month, and the soil moisture observation procedure adopts the oven dry weight method; Research station, on-fixed location point has been chosen Plain, Pingyi website, and observation time is 28 days of every month.Each observes the soil moisture data of website be the 10cm in 2006 and the 20cm degree of depth.

3. soil moisture inversion result and checking

The soil moisture inverting comprises MODIS visible spectral remote sensing inverting and active microwave remote sensing inverting, and concrete grammar and formula are seen the B step, and the assimilation operator makes up flow process and sees Fig. 4.

A.MODIS visible spectral remote sensing inverting:

The present invention is based on the surface soil humidity of TVDI method inverting estimation MODIS data 1km spatial resolution.At first, set up the LST-NDVI space and the relevant limit and the wet frontier juncture system of doing in four season of basin respectively, be used for TVDI and calculate, the average gradient on the wet limit in four LST-NDVI spaces and dried limit and intercept as parameter, have been made up relational expression 2 calculating TVDI coefficients.

$\mathrm{TVDI}=\frac{T_s-(14.898\&times;\mathrm{NDVI}+270.75)}{(-39.550\&times;\mathrm{NDVI}+326.271)-(14.898\&times;\mathrm{NDVI}+270.75)}---\left(36\right)$

Next utilizes the linear relationship that the distributed sites observation data has been set up TVDI and surveyed top layer humidity in the basin; The present invention adopts the correlationship of formula 3 approximate match TVDI and soil moisture (SM), and relation combines original TVDI coefficient estimate that is drawn by the MODIS image to obtain basin topsoil humidity in view of the above.Among the estimation result, because actual TVDI calculated value can exceed 0,1 scope, soil moisture is distributed between 10% to 45% (percent by volume).

SM＝-23.172×TVDI+29.551 (37)

Ultimate analysis MODIS estimation obtain the soil moisture spatial framework, and the continuous soil moisture in 2006 that inverting obtains are estimated that results verify through the measured data of canned paragraph observation website.

Fig. 5 has provided the contrast scatter diagram of MODIS inversion result and website actual observation data, and data come from different MODIS inversion results constantly and corresponding measured value constantly, and SM representes soil moisture among the figure.

Fig. 5 scatter diagram is the result show, MODIS estimation result and measured data correlativity are higher, and root-mean-square error is 0.049m ³/ m ³, precision is higher.Still keep similar correlationship between TVDI and topsoil humidity, so the effectively moistening degree and the dynamic change thereof of image study district topsoil of TVDI method.Utilize the advantage of MODIS data estimation soil moisture be its have higher spatial resolution (～1km); And the time that returns to is short; Help disclosing the dynamic rule of soil moisture; And that shortcoming is that image receives weather effect is bigger; Especially seven, rainy season of August, the soil moisture change information between flush period can't obtain.And through combining microwave data can effectively improve this situation.

B. active microwave remote sensing inverting:

Through the detailed ASAR data humidity estimation steps among the step B in the summary of the invention, the 14 scape study area ASAR WS image invertings that utilized obtain in 2006 have obtained basin topsoil humidity result, have drawn the inversion result of DOY148 day in 2006.DOY148 day is the fixed observer date of all websites, and the 10cm soil moisture data of utilizing a plurality of eyeballs in the study area are verified the estimation result of ASAR image.Measured data estimates that with corresponding ASAR image result's diffusing some distribution is provided by Fig. 6, finds out that easily remote sensing appraising result and measured value have correlativity preferably, and its estimation precision (RMSE) is 0.063m ³/ m ³, mean deviation (MAE) is 0.053m ³/ m ³This estimation result has reached ripe ERS soil moisture estimation precision scope (0.03m ³/ m ³-0.07m ³/ m ³).

4. River Basin Hydrology process simulation and checking

Detailed step and formula through detailed generalization of hydrologic process among the step C in the summary of the invention; Set up hydrological distribution model assimilation platform; And being applied to the hydrologic process simulation in basin, Yihe, the actual observation of using outlet hydrometric station, basin data on flows is day by day tested to the simulation precision of model and is estimated; For the analog capability of the basin topsoil humidity of verification model, a plurality of website soil moisture observational datas in 2006 that utilize the basin, Yihe to obtain have carried out check analysis to the soil moisture analog result of the 10cm and the 20cm degree of depth.

A. footpath flow field simulation and evaluation day by day:

Selecting above basin, Linyi, Yihe to carry out the simulation test of runoff day by day of model, wherein 2001～2003 years was that the rate of modeling is regular, and 2006～2007 years as the modelling verification phase.

Fig. 7 provided respectively the model rate regularly with simulation and the measured discharge graph of checking phase.Can find out that the measured discharge process is totally coincide better with simulation runoff process day by day.

B. soil moisture simulation and precision evaluation:

Basin topsoil humidity analog capability for verification model; Be utilized in a plurality of website soil moisture observational datas in 2006 that the basin, Yihe obtains, the 10cm that modeling in 2006 is obtained and the basin soil moisture result of the 20cm degree of depth have carried out check analysis.

Simulation of topsoil (10cm) humidity and precision evaluation:

Utilize the data of actual observation point in the basin, Yihe in 2006 that the modeling value is verified.

Shown of the time series contrast of observation station 0～10cm deep soil humidity measured value among Fig. 8 and Fig. 9 with the modeling value of corresponding grid.At the canned paragraph point is that the observation time in wheatland location is at interval shorter, therefore more observation information is provided, and the observation in all the other on-fixed locations is carried out on the 28th in every month, to also separate the mapping comparison in the both of these case literary composition.

Find out from Fig. 8 that easily the topsoil humidity analogue value has correlativity preferably with the actual measurement humidity data, especially in the fixing observation station in location.Pingyi fixedly variation tendency and the numerical values recited of the location analogue value is all very approaching with the sequence of observations; The vertical simulation that shows soil moisture model has certain reliability; Satisfy condition down at data precision; Can reach higher simulation precision; This has also guaranteed basic accuracy and reliability for assimilation test generation has the conforming continuous soil moisture data set of physics, and the Pingyi wheatland fixedly square error in on-fixed location, peace city Plain, location is respectively 0.051 and 0.080.All in all; With respect to measured data; There is the trend of underestimating to a certain degree in the modeling value; This mainly is owing to cause the too fast institute of edaphic dehydration to cause to the extraction item calculating of soil horizon moisture is excessive; And in the model root system relevant with transpiration to absorb the water function simulation bigger than normal be primary factor, so the calculating that plant transpiration amount and layering root tissue water status extract still needs further to improve.It is comparatively serious to underestimate trend at the moisture in some soil middle layer, shows that section Root Distribution parameter is unreasonable on the lattice point.The remote sensing soil humidity information will help to improve the situation that simulation is underestimated through assimilating independently.

Simulation of upper layer of soil (20cm) humidity and precision evaluation:

In order further to provide the evaluation of model topsoil humidity analog result, also utilize the 10cm-20cm deep soil humidity data of interior all observation stations of study area in 2006 to contrast among Figure 10 and Figure 11 with the modeling value time series of corresponding grid.

Similar with the edaphic analog result of 10cm, the edaphic humidity analogue value of 20cm is compared with measured data, has shown higher similarity generally.Fixedly the trend fitting of the actual observation data in location and the analogue value is better, though the trend of certain estimation on the low side is arranged on the whole, this layer soil moisture underestimate reason and 10cm soil horizon basically identical.For the on-fixed location, simulation value and observation data are comparatively approaching, and the Pingyi wheatland fixedly square error in on-fixed location, peace city Plain, location is respectively 0.041 and 0.077.

Generally speaking, more than figure has provided the contrast of the analog result and the measured data of soil moisture, has verified that model is used to simulate the feasibility of basin topsoil moisture.From all checking results, there is correlativity preferably between modeling value and the observed reading, for basin soil moisture assimilation research provides the basic accuracy assurance.

5. basin yardstick soil moisture assimilation numerical experimentation and precision evaluation:

The topsoil humidity that the remote sensing observations data inversion is obtained utilizes particle filter to carry out continuous assimilation test as desirable data, and the checking particle filter is used for the feasibility of complex nonlinear hydrological distribution model remotely-sensed data assimilation.With the basin, Yihe is study area, and modeling and assimilation period are 2006, and virtual space resolution is 1000m.In the test, the soil assimilation degree of depth is top layer 10cm, and observation data is the humidity retrieval result of ASAR and MODIS, data, and the humidity simulated time step-length of model is 30 minutes.Model only assimilated calculating and renewal on the same day that remotely-sensed data can be obtained, the soil moisture that obtains after the assimilation is used for next modeling constantly as the River Basin Hydrology state variable accurate estimated value on the same day in the substitution model.Result such as Figure 12.

Figure 13 has provided the top layer 10cm soil moisture result that annual all assimilations in 2006 obtain through model, remote sensing appraising and assimilation constantly.The fate of horizontal ordinate representative assimilation among the figure is arranged according to the chronological order that remote sensing image obtains, and is that threshold value has been rejected the soil moisture extreme value that the minority abnormity point in the remote sensing data causes with 0 and 0.60.Can find out obviously that the assimilation effect of 2 lattice points all very significantly.See from trend; The assimilation result of soil moisture is very high with desirable " true value " degree of agreement; And main difference partly is; When to be that the soil moisture value is too high exceed normal scope in the position that the ideal value of remote sensing observations is undergone mutation; The assimilation result is near the modeling value; The characteristic that has shown algorithm filtering, the relatively stable of topsoil humidity assimilation result kept in the better filtering remote sensing of this characteristic noise effect.Above presentation of results PF is used for topsoil humidity assimilation test has reliability preferably.

Comparatively speaking, the result of remote sensing appraising fluctuation is bigger.And the modeling result compares with remote sensing appraising, and the trend of obviously underestimating is arranged, and the fluctuation of simulated soil humidity variation is comparatively mild, shows simulation humidity line at 70-80 constantly, and 110-130 is June constantly, seems very level and smooth in October and November.Its main cause is, all is in the rare time of precipitation than short-term before reaching during this period of time, and there be over-evaluating to a certain degree in model in the dehydration simulation of topsoil, and its origin cause of formation is analyzed in the humidity simulation precision evaluation of a last chapter; Simultaneously, because when soil moisture was reduced to certain scope, the fall off rate of soil moisture will slow down significantly.It should be noted that most of lattice point (seeing Figure 13) in 81-91 assimilation constantly promptly by the end of June to early August, i.e. basin precipitation more abundant period, modeling and remote sensing observations result's numerical value and trend show all quite similar.Explanation is during precipitation, and the humidity accuracy of simulation is higher than the no precipitation phase possibly.

Figure 14 provided the same day observation data with assimilation result's correlativity relatively, table 5 provided the quantified precision of assimilation effect.

From Figure 14, find out; Though the soil water content of modeling and actual observation has correlativity preferably; To a certain degree underestimate but exist; Occupy its top like 1: 1 diagonal line in the some slip chart that looses; And the result after the particle filter assimilation has obviously suppressed this trend of underestimating, and has positive correlation preferably with the soil moisture of surveying.After the assimilation, the soil moisture precision is greatly improved, and root-mean-square deviation (RMSE) is by 0.081m ³/ m ³Reduce to 0.045m ³/ m ³, mean deviation (MAE) is then by 0.073m ³/ m ³Reduce to 0.036m ³/ m ³

Table 5: simulation and the statistical appraisal index (DOY:2006148) of assimilating the result

Figure 15 has provided the comparison directly perceived of observation, assimilation and analog result.Can find out that assimilation has significantly improved the simulate effect of model, the result more levels off to observed reading.

Claims (9)

1. basin yardstick soil moisture remotely-sensed data assimilation method may further comprise the steps:

A) data are prepared;

B) make up basin soil moisture assimilation observation operator;

C) make up hydrological distribution model assimilation platform;

D) structure is based on the basin soil moisture remotely-sensed data assimilation scheme of hydrological distribution model and particle filter assimilation algorithm.

2. basin according to claim 1 yardstick soil moisture remotely-sensed data assimilation method; It is characterized in that steps A) in, said data prepare to comprise the underlying surface supplemental characteristic, satellite image data of the hydrology, weather data, dem data, vegetation, soil parameters storehouse, soil, basin utilization/cover and soil types and runoff and soil moisture measured data day by day.

3. basin according to claim 2 yardstick soil moisture remotely-sensed data assimilation method is characterized in that step B) in, adopt MODIS image and active microwave ASAR active microwave data to make up basin soil moisture assimilation observation operator.

4. basin according to claim 3 yardstick soil moisture remotely-sensed data assimilation method is characterized in that step C) in, on the frame foundation of hydrological model, make up basin soil moisture remotely-sensed data assimilation platform.

5. basin according to claim 4 yardstick soil moisture remotely-sensed data assimilation method; It is characterized in that step C) in; In the soil moisture transmission; Considered that Richard ' the s equation of gravity and absorption affinity is used to describe vertical unsaturated soil moisture movement of soil profile, numerical solution adopts method of finite difference;

Richard ' s is based on following two equation inferences, unsaturated soil water flux law Buckingham-Darcy equation (formula 1) and continuity equation (formula 2):

The Buckingham-Darcy equation:

$q=-K\left(h\right)(\frac{\&PartialD;h}{\&PartialD;z}+1)---\left(1\right)$

Wherein q is unsaturated water flux density, and h is a holard matric potential, and K (h) is a hydraulic conductivity, and z is the degree of depth of vertical direction;

Continuity equation:

$\frac{\&PartialD;\mathrm{\&theta;}}{\&PartialD;t}+\frac{\&PartialD;q}{\&PartialD;z}=0---\left(2\right)$

Based on the hypothesis of vertical flow, the q of formula 2 is replaced by formula 1, tries to achieve to separate to be one dimension Richard ' s equation expression formula:

$\frac{\&PartialD;\mathrm{\&theta;}}{\&PartialD;t}=\frac{\&PartialD;}{\&PartialD;z}\left(K\left(h\right)\frac{\&PartialD;h}{\&PartialD;z}\right)+\frac{\&PartialD;K\left(h\right)}{\&PartialD;z}---\left(3\right)$

Wherein t is the time, and θ is the soil volumetric(al) moisture content.

6. according to claim 4 or 5 described basin yardstick soil moisture remotely-sensed data assimilation methods; It is characterized in that step C) in; Is to be no more than 10 layers soil horizon with the soil cylinder in basin according to identical layering criteria; Maximum depth of soil is 3.43 meters, thereby the soil cylinder on the lattice point is because the difference of thickness of soil has the different soil numbers of plies; On the basis that obtains the soil layering parameter, utilize the average milliosmolarity down that obtains soil surface to the quantity of precipitation and the runoff yield difference on ground, based on unsaturated soil hydrodynamic force equation, simulate vertical moisture transmission of the moisture state and the interlayer of each layer soil.

7. basin according to claim 6 yardstick soil moisture remotely-sensed data assimilation method; It is characterized in that step D) in; The distributed basin hydrological simulation platform of in utilizing step C, setting up is on the soil water space-time dynamic value analog basis; By particle filter order assimilation algorithm; Merge the surface soil humidity information that satellite remote sensing is obtained; Comprise the surface humidity result that active microwave ENVISAT-ASAR and the inverting of MODIS visible image obtain; And the error of considering simulation and remote-sensing inversion is assimilated calculating with two kinds of information sources; Then the soil moisture data of upgrading are fed back in the distributed simulation platform, constantly simulation obtains having the conforming basin of physics topsoil humidity data collection.

8. basin according to claim 7 yardstick soil moisture remotely-sensed data assimilation method; It is characterized in that step D) in; The calculation process of particle filter assimilation algorithm is following: a, the soil moisture that the simulation of model is obtained are carried out disturbance, produce population, promptly from q (x _k| x _0:k-1, z _1:k) in randomly draw N limited sample, N is a particle number; B, the obtaining constantly of remotely-sensed data, the soil moisture result who uses remote sensing appraising to obtain, based on

Calculate the weight of corresponding particle; C, the weight of particle is done normalization handle; D, resampling utilize importance sampling SIR method that the particle weight is carried out resampling; E, after weight sampling, model state is carried out last update calculation, the soil moisture result after obtaining assimilating.

9. basin according to claim 8 yardstick soil moisture remotely-sensed data assimilation method; It is characterized in that step D) in, utilize the soil moisture estimation result of earth observation satellite TERRA-MODIS and ENVISAT-ASAR that the topsoil humidity variables in the hydrological model is carried out four-dimensional data assimilation.

Images